Reducing Oscillations for Obstacle Avoidance in a Dense Environment Using Deep Reinforcement Learning and Time-Derivative of an Artificial Potential Field

Abstract

Obstacle avoidance plays a crucial role in ensuring the safe path planning of quadrotor unmanned aerial vehicles (QUAVs). In this study, we propose a hierarchical framework for obstacle avoidance, which combines the use of artificial potential field (APF) and deep reinforcement learning (DRL) for training low-level motion controllers. Unlike traditional potential field methods, our approach modifies the state information received by the motion controllers using the outputs of the APF path planner. Specifically, the assumed target position is pushed away from obstacles, resulting in adjustments to the perceived position errors. Additionally, we address path oscillations by incorporating the target’s velocity information, which is calculated based on the time-derivative of the repulsive force. Experimental results have validated the effectiveness of our proposed framework in avoiding collisions with obstacles and reducing oscillations.

1. Introduction

Due to their high flexibility, quadrotor unmanned aerial vehicles (QUAVs) have gained significant popularity in various applications, including parcel delivery [1], precision agriculture [2,3,4], search and rescue [5,6], and surveillance [7]. In these scenarios, the QUAV is typically required to autonomously navigate to a target position. However, in dense environments, unexpected obstacles can obstruct the path and lead to collisions. This task is even more challenging in a multi-agent system [8,9]. Therefore, obstacle avoidance becomes a crucial task for ensuring safe path planning. The QUAV must find an unobstructed path to the target while considering its physical limitations. Conventional approaches to obstacle-free path planning include Dijkstra, probabilistic roadmap (PRM), rapidly-exploring random trees (RRT), and artificial potential field (APF) [10,11,12,13]. Since these approaches do not assume the physical realization of the agent to be of a certain type, they are applicable to QUAV navigation. APF faces a critical challenge in dense environments. When multiple obstacles narrow the passageway to the target, the agent may experience oscillations due to the repulsive force created by the obstacles that conflicts with the attractive force towards the target. As a result, reaching the target smoothly becomes difficult.

Most existing approaches to reducing APF oscillations either rely on second-order optimization theory techniques or employ an escaping mechanism when the agent detects oscillations. However, these approaches are not fully compatible with a hierarchical path planning framework because they directly modify the forces applied to the QUAV without considering its often nonlinear dynamics. Additionally, they commonly attribute oscillations to local minima, which is not always the case. One possible solution to this problem is to utilize global obstacle information to bend an imaginary rubber band path [14]. However, this approach poses challenges for real-time decision-making, as it requires prior knowledge of obstacle positions.

Based on the symmetrical motion controller design in [15], our previous work [16] enabled the utilization of an APF in a hierarchical framework, utilizing only local information provided by sensors mounted on the QUAV. In this framework, an APF path planner generates a state that incorporates modified position errors, while low-level motion controllers are responsible for both position and attitude control. In [17], we introduced the Hessian matrix of the APF as a damping term, which effectively dissipates system energy, and we proved its stability using an energy-based Lyapunov function. In this current work, we utilize the time-derivative of the APF to further reduce oscillations in QUAV path planning within a hierarchical control framework.

Considering that the QUAV dynamic system is nonlinear and strongly coupled, we use deep reinforcement learning (DRL) to train the low-level motion controllers [18]. DRL not only alleviates the labor-intensive process of parameter tuning, as compared to the traditional proportional–integral–derivative (PID) method, but also demonstrates remarkable fitting capability with the aid of deep neural networks. The key contributions of this study are outlined as follows.

• We incorporate the time-derivatives of the potential function to account for the velocities of intermediate target points, thereby resolving the conflict between the fixed-point chasing logic of the DRL motion controller and the potentially erratic movements of the intermediate target points.
• To enhance the control performance of a high-dimensional nonlinear system, we reconfigure the states of the DRL motion controller to eliminate any asymmetry and ensure stability.
• We conducted comparative simulation experiments to validate the effectiveness of the proposed method in reducing oscillations and preventing overshooting.
• The complex dynamics of the QUAV system are considered, which is different from other approaches that either do not discuss agent dynamics or simplify the dynamics as a simple second-order point mass model.

The structure of this paper is outlined as follows. In Section 2, we provide a comprehensive overview of the existing research on obstacle-free path planning. In Section 3, we introduce the necessary background knowledge, including the principles of conventional artificial potential fields, drone dynamics, and deep reinforcement learning. Section 4 elaborates on the methodology that we have employed in this study. We present the experiments conducted in Section 5 to validate the effectiveness of the proposed method in reducing oscillation and overshooting. Section 6 discusses limitations of this work and its applicability to dynamic and irregular environments. Finally, in Section 7, we conclude our findings.

2. Related Work

We identify three major approaches to achieve obstacle-free path planning. Since they do not set restrictions on the type of agents being used, these approaches are applicable to QUAV obstacle avoidance. The first approach is grid search, exemplified by algorithms such as Dijkstra, A*, and D* [10,19,20]. These algorithms are guaranteed to find a viable least-cost path when the environment is finite. However, they suffer from the curse of dimensionality, as the number of vertices to be explored grows exponentially with the dimensionality of the working space. Additionally, the need for higher control resolutions further increases the computational burden. The second approach is sampling-based search, with PRM [11] and RRT [12] being typical examples. For example, Farooq et al. [21] guided the QUAV to avoid dangerous zones in a dynamic environment by computing the PRM. Compared to grid search, these algorithms are capable of quickly producing paths in high-dimensional spaces. However, since sampling-based methods rely on global information, such as environment boundaries, they are not well-suited for real-time control where sensor data are used to generate actions for the next step. The third approach is APF, which assumes that both the target and obstacles generate potential fields that influence the movement of the agent. By leveraging the gradient of the summed potential field, the agent is attracted to the target while being repelled from nearby obstacles. For example, Ma’arif et al. [22] used APF to guide a single QUAV to reach the target and avoid collisions with dynamic obstacles. APF is widely used in QUAV navigation due to its ease of implementation and ability to guide the QUAV using only local information.

Despite its numerous advantages, APF has certain inherent limitations [23]. One recognizable limitation is that the agent can easily become trapped in local minima, where the attraction and repulsion forces cancel each other out. Algorithms designed to address this problem can generally be classified into three categories. The first category is local minimum removal. For example, Kim et al. [24] employed harmonic functions to construct potential fields, allowing for the selection of singularity locations and the elimination of local minima in free space. The second category is local minimum escape. For instance, Park et al. [25] combined APF with simulated annealing, which introduces randomness to the agent’s actions, enabling it to escape from local minima. Wang et al. [26] proposed the Left Turning scheme, which effectively handles U-shaped obstacles and helps the agent escape local minima. Lai et al. [27] proposed a dynamic step-size adjustment method to help the multi-UAV escape the local minima. The third category is local minimum avoidance. For instance, Doria et al. [28] utilized the deterministic annealing approach to expand the repulsive area and shrink the attractive area. This allows the agent to avoid local minima at the beginning, when the potential function is convex due to a high initial temperature. Additionally, Ge et al. [29] addressed a specific case of the local minimum problem known as goals non-reachable with obstacles nearby (GNRON). They considered the relative distance between the agent and the target when designing the repulsive potential function. As the agent approaches the target, this function approaches zero, thereby reducing the repulsive force in the target’s vicinity and overcoming the local minimum problem. Overall, these different approaches highlight the efforts made to overcome the limitations of APF and improve its performance in various scenarios.

APF faces another dilemma, which is the occurrence of oscillations in narrow passages with densely distributed obstacles nearby. Previous research on oscillation reduction is relatively limited and mainly draws inspiration from optimization theory techniques. For instance, Ren et al. [30] proposed the use of Levenberg’s modification of Newton’s method as a solution to oscillation problems. This approach incorporates second-order information by utilizing the Hessian matrix. Additionally, they adjusted the control law to maintain a constant speed, ensuring smooth movement of the agents. Biswas et al. [31] further compared first-order gradient descent methods with two second-order approaches and concluded that Levenberg–Marquardt is more effective in generating smoother trajectories and improving convergence speeds.

The second branch of oscillation reduction algorithms introduces the concept of virtual obstacles or targets. Similar to LME used in tackling local minima, methods belonging to this branch often employ escaping techniques once oscillations are detected. For instance, Zhao et al. [32] enhanced the manipulator’s predictive ability by incorporating dynamic virtual target points and utilized an extreme point jump-out function to escape oscillations. Zhang et al. [33] employed tangent APF to avoid local oscillations and introduced the back virtual obstacle setting strategy-APF algorithm, which enables the agent to return to previous steps and withdraw from concave obstacles. In a rule-based fashion, Zheng et al. [34] specified the condition for adding obstacles, compelling the resultant force to deflect when its angle to the obstacle center is too small. The dynamic step-size adjustment method proposed by Lai et al. [27] is also able to escape the jitter area where a local minimum is encountered, but it does not address oscillations in other cases, such as narrow passageways.

Oscillations can also be mitigated by modifying the repulsive forces in a certain manner. Tran et al. [35] estimated the projection of the repulsive force vector onto the attractive force vector and subtracted this component to prevent the agent from moving in the opposite direction of the attractive force. Martis et al. [36] introduced vortex potential fields to achieve seamless cooperative collision avoidance between mobile agents, and these were validated using Lyapunov stability analysis. For larger obstacles with irregular shapes, Szczepanski [37] combined the benefits of repulsive APF and vortex APF by defining multiple layers around the surface area of the obstacles. This approach surpassed traditional APF and pure vortex APF in terms of path smoothness.

There are also works specifically designed for QUAV path planning and obstacle avoidance. For example, Meradi et al. [38] proposed a nonlinear model predictive control method based on quaternions for QUAVs’ obstacle avoidance. Valencia et al. [39] constructed a QUAV obstacle detection and avoidance system using a monocular camera. Gageik et al. [40] discussed the use of complementary low-cost sensors in QUAV collision avoidance. However, the problem of oscillation reduction in QUAV obstacle-free path planning based on APF is largely under-explored.

With the advancement in DRL technology, motion controllers based on DRL have gained widespread usage in QUAV path planning and obstacle avoidance by virtue of their strong fitting capability. In this context, the APF can be seen as an upper layer that indirectly or directly influences the agent’s current pursuit position, while the low-level motion control task is handled by the DRL agent. For instance, the RL environment may employ convolutional neural networks to receive information about the surrounding potential energy and generate estimated rewards for various actions [41]. Xing et al. [42] combined the enhanced APF method with deep deterministic policy gradient to navigate autonomous underwater vehicles. In our previous research, we developed a hierarchical framework where the APF path planner was utilized to modify the position errors perceived by the DRL motion controllers, effectively altering the target position [16]. In this study, we further address oscillations by considering the velocity of the virtual target point, which has been validated in a collision-free formation control task [17]. Experimental results demonstrate that this improved algorithm reduces oscillations compared to approaches that do not incorporate second-order information, such as velocities.

3. Preliminaries

3.1. Drone Dynamics

The simulation environment is built based on gym-pybullet-drones [43]. In this work, a QUAV with an “X” configuration is considered. Similar to many fixed-wing UAVs [44], it has six degrees of freedom. Its position and Euler angles are defined as $\mathbf{p}={[x,y,z]}^T$ and $\mathbf{\Theta }={[\varphi ,\theta ,\psi ]}^T$, respectively. The Euler angles as shown in Figure 1 are obtained following the intrinsic z-y-x sequence (or, equivalently, extrinsic x-y-z sequence). Therefore, the rotation matrix from the body-fixed frame {${\mathbf{x}}_{\mathrm{b}},{\mathbf{y}}_{\mathrm{b}},{\mathbf{z}}_{\mathrm{b}}$} to the earth frame {${\mathbf{x}}_{\mathrm{e}},{\mathbf{y}}_{\mathrm{e}},{\mathbf{z}}_{\mathrm{e}}$} is calculated by

$$\mathbf{R}=[{\mathbf{x}}_{\mathrm{b}},{\mathbf{y}}_{\mathrm{b}},{\mathbf{z}}_{\mathrm{b}}]=\left[\begin{array}{ccc}{c}_{\theta }{c}_{\psi }& c_{\psi }{s}_{\theta }{s}_{\varphi }-s_{\psi }{c}_{\varphi }& c_{\psi }{s}_{\theta }{c}_{\varphi }+s_{\psi }{s}_{\varphi }\\ c_{\theta }{s}_{\psi }& s_{\psi }{s}_{\theta }{s}_{\varphi }+c_{\psi }{c}_{\varphi }& s_{\psi }{s}_{\theta }{c}_{\varphi }-c_{\psi }{s}_{\varphi }\\ -s_{\theta }& s_{\varphi }{c}_{\theta }& c_{\varphi }{c}_{\theta }\end{array}\right],$$

(1)

where the coordinates of ${\mathbf{x}}_{\mathrm{b}},{\mathbf{y}}_{\mathrm{b}}$, and ${\mathbf{z}}_{\mathrm{b}}$ are in the earth frame, ${\mathbf{z}}_{\mathrm{b}}={[z_{\mathrm{b}}^x,z_{\mathrm{b}}^y,z_{\mathrm{b}}^z]}^T$, $c^{\xi }$ and $s^{\xi }$$(\xi =\varphi ,\theta ,\psi )$ are short for $cos\xi$ and $sin\xi$, respectively, and $\mathbf{R}\in SO\left(3\right)$.

In this work, we use the frame of SE(3) control [45] to avoid using a transformation matrix $\mathbf{W}$ between the Euler angles $\mathbf{\Theta }$ and the angular velocity $\mathit{\omega }={[{\omega }^x,{\omega }^y,{\omega }^z]}^T$, which suffers from singularity as $\theta$ approaches $\frac{\pi }{2}$. This technique overcomes the small angle restriction. The hat map $\widehat{·}:{\mathbb{R}}^3\to \mathfrak{so}\left(3\right)$ is defined as

$$\widehat{\mathit{\alpha }}=\left[\begin{array}{ccc}0& -{\alpha }_3& {\alpha }_2\\ {\alpha }_3& 0& -{\alpha }_1\\ -{\alpha }_2& {\alpha }_1& 0\end{array}\right],$$

(2)

where $\mathit{\alpha }={[{\alpha }_1,{\alpha }_2,{\alpha }_3]}^T$ is a three-dimensional vector. Using the differential property of Lie Group, the QUAV’s dynamic system can be formally described by

$$\begin{array}{cc}\hfill \dot{\mathbf{p}}& =\mathbf{v},\hfill \\ \hfill m\dot{\mathbf{v}}& =\mathbf{RF}-m\mathbf{g},\hfill \\ \hfill \dot{\mathbf{R}}& =\mathbf{R}\widehat{\mathit{\omega }},\hfill \\ \hfill \mathbf{J}\dot{\mathit{\omega }}& =-\mathit{\omega }\times \mathbf{J}\mathit{\omega }+\mathit{\tau },\hfill \end{array}$$

(3)

where the velocity $\mathbf{v}={[v^x,v^y,v^z]}^T$ and the gravity $\mathbf{g}={[0,0,g]}^T$ are in the earth frame, external thrust $\mathbf{F}={[0,0,f]}^T$ and torque $\mathit{\tau }={[{\tau }^{\varphi },{\tau }^{\theta },{\tau }^{\psi }]}^T$ are in the body-fixed frame, m is the mass, and $\mathbf{J}=diag(I^{\varphi },I^{\theta },I^{\psi })$ is the inertia matrix.

The four rotors are controlled by the pulse width modulation (PWM) signal $\mathbf{u}$, and the rotation speed $\mathbf{\Omega }={[{\Omega }_1,{\Omega }_2,{\Omega }_3,{\Omega }_4]}^T$ is described by

$$T_{\mathrm{m}}\dot{\mathbf{\Omega }}=k_{\mathrm{m}}\mathbf{u}+{\mathbf{b}}_{\mathrm{m}}-\mathbf{\Omega },$$

(4)

where $T_{\mathrm{m}},k_{\mathrm{m}}$, and ${\mathbf{b}}_{\mathrm{m}}$ are the motors’ coefficients, which remain constant. The thrust $\mathbf{F}$ and torque $\mathit{\tau }$ are then generated by the four rotors according to

$$\left[\begin{array}{c}f\\ {\tau }^{\varphi }\\ {\tau }^{\theta }\\ {\tau }^{\psi }\end{array}\right]=\left[\begin{array}{cccc}{C}_T& C_T& C_T& C_T\\ -\frac{L}{\sqrt{2}}{C}_T& -\frac{L}{\sqrt{2}}{C}_T& \frac{L}{\sqrt{2}}{C}_T& \frac{L}{\sqrt{2}}{C}_T\\ -\frac{L}{\sqrt{2}}{C}_T& \frac{L}{\sqrt{2}}{C}_T& \frac{L}{\sqrt{2}}{C}_T& -\frac{L}{\sqrt{2}}{C}_T\\ -C_D& C_D& -C_D& C_D\end{array}\right]\left[\begin{array}{c}{\Omega }_1^2\\ {\Omega }_2^2\\ {\Omega }_3^2\\ {\Omega }_4^2\end{array}\right],$$

(5)

where L is the arm length and $C_T,C_D$ are the thrust and torque coefficients, respectively.

The environment contains N obstacles denoted by $\mathcal{O}=[1,2,\cdots ,N]$. The geometric relationship between the QUAV, the target, and obstacle i is shown in Figure 2, in which we define $\mathbf{p}={[x,y,z]}^T$ as the QUAV’s current position, ${\mathbf{p}}_{\mathrm{d}}$ as the target position, and ${\mathbf{o}}_i$ as the closest point of obstacle i to the QUAV. Moreover, ${\mathbf{e}}_i=\mathbf{p}-{\mathbf{o}}_i={[e_i^x,e_i^y,e_i^z]}^T$ and ${\mathbf{e}}_{\mathrm{d}}=\mathbf{p}-{\mathbf{p}}_{\mathrm{d}}={[e^x,e^y,e^z]}^T$ are error vectors to obstacle i and the target, respectively.

3.2. Classical APF

The concept of APF was initially introduced by Khatib [13] with the aim of distributing collision avoidance across different levels of control. This approach involves creating a potential field, similar to an electrostatic field, within the agent’s working space. In this field, the target location is represented as a valley, while obstacles are represented as peaks. By computing the negative gradient of the total potential function, the agent can determine a feasible path to the target while avoiding collisions with obstacles. From another perspective, the target exerts an attractive force on the agent, while the obstacles exert a repulsive force. To ensure the existence of these attractive and repulsive forces, the potential functions must be continuous and differentiable, except at the target position or within the obstacles. The calculation of the attractive and repulsive forces is performed by

$$\begin{array}{cc}\hfill {\mathbf{F}}_{\mathrm{att}}& \hfill =-{\nabla }_{\mathbf{p}}{U}_{\mathrm{att}}\left(\right||{\mathbf{e}}_{\mathrm{d}}\left|\right|),\hfill \\ \hfill {\mathbf{F}}_{{\mathrm{rep}}_i}& \hfill =-{\nabla }_{\mathbf{p}}{U}_{{\mathrm{rep}}_i}\left(\right||{\mathbf{e}}_i\left|\right|),\hfill \end{array}$$

(6)

where $U_{\mathrm{att}}(·)$ and $U_{{\mathrm{rep}}_i}(·)$ are potential functions that receive scalars instead of vectors. This maintains the isotropy of the overall system, since the norms of the APF forces are only determined by the radial distances and not by the relative angles to the attraction or repulsion sources.

Figure 3 depicts the attractive and repulsive forces in the APF algorithm. In this figure, the green star represents the target position, while the orange circles symbolize obstacles with specific regions of influence. The repulsive force generated by obstacle i is denoted as $F_{\mathrm{rep}}$, and the attractive force generated by the target is denoted as $F_{\mathrm{att}}$. The resultant force exerted on the agent is represented as $F_{\mathrm{r}}$. It is worth noting that APF can be applied to both global path planning and local path planning, as illustrated in Figure 3a and Figure 3b, respectively. In the global case, the obstacles have a maximum influence distance ${\rho }_0$, and their positions and shapes are known in advance. However, in the local case, which is employed in this study, the QUAV solely relies on the relative distances and velocities obtained from sensors to calculate the APF forces. The maximum sensing range is denoted as d.

3.3. Deep Reinforcement Learning

RL is a commonly adopted methodology in machine learning. An RL problem is usually characterized as a Markov decision process with the state initialized by a distribution $\mathcal{P}\left(s_1\right)$. In this work, we assume that the Markov property is satisfied. In other words, the transition probability for the next state is only determined by the current state and action, which is formalized as $\mathcal{P}\left(s_{t+1}\right|s_1,a_1,\dots ,s_t,a_t)=\mathcal{P}\left(s_{t+1}\right)|s_t,a_t)$. At time step t, the agent receives the current state $s_t$ and, depending on whether the policy is stochastic or deterministic, takes an action $a_t$ according to its policy ${\pi }_{\alpha }(·|s):\mathcal{S}\to \mathcal{P}\left(\mathcal{A}\right)$ or ${\mu }_{\alpha }\left(s\right):\mathcal{S}\to \mathcal{A}$, where $\alpha$ denotes policy parameters. The environment will then give a scalar reward $r_t\sim \mathcal{R}(s_t,a_t)$ and the next state $s_{t+1}\sim \mathcal{P}(s_t,a_t)$. The aim of the agent is to maximize the expected return $G_t={\sum }_{i=t}^{\infty }{\gamma }^{i-t}{r}_i$, where $\gamma \in (0,1)$ is the discount factor that prevents $G_t$ from approaching infinity and controls how myopic the agent is.

To train the motion controller, we employ twin-delayed deep deterministic policy gradient (TD3), which is a model-free actor-critic DRL algorithm featuring three distinct characteristics. Firstly, it utilizes two critic networks, $Q_{w_1}$ and $Q_{w_2}$, to tackle the issue of action value overestimation. Secondly, it incorporates smooth regularization to minimize the variance of the target update, which means the action is given by Equation (7), ensuring a more stable and reliable training process, which is crucial for achieving optimal performance:

$${\tilde{a}}_{t+1}={\mu }_{{\alpha }^{\prime }}\left(s_{t+1}\right)+ϵ,ϵ\sim \mathrm{clip}(\mathcal{N}(0,\tilde{\sigma }),-c,c),$$

(7)

where ${\mu }_{{\alpha }^{\prime }}$ is the target actor network, $ϵ$ is the Gaussian noise, $\tilde{\sigma }$ is the standard deviation, and c is the noise bound. The loss function for the i-th critic network is represented by

$$L_{\mathrm{TD}3}\left(w_i\right)={\mathbb{E}}_{s_t,a_t,r_t,s_{t+1}\sim D}[{(r_t+\gamma \underset{j=1,2}{min}{Q}_{w_j^{\prime }}(s_{t+1},{\tilde{a}}_{t+1})-Q_{w_i}(s_t,a_t))}^2],$$

(8)

where D is the replay buffer and $Q_{w_j^{\prime }}(j=1,2)$ represent the target critic networks.

The third feature of TD3 is that actor networks are updated only after a fixed number of updates d to the critic. The actor network’s updates follow

$${\nabla }_{\alpha }J\left(\alpha \right)={\mathbb{E}}_{a_t={\mu }_{\alpha }\left(s_t\right)}\left[{\nabla }_{\alpha }{\mu }_{\alpha }\left(s_t\right){\nabla }_{w_1}{Q}_{w_1}(s_t,a_t)\right].$$

(9)

4. Methodology

4.1. Motion Controller Training

As is shown in Figure 4, motion control is divided into two hierarchies: translational control and angular control. Due to the presence of gravity, vertical control differs inherently from horizontal control. Hence, we further subdivide translational control into vertical and horizontal control. As the three angular controllers adhere to the same control logic, and likewise the two horizontal controllers, only three networks are required. The actor networks’ parameters for angular, vertical, and horizontal control are denoted by ${\alpha }_{\mathrm{a}},{\alpha }_{\mathrm{v}},$ and ${\alpha }_{\mathrm{h}}$, respectively. Their states are designed as

$$\begin{array}{cc}\hfill s^{\xi }& ={[e^{\xi },{\omega }^{\xi }]}^T,(\xi =\varphi ,\theta ,\psi )\hfill \\ \hfill s^x& ={[e^x,v^x,\mathrm{sign}\left(e^x\right)e^z,\mathrm{sign}\left(e^x\right)v^z,\overline{\theta },\dot{\overline{\theta }}]}^T,\hfill \\ \hfill s^y& ={[e^y,v^y,\mathrm{sign}\left(e^y\right)e^z,\mathrm{sign}\left(e^y\right)v^z,\overline{\varphi },\dot{\overline{\varphi }}]}^T,\hfill \\ \hfill s^z& ={[e^z,v^z]}^T,\hfill \end{array}$$

(10)

in which $\overline{\theta }$ and $\overline{\varphi }$, variables that characterize the pitch and the roll, are defined as

$$\begin{array}{cc}\hfill \overline{\theta }& \hfill =arctan\frac{z_{\mathrm{b}}^x}{z_{\mathrm{b}}^z}\hfill \\ \hfill \overline{\varphi }& \hfill =arctan\frac{z_{\mathrm{b}}^y}{z_{\mathrm{b}}^z}.\hfill \end{array}$$

(11)

There are a few points to note here. First, $x_{\mathrm{d}},y_{\mathrm{d}},z_{\mathrm{d}}$, and ${\psi }_{\mathrm{d}}$ are given by external sources, while ${\varphi }_{\mathrm{d}}$ and ${\theta }_{\mathrm{d}}$ are calculated from the outputs of the translational controllers. Second, horizontal states not only consider the corresponding position and velocity errors but also take into account height and angular information. Third, angular errors $e_{\xi }$ are defined using the rotation matrix instead of the angular differences, and they are calculated by

$${\mathbf{e}}_R={[e^{\varphi },e^{\theta },e^{\psi }]}^T=\frac{1}{2}{({\mathbf{R}}_{\mathrm{d}}^T\mathbf{R}-{\mathbf{R}}^T{\mathbf{R}}_{\mathrm{d}})}^{\vee },$$

(12)

in which the vee map $\vee :\mathfrak{so}\left(3\right)\to {\mathbb{R}}^3$ is the inverse of the hat map, and ${\mathbf{R}}_{\mathrm{d}}$ is the desired rotation matrix uniquely determined by ${\varphi }_{\mathrm{d}},{\theta }_{\mathrm{d}}$, and ${\psi }_{\mathrm{d}}$.

To obtain ${\varphi }_{\mathrm{d}}$ and ${\theta }_{\mathrm{d}}$, the translational controllers’ networks first generate desired forces along each axis following

$$\begin{array}{cc}\hfill f_{\mathrm{d}}^x& =m{c}_{\mathrm{t}}{\mu }_{{\alpha }_{\mathrm{h}}}\left(s^x\right),\hfill \\ \hfill f_{\mathrm{d}}^y& =m{c}_{\mathrm{t}}{\mu }_{{\alpha }_{\mathrm{h}}}\left(s^y\right),\hfill \\ \hfill f_{\mathrm{d}}^z& =m(c_{\mathrm{t}}{\mu }_{{\alpha }_{\mathrm{v}}}\left(s^z\right)+g),\hfill \end{array}$$

(13)

where the networks’ outputs are bounded by $(-1,1)$, and $c_{\mathrm{t}}$ is the maximum translational acceleration. Then, the desired roll and pitch can be calculated as [15]

$$\begin{array}{cc}\hfill {\varphi }_{\mathrm{d}}& \hfill =arcsin\frac{s_{{\psi }_{\mathrm{d}}}{f}_{\mathrm{d}}^x-c_{{\psi }_{\mathrm{d}}}{f}_{\mathrm{d}}^y}{\left|\right|{[f_{\mathrm{d}}^x,f_{\mathrm{d}}^y,f_{\mathrm{d}}^z]}^T\left|\right|},\hfill \\ \hfill {\theta }_{\mathrm{d}}& \hfill =arctan\frac{c_{{\psi }_{\mathrm{d}}}{f}_{\mathrm{d}}^x+s_{{\psi }_{\mathrm{d}}}{f}_{\mathrm{d}}^y}{f_{\mathrm{d}}^z}.\hfill \end{array}$$

(14)

The ultimate desired thrust and torques are generated following Equation (15), where $c_{\mathrm{r}}$ is the maximum angular acceleration. We can find the desired steady-state rotor speed ${\mathbf{\Omega }}_{\mathrm{ss}}=k_{\mathrm{m}}\mathbf{u}+{\mathbf{b}}_{\mathrm{m}}$ by taking the inverse matrix in Equation (5), and the PWM signal $\mathbf{u}$ will eventually be obtained.

$$\begin{array}{cc}\hfill f_{\mathrm{d}}& =\frac{f_{\mathrm{d}}^z}{c_{\varphi }{c}_{\theta }},\hfill \\ \hfill {\tau }_{\mathrm{d}}^{\xi }& =I^{\xi }{c}_{\mathrm{r}}{\mu }_{{\alpha }_{\mathrm{a}}}\left(s^{\xi }\right).(\xi =\varphi ,\theta ,\psi )\hfill \end{array}$$

(15)

Equation (16) defines the reward function for the angular, vertical, and horizontal reinforcement learning environments. In this equation, $e_t^{\zeta }$ represents the error at the current time, and $\zeta$ indicates the specific controller to which the error is applied. The agent receives a reward when the error approaches zero and is penalized otherwise. One notable advantage of this reward design is the absence of hyperparameters that require manual fine-tuning.

$$r_t^{\zeta }=|e_t^{\zeta }|-|e_{t+1}^{\zeta }|$$

(16)

During the training process, the reinforcement learning environments receive the agent’s action $a_t^{\zeta }$ and output the next state $s_{t+1}^{\zeta }$ as well as the reward $r_t^{\zeta }$ based on the current state $s_t^{\zeta }$. In this work, the pitch and x-axis controllers are chosen as the prototypes for the angular and horizontal controllers, meaning that $\zeta =\theta ,x,z$ during training. The training follows a certain order. First, the roll and vertical controllers are trained independently. When the pitch controller is being trained, we set $f_{\mathrm{d}},{\tau }_{\mathrm{d}}^{\varphi }$, and ${\tau }_{\mathrm{d}}^{\psi }$ as zero, and ${\tau }_{\mathrm{d}}^{\theta }=I^{\theta }{c}_{\mathrm{r}}{\mu }_{{\alpha }_{\mathrm{a}}}\left(s^{\theta }\right)$, while when the vertical controller is being trained, we set $f_{\mathrm{d}}$ as $m(c_{\mathrm{t}}{\mu }_{{\alpha }_{\mathrm{h}}}\left(s^z\right)+g)/c_{\varphi }{c}_{\theta }$ and the desired torques as zero. The x-axis controller is trained only after ${\alpha }_{\mathrm{a}}$ and ${\alpha }_{\mathrm{h}}$ are obtained because $s^x$ contains height and pitch information, which is affected by the vertical and pitch controllers. No obstacles or external disturbances are considered during training.

4.2. State Transformation

The trained motion controllers have the ability to navigate the QUAV in the corresponding dimension. However, when the states in the testing phase are too large, it becomes challenging for the networks to generate reasonable outputs based on the limited training samples. This issue typically arises in translational control, where the position errors are not bounded, except by the simulator. Since it is not feasible to gather information for position errors of arbitrary magnitude within a finite training period, using the policy networks in a large or boundless environment can result in uncontrollable behavior, which is inherently dangerous. To mitigate the impact of large horizontal position errors, we employ a hard-clipping technique on the variables $e^x$ and $e^y$, with a clipping bound of $d_c$. Although no hard-clipping is applied to $e^z$, its range is constrained due to the presence of a floor and a ceiling in the working space of this study.

As observed in [18], the trained DRL controllers will decelerate to avoid overshooting, indicating that there exists a natural velocity upper bound $\tilde{v}$, since ${\mu }_{\alpha }\left(s\right)>0$ when the velocity is smaller than $\tilde{v}$ and ${\mu }_{\alpha }<0$ when it exceeds the bound. Considering that $\tilde{v}$ is sometimes too large to guarantee safety in the testing phase, we manually set an upper bound $\overline{v}$ to the velocity in the testing phase to prevent the QUAV from moving too quickly in a free space. To increase data efficiency, we multiply $e^{v^x}$ and $e^{v^y}$ by $\tilde{v}/\overline{v}$, which can map the horizontal velocity errors from $[0,\overline{v}]$ to $[0,\tilde{v}]$ and therefore force the actual velocity inputs to have a similar scale as the training samples.

In addition to the aforementioned issues, it is necessary to address the problem of anisotropic control performance. As depicted in Figure 5a, we denote the force exerted on the QUAV as ${\mathbf{F}}_1$ when it is located at position A and as ${\mathbf{F}}_2$ when it is at position B. Ideally, if positions A and B are equidistant from the target at a distance r, we expect the two force vectors to have identical magnitudes and angles relative to the corresponding position error vectors. This implies that the control performance should exhibit isotropy. However, the current motion controllers we have obtained do not guarantee isotropy due to two reasons: (1) the x-axis and y-axis controllers independently calculate the force components along each axis based on $s^x$ and $s^y$, and (2) the policy networks are non-linear. To preserve this property, we transform the horizontal coordinate system such that the force vector components are calculated radially and tangentially to the horizontal error vector ${\mathbf{e}}_{\mathrm{h}}={[e^x,e^y,0]}^T$. As shown in Figure 5b, the unit radial and tangential vectors are defined as

$$\begin{array}{ccc}\hfill \mathbf{h}& \hfill =\frac{{\mathbf{e}}_{\mathrm{h}}}{\left|\right|{\mathbf{e}}_{\mathrm{h}}\left|\right|},T\left(\mathbf{h}\right)\hfill & \hfill =\mathbf{h}\times {[0,0,-1]}^T.\end{array}$$

(17)

The horizontal states in the testing phase are so far modified as

$$\begin{array}{cc}\hfill s^{\mathrm{rad}}& ={[e^{\mathrm{rad}},(\tilde{v}/\overline{v}){\mathbf{v}}^T\mathbf{h},\mathrm{sign}\left(e^{\mathrm{rad}}\right)e^z,\mathrm{sign}\left(e^{\mathrm{rad}}\right)v^z,\tilde{\theta },\dot{\tilde{\theta }}]}^T,\hfill \\ \hfill s^{\tan }& ={[e^{\tan },(\tilde{v}/\overline{v}){\mathbf{v}}^{T}T\left(\mathbf{h}\right),\mathrm{sign}\left(e^{\tan }\right)e^z,\mathrm{sign}\left(e^{\tan }\right)v^z,\tilde{\varphi },\dot{\tilde{\varphi }}]}^T,\hfill \end{array}$$

(18)

in which $e^{\mathrm{rad}},e^{\tan },\tilde{\theta }$, and $\tilde{\varphi }$ are defined as

$$\begin{array}{cc}\hfill e^{\mathrm{rad}}& =\mathrm{clip}({\mathbf{e}}_{\mathrm{d}}^T\mathbf{h},-d_c,d_c),\hfill \\ \hfill e^{\tan }& =\mathrm{clip}({\mathbf{e}}_{\mathrm{d}}^{T}T\left(\mathbf{h}\right),-d_c,d_c),\hfill \\ \hfill \tilde{\theta }& =arctan\frac{{\mathbf{z}}_{\mathrm{b}}^T\mathbf{h}}{z_{\mathrm{b}}^z},\hfill \\ \hfill \tilde{\varphi }& =arctan\frac{{\mathbf{z}}_{\mathrm{b}}^{T}T\left(\mathbf{h}\right)}{z_{\mathrm{b}}^z}.\hfill \end{array}$$

(19)

Similarly, the horizontal forces defined in Equation (13) should be modified following

$${[f_{\mathrm{d}}^x,f_{\mathrm{d}}^y,0]}^T=m{c}_{\mathrm{t}}{\mu }_{{\alpha }_{\mathrm{h}}}\left(s^{\mathrm{rad}}\right)\mathbf{h}+m{c}_{\mathrm{t}}{\mu }_{{\alpha }_{\mathrm{h}}}\left(s^{\tan }\right)T\left(\mathbf{h}\right).$$

(20)

It should be noted that, up to now, $e^{\mathrm{rad}}=\sqrt{{\left(e^x\right)}^2+{\left(e^y\right)}^2}\ge 0$ and $e^{\tan }=0$, and there is no reason for the QUAV to move along $T\left(\mathbf{h}\right)$. This is because no obstacles are taken into consideration yet. As we shall see in the next section, they will be further modified by the repulsive potential fields.

4.3. APF for Network Controller

In a hierarchical framework, RL motion controllers are solely responsible for the navigation of the QUAV in open spaces and do not address higher-level tasks such as obstacle avoidance. In this study, we incorporate the use of APF to guide the RL motion controllers by influencing the state inputs, thereby indirectly altering the QUAV’s current pursuit position. Figure 6 provides an illustration of how the target is updated, with the saturated target being the direct result of the clipping operation along the radial direction. Similar to traditional APF methods, any obstacle that falls within the QUAV’s sensing range generates a repulsive force. However, instead of directly applying this force to the agent, it is utilized to further modify the target’s position. This is accomplished by providing the motion controllers with modified position errors in the following manner:

$$\begin{array}{cc}\hfill e^{\mathrm{rad}}& =\mathrm{clip}({\mathbf{e}}_{\mathrm{d}}^T\mathbf{h},0,d_c)-{\mathbf{F}}_{\mathrm{rep}}^T\mathbf{h},\hfill \\ \hfill e^{\tan }& =-{\mathbf{F}}_{\mathrm{rep}}^{T}T\left(\mathbf{h}\right).\hfill \end{array}$$

(21)

It is important to note that distance clipping is performed prior to the modification of errors by ${\mathbf{F}}_{\mathrm{rep}}$, as depicted in Figure 6a. This implies that the radial and tangential position errors, $e^r$ and $e^t$, are not strictly limited by the distance d. This is done to ensure that the repulsive forces exerted by nearby obstacles have a significant impact on the movements of the QUAV. Without this consideration, the actual influence of repulsion on the target would be insignificant, as shown in Figure 6b, particularly when the original target is located far away.

The newly-introduced repulsive forces are dynamic, as the geometric relationship between the QUAV and adjacent obstacles changes. This implies that the position of the target also changes dynamically. Therefore, it is not enough to solely consider the agent’s current velocity in the motion controllers’ states, as the target was assumed to be stationary during the training process. In this study, we calculate the derivative of the repulsive force with respect to time to determine the target’s velocity. Using the chain rule, we can obtain that

$${\mathbf{v}}^{\mathrm{tar}}={\dot{\mathbf{F}}}_{\mathrm{rep}}=-\sum _{i:d_i\le d}{\nabla }_{\mathbf{p}}^2{U}_{{\mathrm{rep}}_i}\left(\right||{\mathbf{e}}_i\left|\right|)\dot{\mathbf{p}}=-\sum _{i:d_i\le d}{\mathbf{H}}_i\mathbf{v},$$

(22)

in which ${\mathbf{v}}^{\mathrm{tar}}$ is the velocity of the dynamic target, $d_i=\sqrt{{\left(e_i^x\right)}^2+{\left(e_i^y\right)}^2}$ is the horizontal distance to the i-th obstacle, and ${\mathbf{H}}_i$ is the Hessian matrix for the i-th obstacle’s repulsive field. A few properties are to be satisfied when designing $U_{{\mathrm{rep}}_i}(·):\left(0,d\right]\to \mathbb{R}$:

• $U_{{\mathrm{rep}}_i}\left(x\right)$ should be monotonically decreasing with respect to x;
• $\left|\right|{\nabla }_{\mathbf{p}}{U}_{{\mathrm{rep}}_i}\left(x\right)\left|\right|$ should be monotonically decreasing with respect to x;
• ${lim}_{x\to 0^{+}}\left|\right|{\nabla }_{\mathbf{p}}{U}_{{\mathrm{rep}}_i}\left(x\right)\left|\right|=+\infty$ and $\left|\right|{\nabla }_{\mathbf{p}}{U}_{{\mathrm{rep}}_i}\left(d\right)=0\left|\right|$;
• $U_{{\mathrm{rep}}_i}\left(x\right)$ should be second-order differentiable.

The first requirement is to ensure that the force is repulsive. The second property states that the influence of the obstacle weakens as the distance increases. The third property implies that the repulsion force approaches infinity at the surface of the obstacle, while it becomes zero when the obstacle is at the boundary of the sensing range. The final requirement guarantees the existence of the Hessian matrix. Since the repulsive APF is only applied through its derivative forms, it is reasonable to directly design ${\mathbf{F}}_{{\mathrm{rep}}_i}$ as

$${\mathbf{F}}_{{\mathrm{rep}}_i}=\left\{\begin{array}{cc}\lambda \left(\frac{1}{d_i^2}-\frac{1}{d^2}\right)\frac{{\mathbf{e}}_i}{\left|\right|{\mathbf{e}}_i\left|\right|}& d_i\le d\\ \mathbf{0}& d_i>d\end{array}\right.,$$

(23)

where $\lambda$ is a manually set factor that scales the intensity of repulsion.

Considering that the calculation of ${\mathbf{H}}_i$ necessitates significant computational resources and is time-consuming, we approximate ${\mathbf{v}}^{\mathrm{tar}}$ by utilizing the time difference of ${\mathbf{F}}_{\mathrm{rep}}$. If the radial component of ${\mathbf{v}}^{\mathrm{tar}}$ points towards the target, we eliminate it, as in such cases, the obstacle poses no threat to the agent. To be more specific, we define a distinct target velocity in the radial direction using the following approach:

$${\mathbf{v}}^{\mathrm{rad}}=\left\{\begin{array}{cc}{\mathbf{v}}^{\mathrm{tar}}& {\mathbf{h}}^T{\mathbf{v}}^{\mathrm{tar}}>0\\ \mathbf{0}& {\mathbf{h}}^T{\mathbf{v}}^{\mathrm{tar}}\le 0\end{array}\right..$$

(24)

This technique prevents $e^r$ from exceeding the bound d and therefore avoids overshooting. Ultimately, the horizontal states in the testing phase are modified as

$$\begin{array}{cc}\hfill s^{\mathrm{rad}}& ={[e^{\mathrm{rad}},(\tilde{v}/\overline{v}){(\mathbf{v}-{\mathbf{v}}^{\mathrm{rad}})}^T\mathbf{h},\mathrm{sign}\left(e^{\mathrm{rad}}\right)e^z,\mathrm{sign}\left(e^{\mathrm{rad}}\right)v^z,\tilde{\theta },\dot{\tilde{\theta }}]}^T,\hfill \\ \hfill s^{\tan }& ={[e^{\tan },(\tilde{v}/\overline{v}){(\mathbf{v}-{\mathbf{v}}^{\mathrm{tar}})}^{T}T\left(\mathbf{h}\right),\mathrm{sign}\left(e^{\tan }\right)e^z,\mathrm{sign}\left(e^{\tan }\right)v^z,\tilde{\varphi },\dot{\tilde{\varphi }}]}^T.\hfill \end{array}$$

(25)

This design enables the seamless integration of RL motion control and APF obstacle avoidance, without requiring any modifications to the underlying structure of the pre-trained motion controllers. The modification of states only occurs during the testing phase. What sets this strategy apart from traditional APF approaches is that, instead of directly applying forces to the QUAV, they are utilized to propel the target. This approach is equivalent to generating a sequence of waypoints that are carefully designed to avoid any potential obstacles.

5. Simulation Experiments

5.1. Training Setup

The QUAV’s dynamics followed Crazyflie 2.0 [46]. All agents were trained using TD3, whose hyperparameters are summarized in

Table 1. TD3 hyperparameter settings.

Symbols	Descriptions	Values
|$\mathcal{D}$|	Replay Buffer Size	$1\times {10}^4$
N	Batch Size	512
$\gamma$	Discount Factor	0.99
$\tau$	Target Update Rate	0.1
$\delta$	Learning Rate	$1\times {10}^{-3}$
$\tilde{\sigma }$	Noise Variance	0.2
c	Noise Clip	0.5
d	Policy Update Delay	2

. The time step was set as 0.02 s to ensure control accuracy, while the upper bounds of translational and angular acceleration were $c_{\mathrm{t}}=5$ m/s² and $c_{\mathrm{r}}=20$ rad/s², respectively.

An MLP with two layers, each consisting of 32 nodes, was utilized to train each controller. The optimizer used in this process was Adam. The targets were randomly sampled from a uniform distribution. Specifically, angular targets were chosen from the range of $(-\frac{\pi }{2},\frac{\pi }{2})$, while translational targets were selected from the range of $[-5,5]$. It is worth noting that the horizontal network was trained only after the angular and vertical networks had been obtained. The training stage of each controller comprised a total of 200 training episodes, with episode lengths set at 100 and 500 for angular and translational control respectively. This ensured that the agent had sufficient time to reach the target before an episode concluded. Throughout the training phase, no obstacles or external disturbances were taken into consideration.

5.2. Testing Setup

We designed three different controllers for all experiments. Controller 1 utilized the default setting, while controller 2 excluded the axial target velocity cancellation operation. On the other hand, controller 3 solely employed APF to influence position errors, without affecting velocity errors.

The testing environment consisted of an open space with a floor and a ceiling that was 5 m high. Additionally, there were several cylindrical obstacles, each standing at a height of 5 m, strategically placed around the center of the space. The positions of these obstacles were determined by uniformly selecting the distances of their axes from the center within the range of $[0,6]$. Furthermore, the radii of the cylinders were uniformly chosen from the range of $[0.9,1.1]$. It was ensured that the minimum distance between any two obstacles was at least $d_{min}$. For each trial, three QUAVs were utilized, which acted as dynamic obstacles and generated APF repulsive forces similar to the fixed obstacles in the environment.

They were initialized at $[8sin\left(\frac{2k\pi }{3}\right),8cos\left(\frac{2k\pi }{3}\right)]$, and the corresponding targets were at $[-8sin\left(\frac{2k\pi }{3}\right),-8cos\left(\frac{2k\pi }{3}\right)]$, where $k=0,1,2$. As the agents moved, proximity sensors mounted on them could detect nearby obstacles. The maximum sensing range d was 2 m, the position error clipping bound $d_c$ was 0.5 m, while the repulsive force scaling factor $\lambda$ was set as 0.18. The manually set velocity upper bound $\overline{v}$ was set as 2 m/s, and through experiments, it was found that the natural velocity upper bound $\tilde{v}=1.5$ m/s.

The agents faced the challenge of reaching the targets on the opposite side without any collisions. The success of navigation was determined by whether each agent entered the 0.1 m proximity of its corresponding target. We created four different environments with specific values for $d_{min}$ and N: $(0.7,15),(0.7,12),(0.9,12)$, and $(0.9,9)$. In each environment, we conducted 1000 trials for every controller setting. Based on these trials, we calculated the success rates and the average time taken for successful trials.

To evaluate path smoothness, we calculated the curvatures of projection of the actual paths on the horizontal plane according to

$$\kappa =\frac{|\dot{x}\ddot{y}-\ddot{x}\dot{y}|}{{({\dot{x}}^2+{\dot{y}}^2)}^{\frac{3}{2}}},$$

(26)

in which the time-derivatives of positions are approximated by the differences between two time steps. It should be emphasized that $\kappa$ may be either positive or negative, depending on the turning direction. Path oscillations are reflected by the occurrence of sharp turns, which correspond to huge absolute curvatures. In other words, if a path contains multiple points with huge absolute curvatures, it is oscillatory. In this work, the curvature threshold for a sharp turn is defined as 1000.

5.3. Experiment Results

5.3.1. General Performance

Figure 7 is a vertical view of the agents’ trajectories as they moved around in the four environments, where the circles represent cylindrical obstacles, the curves in orange, green, and blue illustrate the trajectories, and the stars at their ends are the targets. A smaller $d_{min}$ indicates a more challenging scenario, since the minimum width of the passageway between two obstacles was smaller, while a bigger N means there were more obstacles around the origin. It can be seen that the hierarchical control framework we have proposed is capable of navigating the agents to the corresponding targets smoothly without collisions or detours, even in the densest environment.

To validate the effectiveness of the proposed modifications, we conducted 1000 trials for the three controllers in various environments. The success rates and average reaching times are presented in

Table 2. Success rate (SCR) and average reaching time (ART) of different controllers.

Settings	Controller 1		Controller 2		Controller 3
	SCR	ART (s)	SCR	ART (s)	SCR	ART (s)
(0.7, 15)	0.950	29.545	0.886	27.433	0.364	32.236
(0.7, 12)	0.968	27.194	0.945	25.223	0.675	29.012
(0.9, 12)	0.986	26.000	0.953	24.181	0.950	27.836
(0.9, 9)	0.989	24.529	0.975	22.801	0.968	25.984

. In the environment with low density (0.9, 9), all controllers achieved a high success rate. Among them, controller 2 demonstrated the fastest performance, with an average reaching time of only 22.801 s. This can be attributed to controller 2’s tendency to accelerate more dramatically when the repulsive force contributes to the attractive force. However, as the density increased, controller 2 experienced a more significant decrease in success rate compared to controller 1. This is because the higher acceleration increases the risk of overshooting, which will be discussed in the following section. Furthermore, controller 3, which did not consider any velocity errors, exhibited a severe decline in performance. It not only took more time to reach the targets but also struggled to handle dense environments. This was due to the delay between the assumed target and the actual target, modified by the repulsive APF forces. As a result, oscillatory behaviors were observed. In summary, the introduction of velocity errors effectively mitigated oscillations, while the cancellation of the target velocity radial component enhanced the agents’ ability to address overshooting and the resulting oscillations. However, this came at the expense of slower agent movement.

The method is also applicable to an environment with more QUAVs, as shown in Figure 8. The environment setting is similar to that discussed in the testing setup. In general, the agents are able to reach their targets without collisions. For example, in Figure 8b, where three pairs of agents move head to head, they can still reach the targets safely.

5.3.2. Oscillation Reduction and Overshooting Prevention

Table 3. The average number of sharp turns ($\left|\kappa \right|>1000$) of different controllers.

Settings	Controller 1	Controller 2	Controller 3
(0.7, 15)	79.427	73.855	122.403
(0.7, 12)	54.788	51.731	90.938
(0.9, 12)	41.925	33.658	75.417
(0.9, 9)	30.434	25.136	53.372

lists the average number of points with an absolute curvature greater than 1000 for the three controllers over all successful trials. It can be seen that controller 3 demonstrates the worst performance in terms of path smoothness, which testifies that target velocity information plays an essential role in reducing path oscillations. On the other hand, controllers 1 and 2 are close, with the path of the former being slightly more oscillatory. It should be noted that overshooting is not equivalent to global path smoothness, since it only occurs at certain areas, such as the triangular trap. Therefore, the oscillations caused by overshooting cannot easily be reflected by path curvature.

To give quantitative analysis to the effectiveness of radial target velocity cancellation operation as defined by Equation (24), another experiment setting was used. Concretely, three cylindrical obstacles with a radius of 1.5 m were randomly placed on a ring whose inner and outer radii were 2.3 m and 2.5 m, respectively. The three QUAVs were initialized at $[0.6sin\left(\frac{2k\pi }{3}\right),0.6cos\left(\frac{2k\pi }{3}\right)]$, and the corresponding targets were at $[-8sin\left(\frac{2k\pi }{3}\right),-8cos\left(\frac{2k\pi }{3}\right)]$, where $k=0,1,2$. In this way, they were trapped by the encircled obstacles at first, and the task was to escape the triangular trap. Ideally, the maximum distance between any QUAV and the corresponding target position is 8.6 m, otherwise it is caused by overshooting. In this experiment, overshooting cases are defined as cases where at least one QUAV’s distance to its target exceeded 10 m, as shown in Figure 9a. Furthermore, we define trapped cases where at least one QUAV’s distance to the origin is within 4 m, as shown in Figure 9b.

Table 4. The number of overshooting and trapped cases of different controllers.

	Controller 1	Controller 2	Controller 3
Overshooting Cases	23	145	481
Trapped Cases	24	199	641
SCR	0.969	0.732	0.237
ART (s)	13.735	13.958	18.566

shows the results. It can be seen that controller 1 encountered the fewest overshooting cases, which shows that the velocity cancellation operation helps preventing overshooting. Moreover, since this operation simplifies the obstacle environment by neglecting obstacles behind it, trapped cases are largely avoided as well.

An intuitive comparison of the three controllers can be found in Figure 10. In Figure 10b, the agent represented by the green trajectory experienced oscillation within the triangular area formed by obstacles 1, 2, and 3. The reason behind this oscillation was that the combined force exerted by obstacle 1 and obstacle 2 was directed towards the target and was not counteracted even though they did not pose a threat to the agent. Consequently, the force towards the goal became excessively strong, leading to the agent overshooting in this direction. Furthermore, disregarding velocity information entirely resulted in more frequent oscillations. The reason was that, in this case, the agent assumed the intermediate target to be fixed, while in reality the intermediate target was affected by the relative position between the agent and the obstacles. Oscillations occurred when the intermediate target updated before the agent reached the previous one. This can be observed in Figure 10c, where all agents’ trajectories exhibit a zigzag pattern, especially when maneuvering between obstacles. This increases the risk of collision with nearby obstacles.

5.3.3. The Effect of a Clipping Bound

In the mentioned experiments, we utilized a position error clipping bound of $d_c=0.5$ m. This choice was made based on previous experiments that indicated that it was the safest option. We varied the value of $d_c$ for controller 1, ranging from 0.3 m to 1.0 m with an increment of 0.1 m. A total of 1000 trials were conducted with the environment parameters set at (0.7, 12). Additionally, we observed that, for $d_c$ values of 0.3 m and 0.4 m, the agent’s velocity was relatively low, resulting in a longer time to reach the target. To account for more trials that had the potential to reach the targets, we relaxed the time requirement for success. This adjustment allowed us to consider a broader range of trials. The results of these experiments are presented in

Table 5. Success rate (SCR) and average reaching time (ART) of different controllers.

$\mathit{d}\mathit{c}$	0.3 *	0.4 *	0.5	0.6	0.7	0.8	0.9	1.0
SCR	0.517	0.930	$\mathbf{0}.\mathbf{968}$	0.941	0.851	0.714	0.585	0.568
ART (s)	43.437	34.983	27.194	22.340	19.606	18.705	19.146	$\mathbf{17}.\mathbf{710}$

* The time requirement for success was relaxed from 3000 steps to 5000 steps.

. From the table, it can be observed that the average reaching time is inversely correlated with $d_c$. Furthermore, $d_c=0.5$ yielded the highest success rate. This can be attributed to the fact that, when the agent moved too slowly, it lacked the flexibility to handle imminent collisions. Conversely, when the agent moved too quickly, the outputs of the low-level flight controllers were not sufficiently strong to decelerate the agent in time, thereby increasing the likelihood of collisions.

5.3.4. Comparison to Existing Methods

This section compares the proposed method to two other existing methods in terms of success rate, average reaching time, average absolute curvature, and the number of sharp turns. The first method to be compared uses A* as the path planner and the SE(3) control law proposed in [45] as the motion controller. This method discretizes the environment into cubes with side length 0.2 m, which are also known as voxels. It should be noted that the environment information of this algorithm is known in advance instead of being obtained with sensors. The second method to be compared [47] is a hierarchical DRL method that uses a DQN network to generate waypoints based on the sensor information returned by 16 rangefinders. The waypoints then serve as the target points of the motion controller trained with PPO. The testing environment was a single-angle scenario, where the agent was initialized at [0, 4.5, 2.5], and the target position was [0, −4.5, 2.5]. Twelve cylindrical obstacles with a radius of 0.4 m were randomly placed around the origin. The positions of these obstacles were determined by uniformly selecting the distances of their axes from the center within the range of $[0,3]$ and $d_{min}=0.6$. The results are listed in

Table 6. Comparison of the proposed method to existing methods.

	APF + TD3	A* + SE(3)	DQN + PPO	Controller 3
SCR	$\mathbf{0}.\mathbf{983}$	0.969	0.885	0.761
ART (s)	$\mathbf{15}.\mathbf{312}$	19.458	36.694	18.249
$\left|\kappa \right|$	7.121	$\mathbf{4}.\mathbf{794}$	6.610	11.426
sharp turns	0.651	0.480	$\mathbf{0}.\mathbf{431}$	1.112

, which also includes the performance of controller 3 in this environment. We can see that the proposed method outperformed the other two algorithms in terms of success rate and average reaching time, with the path slightly more oscillatory. Although path oscillation is an inherent problem with APF, it can be reduced to a large extent with the help of velocity information. After addressing this issue, the controller can navigate the QUAV in a more safe and efficient manner.

Figure 11 illustrates the paths of the four methods in the same environment. It can be seen that the proposed method’s path is smooth on a macro level and is far away from obstacles. The path of A* + SE(3) is mostly straight because the waypoints returned by A* are exactly on the centers of the voxels. Figure 11c shows a similar path, in that there are only 16 possible directions to choose from when updating the waypoints. The path of controller 3, on the other hand, is the most oscillatory one.

6. Discussion

6.1. Limitations

Although the designed controller can achieve high success rates with moderate reaching time even in a dense environment, there is still room for improvement. We considered cases as failures when the agents were close to the targets but still outside the 0.1 m neighborhood. This issue is referred to as the GNRON problem, which is discussed in Section 2.

Table 7. The number of GNRON cases for different controllers.

Environments	(0.7, 15)	(0.7, 12)	(0.9, 12)	(0.9, 9)
Controller 1	19	12	10	8
Controller 2	20	10	11	3
Controller 3	6	4	4	1

provides the number of GNRON cases for different controllers, where the success requirement is relaxed to entering the 0.3 m neighborhood of targets. It is evident that the success rates could have been higher if we had employed a strategy that actively addressed this problem. By introducing a damping term in the repulsive forces when the agents are in close proximity to the targets, this problem can be effectively resolved. Additionally, there are other limitations to the proposed method, such as the absence of explicit collision detection, leading to occasional collisions, and the hyperparameter $d_c$, which needs to be manually adjusted.

Another limitation of the proposed method is that it fails to escape local minima created by U-shaped obstacles, which is shown in Figure 12a. Moreover, when the QUAV’s moving direction is almost collinear with the joint force, the path is usually not smooth, as shown in Figure 12b. One possible solution is to combine the proposed method with the vortex potential field method discussed in [36,37], since in this way the agent will move smoothly around the obstacle instead of heading directly towards it and spending a long time adjusting the path.

6.2. Applicability to Dynamic Environments

The proposed method has the potential to be transferred to dynamic environments since collision avoidance behavior is only determined by the repulsive forces updated at every time step. In fact, the preceding experiments have verified this ability because all the QUAVs are dynamic, generating repulsive forces to other QUAVs. As shown in Figure 13, when they move close to each other, the QUAVs are able to take actions to avoid collisions. This behavior resembles the cooperative collision avoidance mentioned in [36]. One significant difference between agents and obstacles is that the former ones usually actively take action to avoid collisions while the latter ones’ behavior may be unpredictable. Therefore, avoiding collisions between moving QUAVs and moving obstacles is more challenging. One possible solution is to combine the current framework with reciprocal velocity obstacles [48] based on the current velocity information of the obstacles.

6.3. Irregular Obstacles

In real applications, the shapes of obstacles may not be as simple as cylinders. Since this work assumes that the QUAV only detects nearby obstacles with sensors, it is impractical to generate a global potential field in advance. There are two simple ways to adapt this framework so as to make it applicable to environments with irregular obstacles. The first method, as shown in Figure 14a, is to replace the surface of the sensed part with smaller cylinders. In this way, the irregular obstacle’s repulsive force equals the joint repulsive force of all cylindrical obstacles. The second method is to assume the part behind the sensed surface as solid and calculate its center of gravity. The distance to the obstacle is obtained as shown in Figure 14b. Some works deal with irregular obstacles in advanced ways. For example, Ge et al. [49] integrated the repulsive force of every part of the obstacle to calculate the ultimate repulsive force. This can model the obstacles more accurately, but it requires the QUAV to have prior global information of the environment. Guo et al. [50] used the closest point on the obstacle’s surface to generate repulsive force, which is the same as our method when the obstacles are cylindrical or round. However, separate discussions are needed for obstacles with corners or edge points. This makes it difficult to achieve stable control using the Hessian matrix considering that its calculation will be even more complex under such circumstances.

7. Conclusions

In this study, we have proposed an efficient hierarchical controller that combines the use of an APF for path planning and obstacle avoidance with DRL motion controllers to generate rotor signals. Our approach differs from traditional APF methods in that the APF force is no longer directly applied to the QUAV but instead used to modify the position errors, effectively pushing the sensed target position. Additionally, we have introduced the calculation of the time-derivative of the APF force and incorporated a second-order term in the velocity error of the DRL controllers’ inputs. This modification has successfully reduced oscillations that occur when an agent is surrounded by multiple obstacles. The oscillation technique not only ensures smoother paths but also significantly improves the robustness of the controller. Furthermore, we have demonstrated that canceling the projection of the target velocity on the attractive force when they are in the same direction is an effective method to prevent overshooting and enhance the success rate. Although this approach may slightly decrease the speed at which the target is reached, it is a reasonable compromise that prioritizes safety, especially in dense environments. We would like to emphasize the importance of the clipping factor $d_c$ as a hyperparameter in our design. It is crucial to find an optimal value for $d_c$, as both excessively small and excessively large values can lead to significant performance degradation. Lastly, we discussed the limitations of the proposed method and provided some possible solutions. The control framework in this work can be applied to a real testing environment considering that it is based on real parameters of a QUAV. Potential challenges include the accumulated errors of the discrete system that generates outputs based on differences, a key point that needs to be addressed in future work.