# Singularity-free Guiding Vector Field for
Robot Navigation

Weijia Yao, Héctor Garcia de Marina, Bohuan Lin, Ming Cao Weijia Yao and Ming
Cao are with Institute of Engineering and Technology (ENTEG), Bohuan Lin is
with Bernoulli Institute of Mathematics, Computer Science and Artificial
Intelligence (BI), all in the University of Groningen, the Netherlands. Hector
Garcia de Marina is a Researcher in the Department of Computer Architecture
and Automatic Control at the Faculty of Physics, Universidad Complutense de
Madrid, 28040 Madrid, Spain. {w.yao, b.lin<EMAIL_ADDRESS><EMAIL_ADDRESS>
###### Abstract

In robot navigation tasks, such as UAV highway traffic monitoring, it is
important for a mobile robot to follow a specified desired path. However, most
of the existing path-following navigation algorithms cannot guarantee global
convergence to desired paths or enable following self-intersected desired
paths due to the existence of singular points where navigation algorithms
return unreliable or even no solutions. One typical example arises in _vector-
field guided path-following (VF-PF) navigation algorithms_. These algorithms
are based on a vector field, and the singular points are exactly where the
vector field diminishes. Conventional VF-PF algorithms generate a vector field
of the same dimensions as those of the space where the desired path lives. In
this paper, we show that it is mathematically _impossible_ for conventional
VF-PF algorithms to achieve global convergence to desired paths that are
_self-intersected_ or even just _simple closed_ (precisely, homeomorphic to
the unit circle). Motivated by this new impossibility result, we propose a
novel method to transform self-intersected or simple closed desired paths to
_non-self-intersected_ and _unbounded_ (precisely, homeomorphic to the real
line) counterparts in a _higher-dimensional_ space. Corresponding to this new
desired path, we construct a singularity-free guiding vector field on a
higher-dimensional space. The integral curves of this new guiding vector field
is thus exploited to enable global convergence to the higher-dimensional
desired path, and therefore the projection of the integral curves on a lower-
dimensional subspace converge to the physical (lower-dimensional) desired
path. Rigorous theoretical analysis is carried out for the theoretical results
using dynamical systems theory. In addition, we show both by theoretical
analysis and numerical simulations that our proposed method is an extension
combining conventional VF-PF algorithms and trajectory tracking algorithms.
Finally, to show the practical value of our proposed approach for complex
engineering systems, we conduct outdoor experiments with a fixed-wing airplane
in windy environment to follow both 2D and 3D desired paths.

##  I Introduction

Several robot navigation tasks, such as highway traffic monitoring, underwater
pipeline inspection, border patrolling and area coverage, require the
fundamental functionality of following a desired path [1], and new
applications are emerging, such as using drones to probe atmospheric phenomena
along prescribed paths [2]. The path-following navigation problem has
attracted the attention from both the robotics community [3, 4, 5, 6, 7] and
the control community [8, 9, 10, 11, 12, 13]. In a path-following navigation
algorithm, the desired path is usually given in the form of a single connected
curve _without temporal information_ , and then robots are guided to converge
to and move along it with sufficient accuracy. Treating the desired path as a
_geometric object_ rather than a time-dependent point, path-following
navigation algorithms sometimes are able to overcome a number of performance
limitations rooted in trajectory tracking [14], e.g., inaccuracy due to
unstable zero dynamics [15] and difficulty to maintain constant tracking speed
[16]. Moreover, comparing to trajectory tracking algorithms, there is separate
interest for the study of path-following navigation algorithms since they are
more suitable for some applications, such as fixed-wing aircraft guidance and
control [4, 17, 18].

Among different path-following navigation algorithms, those using a _guiding
vector field_ have been studied widely [3, 4, 5, 8, 9, 10], and we refer to
these algorithms as _vector-field guided path-following navigation algorithms
(or, VF-PF algorithms)_. The _guidance_ feature of the vector field is
justified as follows: usually robot _kinematics models_ (e.g., single-
integrator and double-integrator models [17, 3]) are considered, and the
guiding vector field, as its name suggests, provides _guidance_ signal inputs
to the model. This is valid based on the common assumption that the _robot-
specific_ inner-loop _dynamics control_ can track these guidance signal inputs
effectively [19, 20, 16]. Thus, one can simply focus on the guidance layer
(i.e., designing a guiding vector field), while considering other control
layers separately. Specifically, the guiding vector field is designed such
that its integral curves are guaranteed to converge to a predefined geometric
desired path. Utilizing the convergence property of the vector field, one can
then derive suitable control laws. It is claimed in [18] that VF-PF algorithms
demonstrate the lowest cross-track error while requiring the least control
effort among several other path-following navigation algorithms. In addition,
[21] shows that VF-PF algorithms achieve better path-following accuracy than
the integral line-of-sight (ILOS) method [12].

Although the VF-PF algorithms are intuitive and easy to implement, the
rigorous analysis remains nontrivial for general desired paths [16, 22, 23,
13]. Significant difficulty in the analysis and application of the VF-PF
algorithms arises when there are singular points111A point where a vector
field becomes zero is called a _singular point_ of the vector field [24, p.
219].  in the vector field (see Fig. 1 and Fig. 2). In such a case, the
convergence of trajectories to the desired path cannot be guaranteed globally,
and the normalization of the vector field at those points is not well-defined
[16, 22, 3, 13]. In [3], it is assumed that these singular points are
repulsive to simplify the analysis, while this assumption is dropped in [16]
for a planar desired path and in [22, 13] for a desired path in 3D. However,
to the best of our knowledge, few efforts have been made on dealing with
singular points directly or on eliminating them effectively. Recently, [4]
presents a simple treatment of the singular point: the robot does not change
its course inside a ball centered at the singular point. Under some
conditions, the Lyapunov function evaluated at the exit point is proved to be
smaller than that at the entry point. However, these conditions are
conservative, since they assume repulsiveness of the singular points, which is
not the case, e.g., for saddle points.

Related to the existence of singular points, one of the challenges for the VF-
PF navigation problem is to follow a self-intersected desired path. Many
existing VF-PF algorithms (e.g., [3, 16, 25, 17, 13]) fail to fulfil this
task. This is rooted in the fact that the vector field degenerates to zero at
the crossing points of a self-intersected desired path, leading to a zero
guidance signal, and thus a robot can get stuck on the desired path (see Fig.
2). Due to the existence of singular points on the desired path, some
effective VF-PF algorithms such as [16, 22, 13] become invalid simply because
the assumptions are violated in this case. In fact, this task is also
challenging for other existing path-following methods, since in the vicinity
of the crossing points, many methods are “ill-defined”. For example, the line-
of-sight (LOS) method [26] is not applicable as there is not a unique
projection point in the vicinity of the intersection of the desired path.
Indeed, many existing path-following navigation algorithms either focus on
simple desired paths such as circles or straight lines [18, 27, 25], or only
deal with desired paths that are sufficiently smooth [16, 22, 3, 13]. One
might retreat to the virtual-target path-following navigation algorithm [28].
In this method, a virtual target has its own dynamics travelling on the
desired path; thus the path-following navigation problem is implicitly
converted to a target tracking problem. Although through this conversion, it
is possible for a robot to follow a self-intersected desired path, this method
is inherently a tracking approach, and thus may inherit the performance
limitations mentioned before, such as limited path-following accuracy.

Figure 1: The _normalized_ vector field [16] for a circle path described by
$\phi(x,y)=x^{2}+y^{2}-1=0$. The red point at the center of the circle is the
singular point of the (un-normalized) vector field. Figure 2: The _normalized_
vector field [16] for a figure “8” path described by
$\phi(x,y)=x^{2}-4y^{2}(1-y^{2})=0$. Three red points are the singular points
of the (un-normalized) vector field.

Another challenging task arising from the VF-PF methods is the description of
the desired path, which is crucial for the derivation of the vector field. For
generality, the desired path is usually determined by the intersection of
several (hyper)surfaces represented by the zero-level sets of some implicit
functions [3, 16, 22, 4, 29, 13]. For planar desired paths, for example, the
implicit function of a star curve might be as complicated as that in [7],
while for desired paths in a higher-dimensional space, it is counter-intuitive
to create (hyper)surfaces such that the intersection is precisely the desired
path, such as a helix. On the other hand, many geometric curves are described
by parametric equations [30] rather than implicit functions. It is possible to
transform the parametric equations to implicit functions and then derive the
vector field, but this might not always be feasible and is computationally
expensive. The restrictive characterization of the desired path limits the
applicability of VF-PF algorithms to some extent.

In this paper, we improve the VF-PF methodology in the sense that we address
the aforementioned three issues: the existence of singular points, the
obstacle of dealing with self-intersected paths, and the difficulty of
representing a generic desired path. Specifically, based on the design of
guiding vector fields in [22], we use an intuitive idea to eliminate singular
points of the vector field so that global convergence to the desired path,
even if self-intersected, is guaranteed. The general idea is to extend the
dimensions of the vector field and eliminate singular points simultaneously.
This procedure naturally leads to a simple way to transform the descriptions
of desired paths from parameterized forms to the intersection of several
(hyper)surfaces, which are required in creating a guiding vector field.

It is important to clarify the terminology used throughout this paper. In many
VF-PF algorithms, the desired path is a geometric object which is not
necessarily parameterized. In a precise mathematical language, we assume that
the desired path is a _one-dimensional connected differential manifold_ [13].
Therefore, we can generally classify desired paths into two categories: those
homeomorphic to the unit circle $\mathbb{S}^{1}$ if they are compact, and
those homeomorphic to the real line $\mathbb{R}$ otherwise [31, Theorem 5.27].
This assumption is not a restriction, since many desired paths in practice,
such as a circle, an ellipse, a Cassini oval, a spline and a straight line,
satisfy this assumption. For ease of exposition, we refer to those desired
paths homeomorphic to the unit circle as _simple closed_ desired paths, and
those homeomorphic to the real line as _unbounded_ desired paths. Note that
_self-intersected_ paths do not satisfy this assumption. Nevertheless, we will
introduce in the sequel how to _transform_ a self-intersected physical desired
path to a non-self-intersected and unbounded higher-dimensional desired path
such that the assumption holds to apply our algorithm.

Here, we summarize the major contributions of our paper:

Firstly, we show by rigorous topological analysis that guiding vector fields
with the same dimension as the desired path (e.g., [22, 16, 7, 5]) cannot
guarantee the global convergence to a simple closed or self-intersected
desired path (see Theorem 1 in Section III). With the dichotomy of convergence
discussed in the paper, this implies that singular points of the vector field
_always exist_ for a _simple closed_ or _self-intersected_ desired path
regardless of which hypersurfaces one uses to characterize the desired path.
This explains why many vector-field guided path-following navigation
algorithms in the literature cannot guarantee global convergence in the
Euclidean space to a simple closed desired path. We note that excluding
singular points is important in practice (e.g., for fixed-wing aircraft
guidance and navigation) since degenerated or pathological solutions of system
dynamics can be safely avoided. Therefore, this topological obstacle is the
primary motivation of the subsequent theoretical development including the
introduction of extended dynamics (see Section IV) and the creation of
singularity-free guiding vector fields (see Section V).

Secondly, due to the aforementioned topological obstruction, we improve the
existing VF-PF algorithms such that all singular points are removed, and
global convergence of trajectories to the desired path is rigorously
guaranteed (see Section IV and Section V). We overcome this topological
obstruction by changing the topology of the desired paths. Specifically, we
transform a physical _simple closed_ or _self-intersected_ desired path to a
new _unbounded_ and _non-self-intersected_ desired path in a _higher-
dimensional space_. We then derive the corresponding guiding vector field on
this higher-dimensional space, which is guaranteed to have _no_ singular
points.

Thirdly, our proposed method to create this new singularity-free guiding
vector field is proved to enjoy several appealing features (see Subsection
V-B). For example, we provide theoretical guarantees for global exponential
convergence of trajectories of system dynamics to the desired path. In
addition, the new system dynamics with the singularity-free guiding vector
field is robust against perturbation, such as noisy position measurements (see
Feature 3 in Subsection V-B). Moreover, using our proposed method, it becomes
straightforward to represent (hyper)surfaces of which the intersection is the
new higher-dimensional desired path, as long as a parametrization of the
physical (lower-dimensional) desired path is available (see Feature 2 in
Subsection V-B).

Last but not least, we successfully conduct experiments using a fixed-wing
aircraft to verify the effectiveness of our proposed VF-PF algorithm in 3D
(see Section VI), in addition to the experiment with an e-puck robot [32] in
our previous preliminary work [33]. This verifies the practical significance
of our proposed method for highly complex autonomous vehicles. We also discuss
and conclude that our proposed VF-PF algorithm combines and extends features
of the conventional VF-PF algorithms and trajectory tracking algorithms (see
Section VII). While we do not claim that our proposed new singularity-free
guiding vector field is always superior than traditional trajectory tracking
algorithms in every application scenario (such as quadcopter attitude
tracking), we emphasize that it significantly improves conventional VF-PF
algorithms by providing a global solution and enabling the path-following
behavior of complicated or unconventional desired paths (e.g., a self-
intersected Lissajous curve). This is imperative and irreplaceable in
applications such as fixed-wing aircraft guidance and navigation where
convergence to and propagation along a desired path from _every_ initial
position is required.

This paper is a significant extension of our preliminary version in [33] in
the following aspects: 1) [33] only considers physical planar desired paths
and the _higher-dimensional_ 3D vector field, but we generalize all the
results to $n$-dimensional; 2) [33] only investigates a _linear transformation
operator_ , but we introduce a more general _transformation operator_ which
can be linear or nonlinear in this paper. Compared with [33], we present five
new results: 1) we prove a topological theorem regarding the feasibility of
global convergence to desired paths, which serves as the primary motivation
for the subsequent development; 2) we provide rigorous guarantees to justify
the features of our proposed approach; 3) we provide experiments with a fixed-
wing aircraft in 3D to verify the practicality of the proposed approach for
complex engineering systems; 4) we elaborate on some important aspects in the
implementation of the proposed approach; 5) we provide new insightful
discussion of our proposed approach and existing algorithms. Although our
previous work in [13] introduces part of the theoretical foundation for this
paper, it does not solve the singularity problem that will be addressed here.

The remainder of this paper is organized as follows. Section II introduces
conventional guiding vector fields for path following. In Section III, a
theorem about the impossibility of global convergence to simple closed or
self-intersected desired paths using the conventional VF-PF navigation
algorithm is elaborated. This is the main motivation for the design of higher-
dimensional guiding vector fields, which will be utilized in Section IV
through _extended dynamics_. Based on the previous sections, the construction
approach of _singularity-free_ guiding vector fields is presented in Section
V. In addition, several appealing features of this method are highlighted in
this section. Then experiments with a fixed-wing aircraft are conducted to
validate the theoretical results in Section VI. Following this, Section VII
discusses how our proposed approach can be viewed as a combined extension of
VF-PF algorithms and trajectory tracking algorithms. Finally, Section VIII
concludes the paper and indicates future work. Some theoretical proofs are
provided in the Appendix.

We present here some notations and basic concepts that are used throughout the
paper.

Notations: Given a positive integer $n$, the _distance_ between a point
$p_{0}\in\mathbb{R}^{n}$ and a nonempty set
$\mathcal{S}\subseteq\mathbb{R}^{n}$ is denoted by
$\mathrm{dist}(p_{0},\mathcal{S})=\mathrm{dist}(\mathcal{S},p_{0}):=\inf\\{||p-p_{0}||:p\in\mathcal{S}\\}$.
The distance between two nonempty sets $\mathcal{A}$ and $\mathcal{B}$ is
$\mathrm{dist}(\mathcal{A},\mathcal{B})=\mathrm{dist}(\mathcal{B},\mathcal{A})=\inf\\{||a-b||:a\in\mathcal{A},b\in\mathcal{B}\\}$.
If $\xi$ is a differentiable function of time $t$, then the derivative of
$\xi$ with respect to time is $\dot{\xi}(t)$. The notation
$\left\lVert\cdot\right\rVert$ denotes the Euclidean norm of a vector or the
matrix two-norm of a matrix. The normalization of a vector $v$ is denoted by
$\hat{v}$; namely, $\hat{v}=v/\left\lVert v\right\rVert$. The transpose of a
vector $v$ is denoted by $v^{\top}$. Suppose there is a function
$\rho:\mathcal{C}\to\mathcal{D}$ from a nonempty set $\mathcal{C}$ to a
nonempty set $\mathcal{D}$. The image of a subset
$\mathcal{F}\subseteq\mathcal{C}$ under $\rho$ is the subset
$\rho(\mathcal{F})\subseteq\mathcal{D}$ defined by
$\rho(\mathcal{F}):=\\{\rho(x)\in\mathcal{D}:x\in\mathcal{F}\\}$. Given two
functions $f:X\to Y,g:Y\to Z$, $g\circ f$ denotes the composition of these two
functions. The notation $:=$ means “defined to be”.

Basic concepts: A trajectory $\xi:[0,+\infty)\to\mathbb{R}^{n}$ asymptotically
converges to a nonempty set $\mathcal{B}\subseteq\mathbb{R}^{n}$ if for any
$\epsilon>0$, there exists $T>0$ such that
$\mathrm{dist}(\xi(t),\mathcal{B})<\epsilon$ for all $t>T$. If the trajectory
can only be maximally prolonged to $t_{*}<\infty$ [34, Chapter 3], then we say
that it converges to the set $\mathcal{B}$ as $t$ approaches $t_{*}$, if for
any $\epsilon>0$, there exists $\delta>0$ such that
$\mathrm{dist}(\xi(t),\mathcal{B})<\epsilon$ for $|t-t_{*}|<\delta$. Two
subsets $U\subseteq\mathbb{R}^{m}$ and $V\subseteq\mathbb{R}^{k}$ are
_homeomorphic_ , denoted by $U\approx V$, if there exists a continuous
bijection $\varphi:U\to V$ of which the inverse is also continuous. The map
$\varphi$ is called a _homeomorphism_ [31].

##  II Guiding Vector Fields for Path Following

In this section, we introduce the vector-field guided path-following (VF-PF)
navigation problem and guiding vector fields. We formally define the VF-PF
navigation problem as follows.

######  Problem 1 (VF-PF navigation problem).

Given a desired path $\mathcal{P}\subseteq\mathbb{R}^{n}$, the VF-PF
navigation problem is to design a continuously differentiable vector field
$\raisebox{2.0pt}{$\chi$}:\mathbb{R}^{n}\to\mathbb{R}^{n}$ for the
differential equation
$\dot{\xi}(t)=\raisebox{2.0pt}{$\chi$}\big{(}\xi(t)\big{)}$ such that the two
conditions below are satisfied:

  1. 1.

There exists a neighborhood $\mathcal{D}\subseteq\mathbb{R}^{n}$ of the
desired path $\mathcal{P}$ in (1) such that for all initial conditions
$\xi(0)\in\mathcal{D}$, the distance $\mathrm{dist}(\xi(t),\mathcal{P})$
between the trajectory $\xi(t)$ and the desired path $\mathcal{P}$ approaches
zero as time $t\to\infty$; that is,
$\lim_{t\to\infty}\mathrm{dist}(\xi(t),\mathcal{P})=0$;

  2. 2.

If a trajectory starts from the desired path, then the trajectory stays on the
path for $t\geq 0$ (i.e., $\xi(0)\in\mathcal{P}\implies\xi(t)\in\mathcal{P}$
for all $t\geq 0$). In addition, the vector field on the desired path is
nonzero (i.e., $0\notin\raisebox{2.0pt}{$\chi$}(\mathcal{P})$).

In this paper, we only investigate the guiding vector field in the Euclidean
space $\mathbb{R}^{n}$, but it can be extended to a general smooth manifold
[35]. For easy understanding of the guiding vector field on $\mathbb{R}^{n}$
in the sequel, we first introduce the one defined on the two-dimensional
Euclidean space $\mathbb{R}^{2}$.

###  II-A Preliminaries on 2D VF-PF Control

In [16], the desired path $\mathcal{P}$ is described by the zero-level set of
an implicit function:

$\mathcal{P}=\\{(x,y)\in\mathbb{R}^{2}:\phi(x,y)=0\\},$ (1)

where $\phi:\mathbb{R}^{2}\to\mathbb{R}$ is twice continuously differentiable.
In this description, $\mathcal{P}$ is a subset of $\mathbb{R}^{2}$. The
description is different from some other works where the desired path is a
parameterized differentiable curve (e.g., [14]); that is, a differentiable map
$f:I\to\mathbb{R}^{n}$ of an open interval $I=(a,b)$ of the real line
$\mathbb{R}$ into $\mathbb{R}^{n}$ [30]. From the definition, we observe that
the mathematical object in (1) is actually the _trace_ of a parameterized
curve $f$ [30], or the image of a mapping $f$. However, this description of
the desired path without any parametrization is common in the field of VF-PF
control [22, 3, 4, 5, 36, 37, 38, 39, 13]. One of the advantages is that the
vector field can be derived directly from the function $\phi(\cdot)$,
independent of the specific parametrization of the desired path. Another
advantage is that we can replace the calculation of the Euclidean
distance222This calculation is generally difficult since one needs to find the
closest point on the desired path to $\xi$; e.g., it is not trivial for even
an ellipse. $\mathrm{dist}(\xi,\mathcal{P})$ between a point
$\xi\in\mathbb{R}^{2}$ and the desired path $\mathcal{P}$ simply by the value
of $|\phi(\xi)|$. For simplicity, rather than referring to $\mathcal{P}$ in
(1) as “the trace of a parameterized curve”, we call $\mathcal{P}$ the desired
path throughout the paper. In fact, one feature of the VF-PF navigation
problem is that the desired path $\mathcal{P}$ is described by a _one-
dimensional connected submanifold_ , so we have the extra freedom of choosing
different parametrizations or analytic expressions for the same desired path.

If the desired path is non-self-intersected, then a valid 2D vector field
$\raisebox{2.0pt}{$\chi$}:\mathbb{R}^{2}\to\mathbb{R}^{2}$ to solve the VF-PF
navigation problem is [16]:

$\raisebox{2.0pt}{$\chi$}(x,y)=E\nabla\phi(x,y)-k\psi\big{(}\phi(x,y)\big{)}\nabla\phi(x,y),$
(2)

where $E\in SO(2)$ is the $90^{\circ}$ rotation matrix333In fact, the matrix
is $-E$ in [16], but we use $E$ for conventional simplicity. This only affects
the direction of the motion (forward or backward) on the desired path.
$\left[\begin{smallmatrix}0&-1\\\ 1&0\end{smallmatrix}\right]$, and
$\psi:\mathbb{R}\to\mathbb{R}$ is a strictly increasing function satisfying
$\psi(0)=0$. For simplicity, one can choose
$\psi\big{(}\phi(x,y)\big{)}=\phi(x,y)$. The first term of the vector field is
“tangential” to the desired path, thus enables a robot to move along the
desired path, while the second term of the vector field is perpendicular to
the first term, helping the robot move closer to the desired path. Therefore,
intuitively, the vector field guides the robot to move towards and along the
desired path at the same time.

After the preliminaries on the guiding vector field on $\mathbb{R}^{2}$ (i.e.,
2D vector field), we are ready to introduce in the next subsection the more
abstract guiding vector field on $\mathbb{R}^{n}$, where $n\geq 2$.

###  II-B General Guiding Vector Field

Based on the description of a desired path
$\mathcal{P}\subseteq\mathbb{R}^{n}$, the corresponding guiding vector field
for path following is derived in this subsection. Some mild assumptions along
with motivation and justification are presented.

####  II-B1 Guiding Vector Field Structure

We introduce the general method of designing a vector field
$\raisebox{2.0pt}{$\chi$}:\mathbb{R}^{n}\to\mathbb{R}^{n}$ corresponding to a
desired path $\mathcal{P}\subseteq\mathbb{R}^{n}$, where $n\geq 2$. The vector
field is a generalization of those in [16, 22], while its structure is the
same as that in [3]. However, we present assumptions that are crucial for this
method and abandon the unnecessary assumption in [3] about the repulsiveness
of the set of singular points.

Suppose a desired path in the $n$-dimensional Euclidean space is described by
the intersection of $(n-1)$ hypersurfaces; that is,

$\mathcal{P}=\\{\xi\in\mathbb{R}^{n}:\phi_{i}(\xi)=0,\;i=1,\dots,n-1\\},$ (3)

where $\phi_{i}:\mathbb{R}^{n}\to\mathbb{R}$, $i=1,\dots,n-1$, are of
differentiability class $C^{2}$. It is naturally assumed that $\mathcal{P}$ in
(3) is _nonempty_ and _connected_. We further require the regularity of the
desired path as shown later in Assumption 1. For better understanding,
$\phi_{i}(\cdot)=0$ can be regarded as $(n-1)$ constraints, resulting in a one
degree-of-freedom desired path.

######  Remark 1.

Topologically, the desired path $\mathcal{P}$ itself is one-dimensional,
independent of the dimensions of the Euclidean space where it lives. However,
with slight abuse of terminology and for convenience, the desired path
$\mathcal{P}$ is called an _$n$ -D (or $n$D) desired path_ if it lives in the
$n$-dimensional Euclidean space $\mathbb{R}^{n}$ and not in any lower-
dimensional subspace $\mathcal{W}\subseteq\mathbb{R}^{n}$ (i.e., the _smallest
subspace_ the desired path lives in). For example, a planar desired path might
be defined in the three-dimensional Euclidean space $\mathbb{R}^{3}$, but we
only consider the two-dimensional subspace
$\mathcal{W}\subseteq\mathbb{R}^{2}$ where it is contained, and it is thus
natural to call it a $2$D (or $2$-D) desired path rather than a $3$D (or
$3$-D) desired path. Sometimes, for simplicity, we refer to a tangent vector
field on the $n$-dimensional Euclidean space $\mathbb{R}^{n}$ as an _$n$
-dimensional vector field_, and we say that _this vector field is
$n$-dimensional_. ∎

The vector field $\raisebox{2.0pt}{$\chi$}:\mathbb{R}^{n}\to\mathbb{R}^{n}$ is
designed as below:

$\raisebox{2.0pt}{$\chi$}=\times(\nabla\phi_{1},\dots,\nabla\phi_{n-1})-\sum_{i=1}^{n-1}k_{i}\phi_{i}\nabla\phi_{i},$
(4)

where $\nabla\phi_{i}$ is the gradient of $\phi_{i}$, $k_{i}>0$ are constant
gains, and
$\times:\mathbb{R}^{n}\times\dots\times\mathbb{R}^{n}\to\mathbb{R}^{n}$ is the
generalized cross product. In particular, let
$p_{i}=(p_{i1},\cdots,p_{in})^{\top}\in\mathbb{R}^{n}$, $i=1,\dots,n-1$, and
$\bm{b_{j}}\in\mathbb{R}^{n}$ be the standard basis column vector with the
$j$th component being $1$ and the other components being $0$. Then

$\displaystyle\times(p_{1},\dots,p_{n-1})=$
$\displaystyle\begin{vmatrix}p_{12}&\cdots&p_{1n}\\\ \vdots&\vdots&\vdots\\\
p_{n-1,2}&\cdots&p_{n-1,n}\end{vmatrix}\bm{b_{1}}+\dots+$
$\displaystyle(-1)^{n-1}\begin{vmatrix}p_{11}&\cdots&p_{1,n-1}\\\
\vdots&\vdots&\vdots\\\ p_{n-1,1}&\cdots&p_{n-1,n-1}\end{vmatrix}\bm{b_{n}}.$

An intuitive _formal expression_ is

$\times(p_{1},\dots,p_{n-1})=\begin{vmatrix}\bm{b_{1}}&\bm{b_{2}}&\cdots&\bm{b_{n}}\\\
p_{11}&p_{12}&\cdots&p_{1n}\\\ \vdots&\vdots&\vdots&\vdots\\\
p_{n-1,1}&p_{n-1,2}&\cdots&p_{n-1,n}\end{vmatrix}.$

In other words, first forming an $(n-1)$ by $n$ matrix by putting vectors
$p_{1}^{\top},\dots,p_{n-1}^{\top}$ row by row, and then the $k$th component
of $\times(p_{1},\dots,p_{n-1})$ is obtained as $(-1)^{k-1}$ multiplying the
determinant of the submatrix deleting the $k$th column. This is similar to the
wedge product in [3]. The generalized cross product has the following property
[40, Proposition 7.2.1]:

######  Lemma 1 (Orthogonality).

It holds that $\times(p_{1},\dots,p_{n-1})$ is orthogonal to each of the
vectors $p_{1},\dots,p_{n-1}$.

######  Remark 2.

Due to Lemma 1, the physical interpretation of the vector field in (4) is
clear. The first term $\times(\nabla\phi_{1},\dots,\nabla\phi_{n-1})$, being
perpendicular to each $\nabla\phi_{i}$, provides a tangential direction to
each surfaces $\phi_{i}(\xi)=0$, and hence “pushes” the robot along the
desired path. The forward or backward direction of movement along the desired
path is determined by the order of the gradient vectors. Thus if the motion
needs to be reversed, it is sufficient to swap any two of these vectors. We
call this term the _propagation term_. The latter term
$-\sum_{i=1}^{n-1}k_{i}\phi_{i}\nabla\phi_{i}$ provides a direction towards
those surfaces, acting as a “pulling force” to the desired path. We call this
term the _converging term_. ∎

To simplify the notations, we define a matrix
$N(\xi)=\big{(}\nabla\phi_{1}(\xi),\cdots,\nabla\phi_{n-1}(\xi)\big{)}\in\mathbb{R}^{n\times(n-1)}$,
a positive definite gain matrix
$K=\mathrm{diag}\\{k_{1},\dots,k_{n-1}\\}\in\mathbb{R}^{(n-1)\times(n-1)}$ and
a $C^{2}$ function $e:\mathbb{R}^{n}\to\mathbb{R}^{n-1}$ by stacking
$\phi_{i}$; that is,

$e(\xi)=\big{(}\phi_{1}(\xi),\cdots,\phi_{n-1}(\xi)\big{)}^{\top}\in\mathbb{R}^{n-1}.$
(5)

In addition, we define $\nabla_{\times}\phi:\mathbb{R}^{n}\to\mathbb{R}^{n}$
by
$\xi\in\mathbb{R}^{n}\mapsto\times\big{(}\nabla\phi_{1}(\xi),\dots,\nabla\phi_{n-1}(\xi)\big{)}$.
Therefore, the vector field (4) can be compactly written as

$\raisebox{2.0pt}{$\chi$}(\xi)=\nabla_{\times}\phi(\xi)-N(\xi)Ke(\xi).$ (6)

Using this notation, the desired path is equivalent to

$\mathcal{P}=\\{\xi\in\mathbb{R}^{n}:e(\xi)=0\\}.$ (7)

We call $e(\xi)$ the _path-following error_ or simply _error_ between the
point $\xi\in\mathbb{R}^{n}$ and the desired path $\mathcal{P}$. An intuitive
example is a 2D circle described by the zero-level sets of
$e(x,y)=\phi_{1}(x,y)=x^{2}+y^{2}-r^{2}$. As the point $(x,y)$ approaches the
circle, the norm of the path-following error $\left\lVert
e(x,y)\right\rVert=|x^{2}+y^{2}-r^{2}|$ decreases. When $\left\lVert
e(x,y)\right\rVert=0$, the point is right on the path. The use of $\left\lVert
e(\cdot)\right\rVert$ is more convenient than that of the distance function
$\mathrm{dist}(\cdot,\mathcal{P})$. However, there are subtle differences
between the norm of the path-following error $\left\lVert
e(\cdot)\right\rVert$ and the distance $\mathrm{dist}(\cdot,\mathcal{P})$. For
example, when the norm of the path-following error converges to zero, the
trajectory may diverge from the desired path (see Example $3$ in [22]). To
avoid this pathological situation, assumptions are proposed in the next
subsection.

######  Remark 3.

The vector field in (4) is basic while effective in VF-PF navigation problems.
Firstly, we note that many existing studies only deal with simple desired
paths such as a circle or a straight line (or a combination of them) [18, 11,
25], but the vector field in (4) is designed for any _generic_ sufficiently
smooth desired path. Secondly, many vector fields in the literature can be
seen as _variants_ of the vector field in (4). One type of variants are
generated by adding $\phi_{i}$-dependent gains to the converging term or (and)
the propagation term in (4) [5, 16, 7]. Another type of variants adds time-
varying gains or an additional time-varying component [3, 17]. Thus, the
vector fields in [5, 16, 7, 41] can be regarded as 2D variants of (4), and
those in [9, 22, 13, 42, 7] as 3D variants of (4). Therefore, the study of the
basic vector field in (4) is of great significance. Note that we do not
consider time-varying gains or components in the vector field as [3, 17] do.
For one thing, this simplifies the structure of the vector field and
facilitates the practical implementation; for another, this clarifies the
topological property of these vector fields as studied in Section III. For
convenience, we refer to these (time-invariant) vector fields in the
literature as _conventional vector fields_. ∎

####  II-B2 Assumptions

To justify using the norm of the path-following error $\left\lVert
e(\cdot)\right\rVert$ instead of $\mathrm{dist}(\cdot,\mathcal{P})$, we need
some assumptions that are easily satisfied in practice. These assumptions are
based on [22, 16], but are extended to $\mathbb{R}^{n}$. To this end, we
define two sets. The _singular set_ consisting of singular points of a vector
field is defined as below:

$\mathcal{C}=\\{\xi\in\mathbb{R}^{n}:\raisebox{2.0pt}{$\chi$}(\xi)=0\\}.$ (8)

Another related set is

$\mathcal{M}=\left\\{\xi\in\mathbb{R}^{n}:N(\xi)Ke(\xi)=0\right\\}.$ (9)

It can be proved that $\mathcal{M}=\mathcal{P}\cup\mathcal{C}$.

######  Lemma 2.

It holds that $\mathcal{M}=\mathcal{P}\cup\mathcal{C}$.

###### Proof.

See Appendix A-A. ∎

Now we are ready to propose the following assumptions.

######  Assumption 1.

There are no singular points on the desired path. More precisely,
$\mathrm{dist}(\mathcal{C},\mathcal{P})>0$.

######  Assumption 2.

In view of (7), as the norm of the path-following error $\left\lVert
e(\xi)\right\rVert$ approaches zero, the trajectory $\xi(t)$ approaches the
desired path $\mathcal{P}$. Similarly, in view of (9), as the “error”
$\left\lVert N(\xi)Ke(\xi)\right\rVert$ approaches zero, the trajectory
$\xi(t)$ approaches the set $\mathcal{M}$.

Assumption 1 leads to the following statement about the topological property
of the desired path, which is an extension of [22].

######  Lemma 3.

The zero vector $\bm{0}\in\mathbb{R}^{n-1}$ is a regular value of the $C^{2}$
function $e$ in (5), and hence the desired path $\mathcal{P}$ is a $C^{2}$
embedded submanifold444See [24, pp. 105, 98-99] for the definitions of a
regular value and an embedded submanifold. in $\mathbb{R}^{n}$.

###### Proof.

See Appendix A-B ∎

######  Remark 4.

Henceforth the “regularity” of the desired path is guaranteed; namely, the
desired path $\mathcal{P}$ is assumed to be a one-dimensional connected
manifold, which can generally be classified into those homeomorphic to the
unit circle if they are compact, and those homeomorphic the real line
otherwise [31, Theorem 5.27]. Thus throughout the paper, we use the notions of
_simple closed desired paths_ and _desired paths homeomorphic to the unit
circle_ interchangeably. The same applies to _unbounded desired paths_ and
_desired paths homeomorphic to the real line_. Note that self-intersected
desired paths do not satisfy Assumption 1, as shown later in Proposition 1,
but we will propose a method in Section V to transform them into unbounded and
non-self-intersected desired paths, which are then homeomorphic to the real
line $\mathbb{R}$. ∎

Assumption 2 is satisfied for many practical cases such as many polynomial
functions and trigonometric functions; examples are demonstrated in [16, 22,
3, 4]. Since the same desired path can be characterized by various choices of
$\phi_{i}(\cdot)$, these assumptions are crucial to exclude some pathological
cases [22, Example 3]. The precise mathematical formulation of Assumption 2 is
presented in Appendix A-C.

##  III Issues on the Global Convergence to Desired Paths

In this section, we show that, under some conditions, it is not possible to
guarantee global convergence to desired paths using the existing VF-PF
algorithms as introduced in Section II. More specifically, given a desired
path $\mathcal{P}\subseteq\mathbb{R}^{n}$ as described in (3), we investigate
solutions (trajectories) of the autonomous ordinary differential equation:

$\dot{\xi}(t)=\raisebox{2.0pt}{$\chi$}\big{(}\xi(t)\big{)},$ (10)

where $\chi$ is defined in (6). We consider the cases of self-intersected
desired paths and simple closed desired paths respectively as follows.

We first show that the crossing points of a self-intersected desired path
$\mathcal{P}$ are singular points of the corresponding vector field $\chi$ in
(4).

######  Proposition 1.

If the desired path $\mathcal{P}$ in (3) is self-intersected, then the
crossing points of the desired path are singular points of the vector field
$\chi$ in (4).

###### Proof.

See Appendix B-A ∎

######  Remark 5.

This proposition shows that $0\in\raisebox{2.0pt}{$\chi$}(\mathcal{P})$ when
$\mathcal{P}$ is a self-intersected desired path, and therefore, the VF-PF
navigation problem (Problem 1) cannot be addressed as the second requirement
about $0\notin\raisebox{2.0pt}{$\chi$}(\mathcal{P})$ is always violated. In
fact, note that Assumption 1 does not hold in this case, but we will propose
in the sequel an approach to transform a self-intersected desired path such
that Assumption 1 holds and the VF-PF problem is solved. ∎

In Fig. 2, for example, the 2D desired path resembling the figure “$8$” is
self-intersected. It can be numerically calculated that the vector field at
the crossing point is zero. This is intuitive in the sense that there is no
“preference” for the vector at this point to point to either the left or right
portion of the desired path, leaving the only option of zero.

Now we consider simple closed desired paths. In the planar case, due to the
Poincáre-Bendixson theorem [34, Corollary 2.1], there is at least one singular
point of the _2D_ vector field in the region enclosed by the simple closed
desired path. Thus we can conclude that global convergence to a simple closed
_planar_ desired path is not possible. However, this conclusion cannot be
trivially generalized to the higher-dimensional case since the Poincáre-
Bendixson theorem is restricted to the planar case. Nevertheless, we can still
reach this conclusion with some topological analysis.

######  Proposition 2.

If an $n$-D desired path $\mathcal{P}\subseteq\mathbb{R}^{n}$ described by (3)
is simple closed, under the dynamics (10) where the guiding vector field
$\raisebox{2.0pt}{$\chi$}:\mathbb{R}^{n}\to\mathbb{R}^{n}$ is in (4), then it
is not possible to guarantee the global convergence of trajectories of (10) to
the desired path $\mathcal{P}$; precisely, the domain of attraction of
$\mathcal{P}$ cannot be $\mathbb{R}^{n}$.

###### Proof.

See Appendix B-B. ∎

Based on Proposition 1 and Proposition 2, we can reach the following key
statement about the impossibility of global convergence to some desired paths.

######  Theorem 1 (Impossibility of global convergence).

If an $n$-D desired path $\mathcal{P}\subseteq\mathbb{R}^{n}$ described by (3)
is _simple closed_ or _self-intersected_ , then it is not possible to
guarantee the global convergence to the desired path with respect to the
dynamics in (10) with the $n$-dimensional guiding vector field $\chi$ in (4);
more precisely, the domain of attraction of $\mathcal{P}$ cannot be
$\mathbb{R}^{n}$.

###### Proof.

If the desired path $\mathcal{P}$ is self-intersected, then by Proposition 1,
there is at least one singular point on the desired path. Obviously, the path-
following problem formulated by Problem 1 cannot be solved. If the desired
path $\mathcal{P}$ is simple closed, then the global convergence to the
desired path is impossible by Proposition 2. ∎

######  Remark 6.

We note that the topological obstacle to global convergence to the desired
path roots in two aspects: i) the geometry of the desired path: being either
simple closed or self-intersected; ii) the time-invariance property of the
vector field. Although we show this topological obstacle only for the basic
vector field in (4), this obstacle also exists for other variants of vector
fields as listed in Remark 3. This is because the two aspects that lead to the
topological obstruction mentioned above also hold for these time-invariant
vector fields. Note that a state-dependent positive scaling (e.g., the
normalization) of vector fields does not affect the topological properties of
interest (i.e., the phase portrait, or the convergent results) [43,
Proposition 1.14]. ∎

To overcome this topological obstacle and satisfy Assumption 1 even for self-
intersected desired paths, we propose a new idea in the sequel to construct
unbounded and non-self-intersected desired paths from the originally simple
closed or self-intersected desired paths by “stretching” them in a higher-
dimensional space. Indeed, such a higher-dimensional desired path will codify
or contain information about the (lower-dimensional) physical desired path.
Based on the proposed higher-dimensional desired paths, we can derive a
guiding vector field on this higher-dimensional space and show that its
singular set is empty. However, to take advantage of the new guiding vector
field, we need to _transform_ (or _project_ in the linear transformation case)
its integral curves into a lower-dimensional subspace that contains the
information of the physical desired path. The details of transformation into
another space will be discussed in Section IV, and the detailed construction
of a singularity-free guiding vector field on a higher-dimensional space will
be presented in Section V.

##  IV Extended Dynamics and Convergence Results

In this section and the subsequent sections, we consider an $m$-dimensional
Euclidean space $\mathbb{R}^{m}$, where $m>n$. The reason is self-evident as
the paper develops, but it is not necessary to bother with this difference
now. To proceed, we introduce the _extended dynamics_ and derive related
convergence results. The extended dynamics relates to a _transformation
operator_ , which is define as follows:

######  Definition 1 (Transformation operator).

A transformation operator is a function
$G_{l}:\mathbb{R}^{m}\to\mathbb{R}^{m}$ which is twice continuously
differentiable and _globally Lipschitz continuous_ with the Lipschitz constant
$l$.

One can observe that the corresponding Jacobian matrix function of a
transformation operator $\bm{\mathrm{D}}G_{l}=\partial G_{l}/\partial
x:\mathbb{R}^{m}\to\mathbb{R}^{m\times m}$ is locally Lipschitz continuous,
where $x$ is the argument of $G_{l}$. The transformation operator is able to
transform a space into another space (or subspace). One example is a _linear
transformation operator_ defined by $G_{l}(x)=Ax$, where $A$ is a nonzero
matrix, called the _matrix representation_[44, Remark 6.1.15] of this
particular linear transformation operator $G_{l}$. Having defined the
transformation operator, we can now introduce the extended dynamics as
follows.

######  Lemma 4 (Extended dynamics).

Let
$\raisebox{2.0pt}{$\chi$}:\mathcal{D}\subseteq\mathbb{R}^{m}\to\mathbb{R}^{m}$
be a vector field that is locally Lipschitz continuous. Given an initial
condition $\xi(0)=\xi_{0}\in\mathcal{D}$, suppose that $\xi(t)$ is the unique
solution to the differential equation
$\dot{\xi}(t)=\raisebox{2.0pt}{$\chi$}\big{(}\xi(t)\big{)}$, then
$\big{(}\xi(t),{\xi}^{\mathrm{trs}}(t)\big{)}\in\mathbb{R}^{2m}$, where
${\xi}^{\mathrm{trs}}(t):=G_{l}\big{(}\xi(t)\big{)}$ and $G_{l}$ is a
transformation operator, is the unique solution to the following initial value
problem:

$\begin{cases}\dot{\xi}(t)=\raisebox{2.0pt}{$\chi$}\big{(}\xi(t)\big{)}&\xi(0)=\xi_{0}\\\
{\dot{\xi}}^{\mathrm{trs}}(t)=\bm{\mathrm{D}}G_{l}\big{(}\xi(t)\big{)}\cdot\raisebox{2.0pt}{$\chi$}\big{(}\xi(t)\big{)}&{\xi}^{\mathrm{trs}}(0)=G_{l}(\xi_{0}).\end{cases}$
(11)

Moreover, if the trajectory $\xi(t)$ asymptotically converges to some set
$\mathcal{A}\neq\emptyset\subseteq\mathbb{R}^{m}$, then
${\xi}^{\mathrm{trs}}(t)$ asymptotically converges to the _transformed set_
${\mathcal{A}}^{\mathrm{trs}}:=G_{l}(\mathcal{A})=\\{p\in\mathbb{R}^{m}:p=G_{l}(q),\;q\in\mathcal{A}\\}.$

###### Proof.

See Appendix C-A. ∎

We call the ordinary differential equation with the initial condition in (11)
_extended dynamics_. Correspondingly, we call
${\xi}^{\mathrm{trs}}(t):=G_{l}\big{(}\xi(t)\big{)}$ the _transformed
solution_ or _transformed trajectory_ of (11). Before presenting the next
corollary related to the VF-PF navigation problem, we first define the
_transformed desired path_ and _transformed singular set_.

######  Definition 2.

The _transformed desired path_ ${\mathcal{P}}^{\mathrm{trs}}$ of
$\mathcal{P}\subseteq\mathbb{R}^{m}$ in (7) and the _transformed singular set_
${\mathcal{C}}^{\mathrm{trs}}$ of $\mathcal{C}\subseteq\mathbb{R}^{m}$ in (8)
are defined below:

$\displaystyle{\mathcal{P}}^{\mathrm{trs}}$
$\displaystyle:=G_{l}(\mathcal{P})=\\{p\in\mathbb{R}^{m}:p=G_{l}(q),\;q\in\mathcal{P}\\}$
(12) $\displaystyle{\mathcal{C}}^{\mathrm{trs}}$
$\displaystyle:=G_{l}(\mathcal{C})=\\{p\in\mathbb{R}^{m}:p=G_{l}(q),\;q\in\mathcal{C}\\}.$
(13)

In some practical applications, it is desirable to scale the vector field to
have a specified constant length. This is useful if a robot takes the vector
field as the control input directly and is required to move at a constant
speed. In this case, the properties of the integral curves of the scaled
vector field are stated in the corollary below. Recall that the solution
$x(t)$ to an initial value problem $\dot{x}=f(x),\;x(0)=x_{0}$, where $f(x)$
is sufficiently smooth, is not always possible to be prolonged to infinity. In
other words, the solution might only be well-defined in a finite time interval
$[0,t^{*})$, where $t^{*}<\infty$ [43]. The time instant $t^{*}$ is called the
maximal prolonged time of the solution.

######  Corollary 1.

Suppose the desired path $\mathcal{P}$ in (7) is unbounded (i.e.,
$\mathcal{P}\approx\mathbb{R}$). Let
$\raisebox{2.0pt}{$\chi$}:\mathcal{D}\subseteq\mathbb{R}^{m}\to\mathbb{R}^{m}$
be the vector field defined in (4). Suppose $\xi(t)$ is the unique solution to
the initial value problem
$\dot{\xi}(t)=s\hat{\raisebox{2.0pt}{$\chi$}}\big{(}\xi(t)\big{)},\;\xi(0)=\xi_{0}\notin\mathcal{C}$,
where $s>0$ is a constant and $\hat{\cdot}$ is the normalization operator.
Consider the following dynamics

$\begin{cases}\dot{\xi}(t)=s\hat{\raisebox{2.0pt}{$\chi$}}\big{(}\xi(t)\big{)}&\xi(0)=\xi_{0}\notin\mathcal{C}\\\
{\dot{\xi}}^{\mathrm{trs}}(t)=\bm{\mathrm{D}}G_{l}\cdot
s\hat{\raisebox{2.0pt}{$\chi$}}\big{(}\xi(t)\big{)}&{\xi}^{\mathrm{trs}}(0)=G_{l}(\xi_{0}),\end{cases}$
(14)

where $G_{l}$ is a transformation operator. Suppose $t^{*}\leq\infty$ is the
maximal prolonged time of the transformed solution ${\xi}^{\mathrm{trs}}(t)$
to (14). Then ${\xi}^{\mathrm{trs}}(t)$ asymptotically converges to the
transformed desired path ${\mathcal{P}}^{\mathrm{trs}}$ in (12) as
$t\rightarrow\infty$ or the transformed singular set
${\mathcal{C}}^{\mathrm{trs}}$ in (13) as $t\rightarrow t^{*}$.

###### Proof.

See Appendix C-B. ∎

######  Remark 7.

Due to the normalization of the vector field in (14), the right-hand side of
the differential equation is not well defined at singular points of the vector
field. Therefore, if the transformed singular set
${\mathcal{C}}^{\mathrm{trs}}$ is bounded, then the maximal interval to which
the transformed trajectory ${\xi}^{\mathrm{trs}}(t)$ can be prolonged is only
finite when the transformed trajectory ${\xi}^{\mathrm{trs}}(t)$ is converging
to ${\mathcal{C}}^{\mathrm{trs}}$. This happens when the initial value
$\xi(0)$ is in the invariant manifold of the singular set $\mathcal{C}$. ∎

The previous lemma states that the _transformed trajectory_ converges to
either the _transformed desired path_ or the _transformed singular set_ for
initial conditions ${\xi}^{\mathrm{trs}}(0)\in\mathbb{R}^{m}\setminus
G_{l}(\mathcal{C})$, while the latter case is undesirable. A preferable
situation is where the (transformed) singular set is empty. Moreover, as
indicated by Theorem 1, to seek for global convergence, the only possibility
is to consider unbounded and non-self-intersected desired paths (i.e.,
$\mathcal{P}\approx\mathbb{R}$). Therefore, we reach the following corollary.

######  Corollary 2 (Global convergence to ${\mathcal{P}}^{\mathrm{trs}}$).

Suppose the desired path $\mathcal{P}$ in (12) is unbounded (i.e.,
$\mathcal{P}\approx\mathbb{R}$). If $\mathcal{C}=\emptyset$ (equivalently,
${\mathcal{C}}^{\mathrm{trs}}=\emptyset$), then the transformed trajectory
${\xi}^{\mathrm{trs}}(t)$ of (14) globally asymptotically converges to the
transformed desired path ${\mathcal{P}}^{\mathrm{trs}}$ as
$t\rightarrow\infty$ in the sense that the initial condition $\xi(0)$ (and
hence ${\xi}^{\mathrm{trs}}(0)$) can be arbitrarily chosen in
$\mathbb{R}^{m}$.

As will be shown later, only the second differential equation of (11) or (14)
is relevant to the physical robotic system. This corollary thus motivates us
to design a (higher-dimensional) vector field such that the singular set is
empty, in which case global convergence to the (transformed) desired path is
guaranteed. In the next section, we will introduce an intuitive idea to
“stretch” a possibly simple closed or self-intersected physical desired path
and create a higher-dimensional singularity-free vector field.

##  V High-dimensional Singularity-free Guiding Vector Field Construction

Following Section III and Section IV, in this section, we explain how to
implicitly construct an unbounded desired path from the physical desired path
(possibly simple closed or self-intersected) together with a higher-
dimensional guiding vector field without any singular points (a.k.a,
singularity-free guiding vector field).

For simplicity, we restrict the transformation operator
$G_{l}:\mathbb{R}^{m}\to\mathbb{R}^{m}$ to a linear one defined by
$G_{l}(x)=P_{a}x$, where $P_{a}\in\mathbb{R}^{m\times m}$ is a nonzero matrix
defined by

$P_{a}=I-\hat{a}\hat{a}^{\top},$ (15)

where $I$ is the identity matrix of suitable dimensions and
$\hat{a}=a/\left\lVert a\right\rVert\in\mathbb{R}^{m}$ is a normalized nonzero
vector. In this case, $G_{l}$ is actually a linear transformation that
projects an arbitrary vector to the hyperplane orthogonal to the given nonzero
vector $a$, and $P_{a}$ is called the _matrix representation_ of $G_{l}$ [44].
One can observe that the linear transformation $G_{l}$ is globally Lipschitz
continuous with the Lipschitz constant $l=\left\lVert P_{a}\right\rVert=1$,
where $\left\lVert\cdot\right\rVert$ is the induced matrix two-norm. In
addition, the Jacobian is simply $\bm{\mathrm{D}}G_{l}=P_{a}$.

Before formulating the problem in the sequel, we define the _coordinate
projection function_ $\pi_{(1,\dots,n)}:\mathbb{R}^{m}\to\mathbb{R}^{n}$ as

$\pi_{(1,\dots,n)}(x_{1},\dots,x_{n},\dots,x_{m})=(x_{1},\dots,x_{n}),$

where $m>n$. In other words, the coordinate projection function
$\pi_{(1,\dots,n)}$ takes only the first $n$ components of an $m$-dimensional
vector and generates a lower-dimensional one.

######  Problem 2.

Given an $n$-D physical desired path555Recall the notion of an $n$-D desired
path in Remark 1. ${\mathcal{P}}^{\mathrm{phy}}$ in $\mathbb{R}^{n}$, we aim
to find an $m$-D desired path ${\mathcal{P}}^{\mathrm{hgh}}$ in
$\mathbb{R}^{m}$, where $m>n$, which satisfies the following conditions:

  1. 1.

there exist functions $\phi_{i}(\cdot)$, $i=1,\dots,m-1$, such that
${\mathcal{P}}^{\mathrm{hgh}}$ is described by (3);

  2. 2.

the singular set ${\mathcal{C}}^{\mathrm{hgh}}$ of the _higher-dimensional_
vector field
${\raisebox{2.0pt}{$\chi$}}^{\mathrm{hgh}}:\mathbb{R}^{m}\to\mathbb{R}^{m}$ in
(4) corresponding to ${\mathcal{P}}^{\mathrm{hgh}}$ is empty;

  3. 3.

there exists a transformation operator $G_{l}:\mathbb{R}^{m}\to\mathbb{R}^{m}$
such that
$\pi_{(1,\dots,n)}({\mathcal{P}}^{\mathrm{trs}})={\mathcal{P}}^{\mathrm{phy}}$,
where the transformed desired path
${\mathcal{P}}^{\mathrm{trs}}=G_{l}({\mathcal{P}}^{\mathrm{hgh}})$.

######  Remark 8.

It is important to distinguish among the _physical desired path_
${\mathcal{P}}^{\mathrm{phy}}$, the _higher-dimensional desired path_
${\mathcal{P}}^{\mathrm{hgh}}$ and the _transformed desired path_
${\mathcal{P}}^{\mathrm{trs}}$. A major difference is the dimensions of their
ambient space; that is, ${\mathcal{P}}^{\mathrm{phy}}\subseteq\mathbb{R}^{n}$,
while
${\mathcal{P}}^{\mathrm{hgh}},{\mathcal{P}}^{\mathrm{trs}}\subseteq\mathbb{R}^{m}$
and $m>n$. Although the higher-dimensional desired path
${\mathcal{P}}^{\mathrm{hgh}}$ and the transformed desired path
${\mathcal{P}}^{\mathrm{trs}}$ both live in $\mathbb{R}^{m}$, the transformed
desired path ${\mathcal{P}}^{\mathrm{trs}}$ lives in a subspace
$\mathcal{W}\subseteq\mathbb{R}^{m}$ _probably_ with $\dim(\mathcal{W})<m$
since ${\mathcal{P}}^{\mathrm{trs}}=G_{l}({\mathcal{P}}^{\mathrm{hgh}})$.
Indeed, for the case of a linear transformation operator in (15), the
transformed desired path
${\mathcal{P}}^{\mathrm{trs}}=P_{a}({\mathcal{P}}^{\mathrm{hgh}})$ lives in
the orthogonal complement subspace $\mathcal{W}$ of the linear space spanned
by the vector $a$ (i.e., $\mathrm{span}\\{a\\}$), and
$\dim(\mathcal{W})=m-1<m$. ∎

Next, we propose the solution to Problem 2 in Subsection V-A. Having found the
higher-dimensional desired path ${\mathcal{P}}^{\mathrm{hgh}}$, then we can
directly derive the corresponding vector field
${\raisebox{2.0pt}{$\chi$}}^{\mathrm{hgh}}$ defined on $\mathbb{R}^{m}$ by
(4). Some features of the approach illustrated in Subsection V-A are
highlighted in Subsection V-B.

###  V-A Construction of a Singularity-free Guiding Vector Field

Suppose an $n$-D physical desired path ${\mathcal{P}}^{\mathrm{phy}}$ is
parameterized by

$x_{1}=f_{1}(w),\,\dots,\,x_{n}=f_{n}(w),$ (16)

where $w\in\mathbb{R}$ is the parameter of the desired path and $f_{i}\in
C^{2},\,i=1,\dots,n$. We can simply let

$\phi_{1}(\xi)=x_{1}-f_{1}(w),\;\dots,\;\phi_{n}(\xi)=x_{n}-f_{n}(w),$ (17)

where $\xi=(x_{1},\dots,x_{n},w)$ has an additional coordinate $w$ now and is
an $m$-dimensional vector, where $m=n+1$. So the $m$-D desired path is

${\mathcal{P}}^{\mathrm{hgh}}=\\{\xi=(x_{1},\dots,x_{n},w)\in\mathbb{R}^{m}:\phi_{i}(\xi)=0,\,i=1,\dots,n\\}.$

(18)

Thus the first requirement of Problem 2 is met. Intuitively, the new higher-
dimensional desired path ${\mathcal{P}}^{\mathrm{hgh}}$ is obtained by
“stretching” the $n$-D desired path ${\mathcal{P}}^{\mathrm{phy}}$ along the
additional virtual $w$-axis. From the higher-dimensional desired path
${\mathcal{P}}^{\mathrm{hgh}}\subseteq\mathbb{R}^{m}$ in (18), we obtain the
corresponding guiding vector field on the higher-dimensional space
$\mathbb{R}^{m}$ by (4):

${\raisebox{2.0pt}{$\chi$}}^{\mathrm{hgh}}=\nabla_{\times}\phi-\sum_{i=1}^{n}k_{i}\phi_{i}\nabla\phi_{i}.$

It can be calculated that
$\nabla\phi_{i}=\big{(}0,\dots,1,\dots,-f^{\prime}_{i}(w)\big{)}^{\top}$ for
$i=1,\dots,n$, where
$f^{\prime}_{i}(w):=\frac{\mathrm{d}f_{i}(w)}{\mathrm{d}w}$ and $1$ is the
$i$-th component of the gradient vector. Therefore,

$\nabla_{\times}\phi=\times(\nabla\phi_{1},\dots,\nabla\phi_{n})=(-1)^{n}\begin{bmatrix}f^{\prime}_{1}(w)\\\
\vdots\\\ f^{\prime}_{n}(w)\\\
1\end{bmatrix}\in\mathbb{R}^{m}=\mathbb{R}^{n+1}.$

It is interesting that the $m$-th coordinate of this vector is a constant
$(-1)^{n}$ regardless of the specific parametric form of the desired path.
This means that $\left\lVert\nabla_{\times}\phi(\xi)\right\rVert\neq 0$ for
$\xi\in\mathbb{R}^{m}$ globally. From Lemma 1, we know that the propagation
term $\nabla_{\times}\phi$ of the vector field is always linearly independent
from the converging term $\sum_{i=1}^{n}k_{i}\phi_{i}\nabla\phi_{i}$ unless
they are zero vectors. However, as we have shown that
$\left\lVert\nabla_{\times}\phi\right\rVert\neq 0$ in $\mathbb{R}^{m}$
globally, this reveals the appealing property that the vector field
${\raisebox{2.0pt}{$\chi$}}^{\mathrm{hgh}}(\xi)\neq 0$ for any point
$\xi\in\mathbb{R}^{m}$, implying that there are no singular points in the
higher-dimensional space $\mathbb{R}^{m}$; i.e.,
${\mathcal{C}}^{\mathrm{hgh}}=\emptyset$. Thus, the second requirement of
Problem 2 (as well as a related condition in Corollary 2) is satisfied.

To let the third requirement of Problem 2 be satisfied, we retreat to a linear
transformation operator with a matrix representation $P_{a}$. One of the
simplest linear transformation operators corresponds to
$a=\bm{b_{n+1}}\in\mathbb{R}^{m}$, which is a standard basis column vector
with the $(n+1)$-th component being $1$ and the other components being $0$.
This can be used to transform an $m$-dimensional space to an $n$-dimensional
subspace by “zeroing” the last coordinate. Specifically, we let
$a=\bm{b_{n+1}}$, then the matrix representation of the linear transformation
operator is $P_{a}=\left[\begin{smallmatrix}I_{n\times n}&\bm{0}\\\
\bm{0}&0\end{smallmatrix}\right]$, where $\bm{0}$ are zero vectors that are of
suitable dimensions. Observe that the $n$-D desired path
${\mathcal{P}}^{\mathrm{phy}}\subseteq\mathbb{R}^{n}$ is the orthogonal
projection of the higher-dimensional desired path
${\mathcal{P}}^{\mathrm{hgh}}\subseteq\mathbb{R}^{m}$ on the plane where
$w=0$; that is,

$\pi_{(1,\dots,n)}({\mathcal{P}}^{\mathrm{hgh}})=\pi_{(1,\dots,n)}({\mathcal{P}}^{\mathrm{trs}})={\mathcal{P}}^{\mathrm{phy}}.$

Therefore, the third requirement of Problem 2 is also satisfied. By the
construction in (18), the higher-dimensional desired path
${\mathcal{P}}^{\mathrm{hgh}}\subseteq\mathbb{R}^{m}$ satisfying all the
conditions in Problem 2 is thus found. Ultimately, we can take advantage of
the new “well-behaved” guiding vector field
${\raisebox{2.0pt}{$\chi$}}^{\mathrm{hgh}}\in\mathbb{R}^{m}$ derived from
${\mathcal{P}}^{\mathrm{hgh}}\subseteq\mathbb{R}^{m}$ as mentioned above. This
result is formally stated in the following theorem.

######  Theorem 2.

Suppose an $n$-D physical desired path
${\mathcal{P}}^{\mathrm{phy}}\subseteq\mathbb{R}^{n}$ is parameterized by
(16). If $\phi_{1},\dots,\phi_{n}$ are chosen as in (17), then there are no
singular points in the corresponding guiding vector field
${\raisebox{2.0pt}{$\chi$}}^{\mathrm{hgh}}:\mathbb{R}^{n+1}\to\mathbb{R}^{n+1}$
defined on the $(n+1)$-dimensional space $\mathbb{R}^{n+1}$. Let
$a=\bm{b_{n+1}}$ for the linear transformation operator $P_{a}$. Suppose the
transformed trajectory of the _extended dynamics_ (14) is
${\xi}^{\mathrm{trs}}(t):=\big{(}x_{1}(t),\dots,x_{n}(t),w(t)\big{)}^{\top}$.
Then the _projected transformed trajectory_

${\xi}^{\mathrm{prj}}(t):=\pi_{(1,\dots,n)}\big{(}{\xi}^{\mathrm{trs}}(t)\big{)}=\big{(}x_{1}(t),\dots,x_{n}(t)\big{)}^{\top}$

globally asymptotically converges to the physical desired path
${\mathcal{P}}^{\mathrm{phy}}$ as $t\to\infty$.

###### Proof.

By (4) and (17), the guiding vector field on the $(n+1)$-dimensional space
$\mathbb{R}^{n+1}$ is

${\raisebox{2.0pt}{$\chi$}}^{\mathrm{hgh}}(x_{1},\dots,x_{n},w)=\begin{bmatrix}(-1)^{n}f_{1}^{\prime}(w)-k_{1}\phi_{1}\\\
\vdots\\\ (-1)^{n}f_{n}^{\prime}(w)-k_{n}\phi_{n}\\\
(-1)^{n}+\sum_{i=1}^{n}k_{i}\phi_{i}f_{i}^{\prime}(w)\end{bmatrix}.$ (19)

As discussed before, the singular set
${\mathcal{C}}^{\mathrm{hgh}}=\emptyset$. According to Corollary 2,
${\xi}^{\mathrm{trs}}(t)$ globally asymptotically converges to the transformed
desired path
${\mathcal{P}}^{\mathrm{trs}}=G_{l}({\mathcal{P}}^{\mathrm{hgh}})=P_{a}({\mathcal{P}}^{\mathrm{hgh}})$
as $t\rightarrow\infty$. Since
$a^{\top}{\xi}^{\mathrm{trs}}=a^{\top}P_{a}\xi=0$, the $(n+1)$-th coordinate
$w(t)$ of the transformed trajectory ${\xi}^{\mathrm{trs}}(t)$ is equal to
$0$, meaning that the transformed trajectory ${\xi}^{\mathrm{trs}}(t)$ lies on
the subspace
$\mathcal{W}:=\\{(x_{1},\dots,x_{n+1})\in\mathbb{R}^{n}:x_{n+1}=0\\}$.
Therefore, the projected transformed trajectory
${\xi}^{\mathrm{prj}}(t)=\pi_{(1,\dots,n)}\big{(}{\xi}^{\mathrm{trs}}(t)\big{)}$
globally asymptotically converges to the physical desired path
${\mathcal{P}}^{\mathrm{phy}}$. ∎

######  Remark 9.

Note that the proof of convergence to the physical desired path
${\mathcal{P}}^{\mathrm{phy}}$ is indirect. The norm of the path-following
error $\left\lVert
e(\cdot)\right\rVert=\left\lVert\big{(}\phi_{1}(\cdot),\dots,\phi_{n}(\cdot)\big{)}\right\rVert$
captures the distance to the higher-dimensional desired path
${\mathcal{P}}^{\mathrm{hgh}}$, taking into account the additional coordinate
$w$ as well. It is shown first that in the higher-dimensional space
$\mathbb{R}^{n+1}$, the norm of the path-following error $\left\lVert
e(\cdot)\right\rVert$ approaches zero asymptotically. Then the convergence to
the transformed desired path ${\mathcal{P}}^{\mathrm{trs}}$ is obtained from
Corollary 1 (or Corollary 2). Due to the special choice of the linear
transformation operator $P_{a}$, where $a=\bm{b_{n+1}}$, the transformed
desired path ${\mathcal{P}}^{\mathrm{trs}}$ is “almost” the same as the
physical desired path ${\mathcal{P}}^{\mathrm{phy}}$, except that it has an
additional but constant coordinate $w(t)\equiv 0$. ∎

We have shown that, by extending the vector field from $\mathbb{R}^{n}$ to
$\mathbb{R}^{n+1}$, the new guiding vector field does not have any singular
points. Therefore, by using the extended dynamics, the convergence to the
physical desired path is guaranteed globally. When $n>3$, this case
corresponds to some configuration spaces, such as the robot arm joint space in
a smooth manifold embedded in $\mathbb{R}^{n}$. See [35] for more details.

###  V-B Features of the Approach

There are several intriguing features of our proposed approach discussed above
in Section V-A. These features are summarized below, and the corresponding
theoretical guarantees are presented in Appendix D. For ease of narration and
without loss of generality, we take the case of a 2D physical desired path
${\mathcal{P}}^{\mathrm{phy}}\subseteq\mathbb{R}^{2}$ for discussion (i.e.,
$n=2$).

######  Feature 1.

The corresponding higher-dimensional desired path
${\mathcal{P}}^{\mathrm{hgh}}=\\{\xi\in\mathbb{R}^{2+1}:\phi_{1}(\xi)=0,\phi_{2}(\xi)=0\\}$
is not self-intersected. This is due to the fact that a crossing point must be
a singular point (see Proposition 1), but we have shown that there are no
singular points in the higher-dimensional guiding vector field. In fact, the
parameter of the desired path $w$ in (16) is implicitly transformed to an
additional coordinate of the higher-dimensional desired path. Thus the
physical planar desired path ${\mathcal{P}}^{\mathrm{phy}}$ is “cut” and
“stretched” into the three-dimensional Euclidean space, and becomes unbounded
and non-self-intersected along the additional dimension. The significance of
this feature is that even a self-intersected physical desired path
${\mathcal{P}}^{\mathrm{phy}}$ described by (16) can be successfully followed
by using the new singularity-free guiding vector field, which in fact
corresponds to a non-self-intersected “stretched” desired path
${\mathcal{P}}^{\mathrm{hgh}}$.

######  Feature 2.

This approach facilitates the expression of hypersurfaces characterized by
implicit functions $\phi_{i}$. Usually, a parameterized form of the desired
path is more readily available than the hypersurfaces of which the
intersection is the desired path. Therefore, given the parameterized form in
(16), we do not need to convert them into $\phi(x,y)=0$ and derive the
corresponding 2D vector field. Instead, by simply defining two $\phi$
functions as in (17), we obtain a singularity-free vector field
${\raisebox{2.0pt}{$\chi$}}^{\mathrm{hgh}}$ defined on $\mathbb{R}^{2+1}$.

######  Feature 3.

One only needs to examine the boundedness of $|f^{\prime}_{i}(z)|$, $i=1,2$,
in the vicinity of the higher-dimensional desired path
${\mathcal{P}}^{\mathrm{hgh}}$ to guarantee both the property of local
exponential vanishing of the norm of the path-following error $\left\lVert
e\right\rVert$ and the property of robustness against disturbance of the
system dynamics (10), while these properties usually require more conditions
to be satisfied for general vector fields [13]. The theoretical guarantees are
shown in Proposition 4 and Proposition 5 in Appendix D respectively.

######  Feature 4.

Only Assumption 2 is required. Since the new guiding vector field does not
have any singular points, the other assumption, Assumption 1, is vacuously
true. This is independent of the specific parametrizations of the desired path
in (16). The intuitive Assumption 2 holds for many practical examples; thus
one might safely ignore it in practice [13, 22, 16].

######  Feature 5.

The additional virtual coordinate can be used to realize scalable distributed
multi-robot coordinated path-following navigation by adding a consensus term
and guarantee collision avoidance by using a safety barrier certificate [45].

##  VI Experiments with an autonomous aircraft

In this section, we will demonstrate the effectiveness of our path tracking
approach with an autonomous fixed-wing aircraft. In particular, we are going
to verify the tracking of both 2D and 3D self-intersected desired paths. All
the related software has been developed within the open-source project for
autopilots _Paparazzi_ [46, 23, 6]. In fact, the codes are ready for others to
implement their trajectories and only the corresponding parametric equations
are
needed666https://github.com/noether/paparazzi/tree/gvf_advanced/sw/airborne/modules/guidance/gvf_parametric..

###  VI-A The autonomous aircraft and airfield

For the experiments, we use one off-the-shelf 1.2m delta-wing Opterra as shown
in Fig. 3. Two _elevons_ actuate the aircraft at the wings and one motor acts
in pushing the configuration. The vehicle’s electronics consists of the
autopilot _Apogee_ 777https://wiki.paparazziuav.org/wiki/Apogee/v1.00., an
Ublox GPS receptor, a Futaba receiver, and a X-Bee S1 radio modem. The
Apogee’s core is an STM32F4 microcontroller where our algorithm runs with a
fixed frequency of 50Hz, and all the relevant data are logged in an SD card at
100Hz. The ground segment consists of a standard laptop with another X-Bee S1
radio modem to monitor the telemetry and a Futaba transmitter in case of
taking over manual control of the vehicle.

Figure 3: Autonomous Opterra 1.2m equipped with Paparazzi’s Apogee autopilot.
The airframe is built by E-Flite / Horizon Hobby company.

The flights took place on July 18, 2020, in Ciudad Real (Spain) with GPS
coordinates $(39.184535,-4.020797)$ degrees. The weather forecast reported 36
Celsius degrees and a South wind of 14 km/h.

###  VI-B Aircraft’s guidance system design

We employ a decoupled vertical and horizontal model for setting the aircraft’s
guiding reference signals. In particular, accounting for the non-holonomic
lateral constraint of the aircraft, we consider the following unicycle model

$\dot{x}=v\cos\theta\quad\dot{y}=v\sin\theta\quad\dot{\theta}=u_{\theta}\quad\dot{z}=u_{z},$
(20)

where $(x,y,z)$ is the 3D position, $\theta$ is the heading angle on the XY
plane, $v$ is the ground speed, $u_{\theta}$ is the angular velocity
control/guiding signal to change the heading, and $u_{z}$ is the guiding
signal for the climbing velocity. In this section, we are going to show how to
design the guiding signals $u_{\theta}$ and $u_{z}$, which are injected into
the control system of the aircraft that deals with the non-trivial couplings
of the lateral and vertical modes. In particular, $u_{\theta}$ is tracked by
banking the aircraft depending on the current speed $v$ and _pitch_ angle to
achieve a _coordinated turn_
888https://en.wikipedia.org/wiki/Coordinated_flight., and $u_{z}$ is tracked
by controlling the pitch angle and the propulsion to vary the lift and the
vertical component of the pushed force coming from the propeller999We leave
the reader to check the details of the employed low-level controllers at
http://wiki.paparazziuav.org/wiki/Control_Loops.. The experimental results
will show the compatibility of our path tracking algorithm together with the
model (20) and the low-level controller employed in Paparazzi for a fixed-wing
aircraft.

We note that the wind has a noticeable impact on the ground speed of the
aircraft. Nevertheless, as the experimental results indicate, such a wind
speed does not impact the intended performance of the algorithm. In practice,
we consider $\theta$ as the heading angle (given by the velocity vector), not
the attitude _yaw_ angle. If there is no wind, both angles are the same in our
setup. When we consider the heading instead of the yaw for the unicycle model
(20), the aircraft ends up compensating the lateral wind by _crabbing_ so that
aerodynamic angle _sideslip_ is almost zero.101010Crabbing happens when the
inertial velocity makes an angle with the nose heading due to wind. _Slipping_
happens when the aerodynamic velocity vector makes an angle (sideslip) with
the body ZX plane. Slipping is (almost) always undesirable because it degrades
aerodynamic performance. Crabbing is not an issue for the aircraft.

For tracking 3D paths (including 2D paths at constant altitude), we will
employ a higher-dimensional 4D vector field. The generalized 4D velocity
vector of the aircraft is defined as
$\dot{\xi}=(\dot{x},\dot{y},\dot{z},\dot{w})^{\top}$, where
$(\dot{x},\dot{y})$ is the actual ground velocity of the aircraft, $\dot{z}$
is the vertical speed, and $\dot{w}$ is the velocity in the additional
coordinate to be determined.

Now we present the control algorithm design; that is, the design of
$u_{\theta}$ and $u_{z}$ in (20) with the following proposition:

######  Proposition 3.

Suppose the 3D physical desired path
${\mathcal{P}}^{\mathrm{phy}}\subseteq\mathbb{R}^{3}$ to follow is
parameterized by (16). Then a corresponding 4D vector field
$\raisebox{2.0pt}{$\chi$}:\mathbb{R}^{4}\to\mathbb{R}^{4}$ can be constructed
by Theorem 2. Assume that the vector field satisfies
$\raisebox{2.0pt}{$\chi$}_{1}(\xi)^{2}+\raisebox{2.0pt}{$\chi$}_{2}(\xi)^{2}\neq
0$ for $\xi\in\mathbb{R}^{4}$, where $\raisebox{2.0pt}{$\chi$}_{i}$ denotes
the $i$-th entry of $\chi$. Consider the model (20), and let the dynamics of
the additional coordinate $w$ be

$\dot{w}=\frac{v\raisebox{2.0pt}{$\chi$}_{4}}{\sqrt{\raisebox{2.0pt}{$\chi$}_{1}^{2}+\raisebox{2.0pt}{$\chi$}_{2}^{2}}}.$
(21)

Let the angular velocity control input $u_{\theta}$ and the climbing velocity
input $u_{z}$ be

$\displaystyle u_{\theta}$
$\displaystyle=\underbrace{\left(\frac{-1}{\left\lVert\raisebox{2.0pt}{$\chi$}^{p}\right\rVert}\hat{\raisebox{2.0pt}{$\chi$}^{p}}^{\top}EJ(\raisebox{2.0pt}{$\chi$}^{p})\dot{\xi}\right)}_{:=\dot{\theta}_{d}}-k_{\theta}\hat{h}^{\top}E\hat{\raisebox{2.0pt}{$\chi$}^{p}},$
(22a) $\displaystyle u_{z}$
$\displaystyle=\frac{v\raisebox{2.0pt}{$\chi$}_{3}}{\sqrt{\raisebox{2.0pt}{$\chi$}_{1}^{2}+\raisebox{2.0pt}{$\chi$}_{2}^{2}}},$
(22b)

where $k_{\theta}>0$ is a gain constant, $h=(\cos\theta,\sin\theta)^{\top}$,
$E=\left[\begin{smallmatrix}0&-1\\\ 1&0\end{smallmatrix}\right]$,
$\raisebox{2.0pt}{$\chi$}^{p}=(\hat{\raisebox{2.0pt}{$\chi$}}_{1},\hat{\raisebox{2.0pt}{$\chi$}}_{2})^{\top}$
and $J(\raisebox{2.0pt}{$\chi$}^{p})$ is the Jacobian matrix of
$\raisebox{2.0pt}{$\chi$}^{p}$ with respect to the generalized position
$\xi=(x,y,z,w)$ and $\dot{\xi}=(\dot{x},\dot{y},\dot{z},\dot{w})^{\top}$ is
the generalized velocity. Let the angle difference directed from
$\hat{\raisebox{2.0pt}{$\chi$}^{p}}$ to $\hat{h}$ be denoted by
$\beta\in(-\pi,\pi]$. If the initial angle difference satisfies
$\beta(0)\in(-\pi,\pi)$, then it will vanish asymptotically (i.e.,
$\beta(t)\to 0$). Furthermore, the actual robot trajectory $(x(t),y(t),z(t))$
will converge to the physical desired path ${\mathcal{P}}^{\mathrm{phy}}$
asymptotically as $t\to\infty$.

###### Proof.

Let

$\raisebox{2.0pt}{$\chi$}^{\prime}:=\frac{1}{\sqrt{\raisebox{2.0pt}{$\chi$}_{1}^{2}+\raisebox{2.0pt}{$\chi$}_{2}^{2}}}\raisebox{2.0pt}{$\chi$}$

be the scaled 4D vector field. We aim to let the generalized robot velocity
$\dot{\xi}=(\dot{x},\dot{y},\dot{z},\dot{w})^{\top}$ eventually align with and
point towards the same direction as the scaled vector field. Specifically, let
the orientation error be defined by

$e_{ori}(t)=\dot{\xi}-v\raisebox{2.0pt}{$\chi$}^{\prime}\stackrel{{\scriptstyle\eqref{eq_u_w},\eqref{eq_u_z}}}{{=}}v\begin{bmatrix}\cos\theta-\raisebox{2.0pt}{$\chi$}^{\prime}_{1}\\\
\sin\theta-\raisebox{2.0pt}{$\chi$}^{\prime}_{2}\\\ 0\\\
0\end{bmatrix}=\begin{bmatrix}h-g\\\ \bm{0}\end{bmatrix}\in\mathbb{R}^{4},$

where $h=(\cos\theta,\sin\theta)^{\top}$ and
$g=(\raisebox{2.0pt}{$\chi$}^{\prime}_{1},\raisebox{2.0pt}{$\chi$}^{\prime}_{2})^{\top}$.
It is obvious that $e_{ori}\to 0$ if and only if $h-g\to 0$. Therefore, it
suffices to show that the orientation of $h$ asymptotically aligns with that
of $g$. Note that

$\hat{\raisebox{2.0pt}{$\chi$}^{p}}=\frac{1}{\sqrt{\hat{\raisebox{2.0pt}{$\chi$}}_{1}^{2}+\hat{\raisebox{2.0pt}{$\chi$}}_{2}^{2}}}\begin{bmatrix}\hat{\raisebox{2.0pt}{$\chi$}}_{1}\\\
\hat{\raisebox{2.0pt}{$\chi$}}_{2}\end{bmatrix}=\frac{1}{\sqrt{\raisebox{2.0pt}{$\chi$}_{1}^{2}+\raisebox{2.0pt}{$\chi$}_{2}^{2}}}\begin{bmatrix}\raisebox{2.0pt}{$\chi$}_{1}\\\
\raisebox{2.0pt}{$\chi$}_{2}\end{bmatrix}=g$

and $\hat{h}=h$. Therefore, we can define a new orientation error as
$e_{ori}^{\prime}:=\hat{h}-\hat{\raisebox{2.0pt}{$\chi$}^{p}}\in\mathbb{R}^{2}$.
Choose the Lyapunov function candidate
$V=1/2\,e^{\prime\top}_{ori}e^{\prime}_{ori}$ and its time derivative is

$\begin{split}\dot{V}=\dot{e}^{\prime\top}_{ori}e^{\prime}_{ori}&=(\dot{\theta}E\hat{h}-\dot{\theta}_{d}E\hat{\raisebox{2.0pt}{$\chi$}^{p}})^{\top}(\hat{h}-\hat{\raisebox{2.0pt}{$\chi$}^{p}})\\\
&=(\dot{\theta}-\dot{\theta}_{d})\hat{h}^{\top}E\hat{\raisebox{2.0pt}{$\chi$}^{p}}\\\
&\overset{\eqref{eq_u_theta}}{=}-k_{\theta}(\hat{h}^{\top}E\hat{\raisebox{2.0pt}{$\chi$}^{p}})^{2}\leq
0.\end{split}$ (23)

The second equation makes use of the identities:
$\dot{\hat{h}}=\dot{\theta}E\hat{h}$ and
$\dot{\hat{\raisebox{2.0pt}{$\chi$}^{p}}}=\dot{\theta}_{d}E\hat{\raisebox{2.0pt}{$\chi$}^{p}}$,
where $\dot{\theta}_{d}$ is defined in (22a). The third equation is derived by
exploiting the facts that $E^{\top}=-E$ and $a^{\top}Ea=0$ for any vector
$a\in\mathbb{R}^{2}$. Note that $\dot{V}=0$ if and only if the angle
difference between $\hat{h}$ and $\hat{\raisebox{2.0pt}{$\chi$}^{p}}$ is
$\beta=0$ or $\beta=\pi$. Since it is assumed that the initial angle
difference $\beta(t=0)\neq\pi$, it follows that $\dot{V}(t=0)<0$, and thus
there exists a sufficiently small $\epsilon>0$ such that
$V(t=\epsilon)<V(t=0)$. It can be shown by contradiction that $|\beta(t)|$ is
monotonically decreasing with respect to time $t$111111Suppose there exist
$0<t_{1}<t_{2}$ such that $|\beta(t_{1})|<|\beta(t_{2})|$. It can be
calculated that $V(t)=1-\cos\beta(t)$, and thus $V(t_{1})<V(t_{2})$,
contradicting the decreasing property of $\dot{V}$. Thus $|\beta(t)|$ is
indeed monotonically decreasing.. By (23), one observes that $|\beta(t)|$ and
$V(t)$ tends to $0$, implying that the generalized velocity $\dot{\xi}$ will
converge asymptotically to the scaled vector field
$v\raisebox{2.0pt}{$\chi$}^{\prime}$. Note that the integral curves of the
state-dependent positive scaled vector field
$\raisebox{2.0pt}{$\chi$}^{\prime}$ has the same convergence results as those
for the original vector field $\chi$ [43, Proposition 1.14]. Therefore, the
generalized trajectory $(x(t),y(t),z(t),w(t))$ will converge to the higher-
dimensional desired path ${\mathcal{P}}^{\mathrm{hgh}}$ in (18). From Theorem
2, the actual robot trajectory (i.e., the projected transformed trajectory)
$(x(t),y(t),z(t))$ will converge to the physical desired path
${\mathcal{P}}^{\mathrm{phy}}$ asymptotically as $t\to\infty$. ∎

$J(\raisebox{2.0pt}{$\chi$})=L\begin{bmatrix}-k_{1}L&0&0&-f_{1}^{\prime\prime}(\beta
w)L^{2}\beta^{2}+k_{1}\beta Lf_{1}^{\prime}(\beta w)\\\
0&-k_{2}L&0&-f_{2}^{\prime\prime}(\beta w)L^{2}\beta^{2}+k_{2}\beta
Lf_{2}^{\prime}(\beta w)\\\ 0&0&-k_{3}L&-f_{3}^{\prime\prime}(\beta
w)L^{2}\beta^{2}+k_{3}\beta Lf_{3}^{\prime}(\beta w)\\\ k_{1}\beta
Lf_{1}^{\prime}(\beta w)&k_{2}\beta Lf_{2}^{\prime}(\beta w)&k_{3}\beta
Lf_{3}^{\prime}(\beta
w)&\beta^{2}\sum_{i=1}^{3}\left[k_{i}\phi_{i}f_{i}^{\prime\prime}(\beta
w)-k_{i}Lf_{i}^{\prime 2}(\beta w)\right]\end{bmatrix}$ (24)

We set our aircraft to fly at constant airspeed (around $12m/s$) while flying
at constant altitude; therefore, we have a bounded speed $v$ (estimated
onboard with an Inertial Navigation System) when we account for the wind. For
tracking 3D paths, the aircraft will nose down or change the propeller’s
r.p.m.; nevertheless, the airspeed is also bounded between $9m/s$ and $16m/s$.
Note that both ground and airspeed are not control/guiding signals; therefore,
we do not face any saturation problems regarding these variables.

###  VI-C Accommodating the guidance to the aircraft’s dynamics

An arbitrary function $\phi_{i}(\cdot)$ in (17), which depends on a specific
parametrization $f_{i}(\cdot)$, may result in a highly sensitive coordinate
$w$. This can lead to considerable vibrations of the guidance signals, due to
noisy sensor readings or disturbances of the position, that cannot be tracked
effectively by the aircraft.

We propose two approaches, which can be combined to mitigate this practical
effect. The first one is to re-parameterize the equations for the 3D desired
path ${\mathcal{P}}^{\mathrm{phy}}$; this does not affect the convergence
result. Suppose ${\mathcal{P}}^{\mathrm{phy}}$ is re-parameterized by

$x=f_{1}\big{(}g(w)\big{)},\quad y=f_{2}\big{(}g(w)\big{)},\quad
z=f_{3}\big{(}g(w)\big{)},$

where $g:\mathbb{R}\to\mathbb{R}$ is a smooth bijection with nonzero
derivative (i.e., $\frac{\mathrm{d}g}{\mathrm{d}w}(w)\neq 0$ for all
$w\in\mathbb{R}$). A simple example of $g$ is $g(w)=\beta w$, where $\beta$ is
a positive constant. This is adopted for the experiments. Let $\phi_{1}$,
$\phi_{2}$, and $\phi_{3}$ be chosen as in (17), then the first term of the
higher-dimensional vector field becomes (for simplicity, the arguments are
omitted)

$\times(\nabla\phi_{1},\nabla\phi_{2},\nabla\phi_{3})=-\begin{bmatrix}\frac{\mathrm{d}f_{1}}{\mathrm{d}g}\cdot\frac{\mathrm{d}g}{\mathrm{d}w}\\\
\frac{\mathrm{d}f_{2}}{\mathrm{d}g}\cdot\frac{\mathrm{d}g}{\mathrm{d}w}\\\
\frac{\mathrm{d}f_{3}}{\mathrm{d}g}\cdot\frac{\mathrm{d}g}{\mathrm{d}w}\\\
1\end{bmatrix}.$ (25)

To reduce the effect of the “virtual speed” from the fourth coordinate of
(25), the “gain” $\frac{\mathrm{d}g}{\mathrm{d}w}$ can be chosen large such
that
$(\frac{\mathrm{d}f_{1}}{\mathrm{d}g}\cdot\frac{\mathrm{d}g}{\mathrm{d}w})^{2}+(\frac{\mathrm{d}f_{2}}{\mathrm{d}g}\cdot\frac{\mathrm{d}g}{\mathrm{d}w})^{2}+(\frac{\mathrm{d}f_{3}}{\mathrm{d}g}\cdot\frac{\mathrm{d}g}{\mathrm{d}w})^{2}\gg
1$, which implies that
$\left\lVert\nabla\phi_{1}\times\nabla\phi_{2}\times\nabla\phi_{3}\right\rVert\approx\left|\frac{\mathrm{d}g}{\mathrm{d}w}\right|\sqrt{(\frac{\mathrm{d}f_{1}}{\mathrm{d}g})^{2}+(\frac{\mathrm{d}f_{2}}{\mathrm{d}g})^{2}+(\frac{\mathrm{d}f_{3}}{\mathrm{d}g})^{2}}$.
However, from the analytic expression of the vector field

$\raisebox{2.0pt}{$\chi$}=\begin{bmatrix}-\frac{\mathrm{d}g}{\mathrm{d}w}\cdot\frac{\mathrm{d}f_{1}}{\mathrm{d}g}-k_{1}\phi_{1}\\\
-\frac{\mathrm{d}g}{\mathrm{d}w}\cdot\frac{\mathrm{d}f_{2}}{\mathrm{d}g}-k_{2}\phi_{2}\\\
-\frac{\mathrm{d}g}{\mathrm{d}w}\cdot\frac{\mathrm{d}f_{3}}{\mathrm{d}g}-k_{3}\phi_{1}\\\
-1+\frac{\mathrm{d}g}{\mathrm{d}w}\left(k_{1}\phi_{1}\frac{\mathrm{d}f_{1}}{\mathrm{d}g}+k_{2}\phi_{2}\frac{\mathrm{d}f_{2}}{\mathrm{d}g}+k_{3}\phi_{3}\frac{\mathrm{d}f_{3}}{\mathrm{d}g}\right)\end{bmatrix},$

one observes that, when $\left\lVert(\phi_{1},\phi_{2},\phi_{3})\right\rVert$
is large, (i.e., the aircraft is far from the desired path), the additional
coordinate of the vector has also been enlarged approximately by a factor of
$\frac{\mathrm{d}g}{\mathrm{d}w}$. Thus, the “gain”
$\lvert\frac{\mathrm{d}g}{\mathrm{d}w}\rvert$ should not be chosen too large.

The second approach is to scale down the functions $\phi_{i}$. That is, the
equations (17) are changed to

$\phi_{i}(x,y,z,w)=L(x_{i}-f_{i}(w)),i=1,2,3,$

where $L\in(0,1)$. The corresponding 3D vector field is thus changed to

$\tilde{\raisebox{2.0pt}{$\chi$}}=L\begin{bmatrix}-L^{2}\frac{\mathrm{d}f_{1}}{\mathrm{d}w}-k_{1}\phi_{1}\\\
-L^{2}\frac{\mathrm{d}f_{2}}{\mathrm{d}w}-k_{2}\phi_{2}\\\
-L^{2}\frac{\mathrm{d}f_{3}}{\mathrm{d}w}-k_{3}\phi_{3}\\\
-L^{2}+k_{1}\phi_{1}\frac{\mathrm{d}f_{1}}{\mathrm{d}w}+k_{2}\phi_{2}\frac{\mathrm{d}f_{1}}{\mathrm{d}w}+k_{3}\phi_{3}\frac{\mathrm{d}f_{1}}{\mathrm{d}w}\end{bmatrix}.$

The new guiding vector field is scaled down; thus, it helps to lower the
sensitivity of the additional coordinate $w$.

###  VI-D The 2D trefoil curve

We start with following a 2D self-intersected desired path, the trefoil curve,
at a constant altitude $z_{o}=50m$ over the ground level. The parametric
equations of the trefoil curve are given by

$\displaystyle f_{1}(w)$ $\displaystyle=\cos(\beta w\,\omega_{1})(a\cos(\beta
w\,\omega_{2})+b)$ $\displaystyle f_{2}(w)$ $\displaystyle=\sin(\beta
w\,\omega_{1})(a\cos(\beta w\,\omega_{2})+b)$ $\displaystyle f_{3}(w)$
$\displaystyle=0,$

where we have set $\beta=\frac{\mathrm{d}g}{\mathrm{d}w}=0.45$ (the “gain”
introduced in Section VI-C), $\omega_{1}=0.02$, $\omega_{2}=0.03$, $a=80$, and
$b=160$. In order to fit into the available flying space, these parametric
equations have been rotated by $\alpha$ and shifted adequately by
$(x_{o},y_{o})$ in the autopilot, i.e.,

$\displaystyle f_{1}^{*}(w)$
$\displaystyle=\cos(\alpha)f_{1}(w)-\sin(\alpha)f_{2}(w)+x_{o}$ $\displaystyle
f_{2}^{*}(w)$ $\displaystyle=\sin(\alpha)f_{1}(w)+\cos(\alpha)f_{2}(w)+y_{o}$
$\displaystyle f_{3}^{*}(w)$ $\displaystyle=f_{3}(w)+z_{o}.$

Note that the same affine transformation must be done for both
$f_{i}^{\prime}$ and $f_{i}^{\prime\prime}$ (needed for the Jacobian of $\chi$
as we will see shortly). In particular, for the presented flight, we have set
$x_{o}=79$, $y_{o}=-68.10$ and $z_{o}=50$ meters and $\alpha=0$. We set the
scaling factor $L=0.1$ for the construction of $\tilde{\phi}_{i}$ as in
Section VI-C, and we have chosen the gains $k_{1}=k_{2}=k_{3}=0.002$. We
finally set $k_{\theta}=1$ for the control/guidance signal $u_{\theta}$ in
Proposition 3.

Note that for computing all the control signals (22), we need $f_{i}(w)$ and
their derivatives $f_{i}^{\prime}(w)$ and $f_{i}^{\prime\prime}(w)$ with
respect to $w$. For the sake of completeness, we provide the Jacobian
$J(\raisebox{2.0pt}{$\chi$}^{p})$ in (22a) which is given by

$J(\raisebox{2.0pt}{$\chi$}^{p})=FJ(\hat{\raisebox{2.0pt}{$\chi$}})=F(I-\hat{\raisebox{2.0pt}{$\chi$}}\hat{\raisebox{2.0pt}{$\chi$}}^{\top})J(\raisebox{2.0pt}{$\chi$})/\left\lVert\raisebox{2.0pt}{$\chi$}\right\rVert,$
(26)

where $F=\left[\begin{smallmatrix}1&0&0&0\\\ 0&1&0&0\end{smallmatrix}\right]$,
and $J(\raisebox{2.0pt}{$\chi$})$ is shown in (24). We show the flight results
for the trefoil curve in Figure 4.

Figure 4: Flight results of the trefoil path tracking at a constant altitude.
On the low-right side, we show a screenshot of the ground control station,
where the green circle denotes the _stand-by_ trajectory just before starting
the tracking of the trefoil. On top, we show the first minute of the flight
once the tracking of the trefoil started. All the tracking is done at the
constant altitude of $50$ meters, and the evolution of the coordinate $z$ is
thus omitted. We can see how the blue dot, i.e., $(f_{1}(w),f_{2}(w))$ goes
and _waits_ for the aircraft at time $t=321$. Afterward, the vehicle converges
to the desired path as the first two low-left plots indicate with $\phi_{1}$
and $\phi_{2}$ staying around $0$ meters. The third plot shows the evolution
of the virtual coordinate $w$. We can see how $w$ does not follow a constant
growth, but it varies since it is in closed-loop with the position of the
aircraft to facilitate the convergence to the path. Coming back to the ground
control station, we can see in blue color the described 2D trajectory of the
aircraft, showing one period of the trefoil.

###  VI-E The 3D Lissajous curve

We consider the 3D Lissajous curve described as below:

$\displaystyle f_{1}(w)$ $\displaystyle=c_{x}\cos(\beta w\,\omega_{x}+d_{x})$
$\displaystyle f_{2}(w)$ $\displaystyle=c_{y}\cos(\beta w\,\omega_{y}+d_{y})$
$\displaystyle f_{3}(w)$ $\displaystyle=c_{z}\cos(\beta w\,\omega_{z}+d_{z}),$

where we have set $\beta=\frac{\mathrm{d}g}{\mathrm{d}w}=0.01$,
$\omega_{x}=1$, $\omega_{y}=\omega_{z}=2$, $c_{x}=c_{y}=225$, $c_{z}=-20$,
$d_{x}=d_{z}=0$, and $d_{y}=\pi/2$. This selection of parameters gives us an
_eight-shaped_ desired path that is bent along the vertical axis. As with the
trefoil curve, we have added an affine transformation of
$f_{i}(w),f_{i}^{\prime}(w)$ and $f_{i}^{\prime\prime}(w)$ in the autopilot to
fit the Lissajous curve into the available flying space. In particular, we
have set $x_{o}=79,y_{o}=-68.10,z_{o}=50,\alpha=0.66$). Finally, for the
construction of $\tilde{\phi}_{i}$, we have chosen $L=0.1$,
$k_{1}=k_{2}=0.002$ and $k_{3}=0.0025$. We finally set $k_{\theta}=1$ for the
control/guiding signal $u_{\theta}$ in Proposition 3. We show the flight
results in in Figure 5.

Figure 5: Flight results of the 3D Lissajous path tracking. On the low-right
side, we show a screenshot of the ground control station, where the green
circle denotes the _stand-by_ trajectory just before starting the tracking of
the Lissajous figure. In the 3D plots on top, we show the first seconds of the
flight once the tracking of the Lissajous figure started. We can see how the
blue dot, i.e., $(f_{1}(w),f_{2}(w),f_{3}(w))$ travels quickly and _waits_ for
the aircraft at time $t=204$. This quick traveling of $w$ can be seen at the
beginning of the fourth plot in the bottom row. Afterward, the vehicle stays
on the desired path as the three bottom-left plots indicate with
$\phi_{1},\phi_{2}$, and $\phi_{3}$ staying around $0$ meters. The aircraft
has been trimmed to fly at a constant altitude. However, this path requires
the vehicle to track a sinusoidal ascending/descending trajectory, and any
disturbance (e.g., changing wind) makes the aircraft sensitive to track a
climbing/descending speed different from zero accurately. In addition, the
Lissajous curve demands aggressive turnings a bit beyond the capabilities of
the aircraft, in particular, when the aircraft is descending and achieving
maximum speed. We can see how $w$ does not follow a constant growth, but it
varies since it is in closed-loop with the position of the aircraft to
facilitate the convergence to the path. Coming back to the ground control
station, we can see in blue color the described 2D (projected) trajectory of
the aircraft, showing one period of the Lissajous figure. In particular, we
can see that the aircraft passes by the _middle waypoint_ corresponding to the
highest point of the desired track.

##  VII Discussion
Path following or trajectory tracking?

In this section, we show that our proposed higher-dimensional VF-PF algorithm
is an extension that combines elements from both conventional VF-PF algorithms
(e.g., [16], see Remark 3) and trajectory tracking algorithms (e.g., [47, p.
506]). While our generated guiding vector field is the standard output for the
path-following approach, we will argue that our algorithm can also be seen as
a fair extension of a trajectory tracking approach. Therefore, our algorithm,
to some extent, combines and extends elements from both approaches. For ease
of explanation and without loss of generality, we restrict our focus to a
physical planar desired path in $\mathbb{R}^{2}$; that is,
${\mathcal{P}}^{\mathrm{phy}}\subseteq\mathbb{R}^{2}$.

Compared to trajectory tracking algorithms, a similarity exists in the sense
that the additional coordinate $w$ in the proposed VF-PF algorithms acts like
the time variable in trajectory tracking algorithms. However, our approach is
an extension in the sense that the time-like variable is in fact state-
dependent. In trajectory tracking algorithms, a desired trajectory
$\big{(}x_{d}(t),y_{d}(t)\big{)}$ is given. Then, at any time instant $t$, the
algorithm aims to decrease the distance to the _desired trajectory point_
$\big{(}x_{d}(t),y_{d}(t)\big{)}$, which moves as time $t$ advances. Note that
the dynamics of the _desired trajectory point_
$\big{(}x_{d}(t),y_{d}(t)\big{)}$ is _open-loop_ in the sense that it does not
depend on the current states of the robot, but only depends on time $t$. From
(17), if we let $\phi_{i}=0$, $i=1,2$, then we may call the point
$\big{(}f_{1}(w(\xi(t))),f_{2}(w(\xi(t)))\big{)}$ the _guiding point_ , since
it always stays on the desired path and may be regarded as the counterpart of
the _desired trajectory point_ in trajectory tracking algorithms. But as we
will show later, the _guiding point_ is essentially different from the
_desired trajectory point_. Note that the _guiding point_
$\big{(}f_{1}(w(\xi(t))),f_{2}(w(\xi(t)))\big{)}$ in our VF-PF algorithm
depends on the evolution of the additional coordinate $w\big{(}\xi(t)\big{)}$,
of which the dynamics is state-dependent as shown in (21). This might be
roughly regarded as a _closed-loop_ version of the _desired trajectory point_.
An intuitive consequence of this difference is that the _desired trajectory
point_ $\big{(}x_{d}(t),y_{d}(t)\big{)}$ in trajectory tracking algorithms
always moves _unidirectionally_ along the desired trajectory as $t$
monotonically increases, while the _guiding point_ can move _bidirectionally_
along the desired path, subject to the current state (i.e., the robot
position). In fact, when the initial position of the _guiding point_
$\big{(}f_{1}(w(\xi(0))),f_{2}(w(\xi(0)))\big{)}$ is far from the initial
position of the robot, the _guiding point_ “proactively” moves towards the
robot along the desired path to accelerate the path-following process. This
feature, along with better robustness against perturbation in some cases, are
experimentally studied in our previous work [33, Section VII]. To illustrate
this closed-loop feature more intuitively, after the robot has successfully
followed the desired path, we manually move the robot far away from the
desired path and keep it stationary (to mimic the situation of erroneous
localization and operation failure of the robot). As is clear in the
supplementary video121212https://www.youtube.com/watch?v=jxWPsm0g-Ro.,
although the robot is kept stationary, the guiding point
$\big{(}f_{1}(w(\xi(t))),f_{2}(w(\xi(t)))\big{)}$ can still move in the
reverse direction to approach the robot along the desired path such that the
path-following error decreases, and the guiding point eventually stops at some
place on the desired path. After that, the guiding point does not move until
the robot is released to move again.

In existing VF-PF algorithms, a two-dimensional vector field on
$\mathbb{R}^{2}$ is created for guiding the robot movement (see Remark 3).
However, as we aim to create a higher-dimensional (i.e., three-dimensional)
vector field, our approach can be roughly regarded as utilizing an infinite
number of layers of projected two-dimensional vector fields, and thus might be
seen as a _dynamic_ two-dimensional vector field. The dynamic property is due
to the dynamics of the additional coordinate $w$. For example, consider a
circular desired path parameterized by

$x=f_{1}(w)=\cos(4w)\quad y=f_{2}(w)=\sin(4w),$

where $w\in\mathbb{R}$ is the parameter. In conventional VF-PF algorithms, a
2D vector field can be created, as shown in Fig. 1, but there exists a
singular point at the center of the circle. Nevertheless, using our approach,
we can generate a singularity-free 3D vector field, as illustrated in Fig. 6.
For clarity of visualization, we plot the 3D vectors, which originate from
three planes where the $w$ values are $0,0.6$, and $1.4$, respectively. For
each value of the additional coordinate $w$, we can obtain a projected 2D
vector field, as shown in Fig. 7. Therefore, we can observe that these 2D
vector fields change dynamically as $w$ varies. As a result of the dynamics of
$w$, the _guiding point_ $\big{(}f_{1}(w(\xi(t))),f_{2}(w(\xi(t)))\big{)}$
moves along the 2D desired path (not necessarily uni-directionally). Again, we
note that this point is not the same as the _desired trajectory point_ in
trajectory tracking algorithms since the integral curves of the 2D vector
field do not converge to this point, as can be seen graphically from Fig. 7 or
analytically from the expression of the vector field in (6): the second term
leads to convergence to the _guiding point_ , while the first term “deviates”
this convergence, since it controls the propagation along the _higher-
dimensional_ desired path.

In many existing VF-PF algorithms, the desired path is usually not
parameterized but is described by the intersection of (hyper)surfaces, while
the latter case might be restrictive in describing more complicated desired
paths. However, our approach enables the possibility to use a parameterized
desired path directly in the design of a higher-dimensional vector field. Our
approach thus extends the flexibility of conventional VF-PF algorithms. The
desired path can now be described by either the intersection of
(hyper)surfaces or parameterized functions. In the latter case, the parametric
equations can be easily converted using (16), (17) and (18) and leads to a
higher-dimensional desired path and singularity-free guiding vector field.
Theoretically, the parametrization is not instrumental, since it is only
utilized to derive the expressions of functions $\phi_{i}$, of which the zero-
level sets are interpreted as (hyper)surfaces. The subsequent derivation of
the vector field is based on $\phi_{i}$, independent of the specific
parametrization of the desired path.

Figure 6: Three layers of the 3D vector field corresponding to a circle. The
solid line is the 2D desired path while the dashed line is the corresponding
3D (unbounded) desired path. Three layers of the 3D vector field evaluated at
$w=0,0.6,1.4$ respectively are illustrated.  Figure 7: The projected 2D vector
field corresponding to $w=0,0.6,1.4$ respectively. The solid line is the
projected 2D desired path. The solid dots represent the guiding point
$\big{(}cos(4w),sin(4w)\big{)}$.

##  VIII Conclusion and Future work

In this paper, we first show that the integral curves of a time-invariant
continuously differential vector field as in (4) cannot guarantee global
converge to desired paths which are simple closed (i.e., homeomorphic to the
unit circle) or self-intersected. Motivated by this general topological
result, we propose a novel approach to create unbounded desired paths from
simple closed or self-intersected ones, and construct a singularity-free
higher-dimensional guiding vector field. One of the advantages of this
approach is that global convergence to the desired paths, which can be even
self-intersected, is now rigorously guaranteed. This is achieved by the
introduction of a transformation operator and the extended dynamics. Another
advantage is that, given a parameterized desired path, we can easily describe
the (hyper)surfaces as the zero-level set of some implicit functions, and then
the proposed vector field on a higher-dimensional space can be directly
constructed. This increases the applicability of conventional VF-PF
algorithms. In addition, we highlight five features of our approach, with
rigorous theoretic guarantees. We also show that our approach is a combined
extension of both conventional VF-PF algorithms and trajectory tracking
algorithms. Finally, we conduct outdoor experiments with a fixed-wing aircraft
under wind perturbation to validate the theoretical results and demonstrate
the practical effectiveness for complex robotic systems.

Due to the additional coordinates of the vector field, it is difficult for a
robot to follow the desired path at a constant speed. This may be solved by
replacing the normalization of the vector field
$\hat{\raisebox{2.0pt}{$\chi$}}$ by the “partial normalization”; i.e.,
$\tilde{\raisebox{2.0pt}{$\chi$}}:=\raisebox{2.0pt}{$\chi$}/\left\lVert\raisebox{2.0pt}{$\chi$}_{p}\right\rVert$,
where $\raisebox{2.0pt}{$\chi$}_{p}$ is the vector composed of only the first
$(n-1)$ components of $\chi$ (the $n$-th component is the additional
coordinate). It is our future work to use the singularity-free guiding vector
field to solve collision avoidance problems; preliminary results are reported
in [41]. We are also interested in determining the parameterization of the
desired path to achieve the optimal path-following accuracy.

## Acknowledgments

Weijia Yao and Bohuan Lin are funded by the China Scholarship Council. The
work of Héctor Garcia de Marina is supported by the grant Atraccion de Talento
2019-T2/TIC-13503 from the Government of the Autonomous Community of Madrid.
The work of Hector is also partially supported by the Spanish Ministry of
Science and Innovation under research Grant RTI2018-098962-B-C21. The work of
Cao was supported in part by the European Research Council (ERC-CoG-771687)
and the Netherlands Organization for Scientific Research (NWO-vidi-14134). We
thank Murat Bronz and Gautier Hattenberger from the _Paparazzi_ project for
their assistance and feedback. We also thank Yangguang Yu and Kyle Miller for
the technical discussions.

##  Appendix A Proofs and Supplementary Material for Section II

###  A-A Proof of Lemma 2

###### Proof.

First, it is easy to see that for any point
$\xi\in\mathcal{P}\cup\mathcal{C}$, we have $\xi\in\mathcal{M}$, thus
$\mathcal{P}\cup\mathcal{C}\subseteq\mathcal{M}$. Second, for any point
$\xi^{\prime}\in\mathcal{M}$, it follows that
$N(\xi^{\prime})Ke(\xi^{\prime})=\sum_{i=1}^{n-1}k_{i}\phi_{i}(\xi^{\prime})\nabla\phi_{i}(\xi^{\prime})=0$.
If $e(\xi^{\prime})=0$, then $\xi^{\prime}\in\mathcal{P}$. If
$e(\xi^{\prime})\neq 0$, then the former equation implies that
$\nabla\phi_{i}(\xi^{\prime}),\;i=1,\dots,n-1$, are linearly dependent
(recalling that $k_{i}>0$); hence the first term of the vector field becomes
zero (i.e., $\nabla_{\times}\phi(\xi)=0$). Since $\xi^{\prime}\in\mathcal{M}$,
the second term of the vector field is also zero, thus the vector field
$\raisebox{2.0pt}{$\chi$}(\xi^{\prime})=0$ and $\xi^{\prime}\in\mathcal{C}$.
The reasoning shows that $\mathcal{M}\subseteq\mathcal{P}\cup\mathcal{C}$.
Combining $\mathcal{P}\cup\mathcal{C}\subseteq\mathcal{M}$ and
$\mathcal{M}\subseteq\mathcal{P}\cup\mathcal{C}$, it is indeed true that
$\mathcal{M}=\mathcal{P}\cup\mathcal{C}$. ∎

###  A-B Proof of Lemma 3

###### Proof.

By the definition of the desired path $\mathcal{P}$ in (7), it is obvious that
the desired path is the inverse image of $e(\bm{0})$; that is,
$\mathcal{P}=e^{-1}(\bm{0})$. Under Assumption 1, for any $\xi\in\mathcal{P}$,
the Jacobian matrix $\mathrm{d}\,e/\mathrm{d}\,\xi=N^{\top}(\xi)$ is of full
row rank. Thus the result follows directly from [24, Corrollary 5.14]. ∎

###  A-C Mathematical Formulation of Assumption 2

To keep the paper self-contained, we restate Assumption 2 using precise
mathematical descriptions. The first part is: for any given constant
$\kappa>0$, it follows that
$\inf\\{||e(\xi)||:\mathrm{dist}(\xi,\mathcal{P})\geq\kappa\\}>0$. Under this
assumption, it can be proved that as the norm of the path-following error
$\left\lVert e(\cdot)\right\rVert$ vanishes, the distance to the desired path
$\mathrm{dist}(\cdot,\mathcal{P})$ also vanishes; that is,
$\lim\limits_{t\rightarrow\infty}||e\big{(}p(t)\big{)}||=0\implies\lim\limits_{t\rightarrow\infty}\mathrm{dist}(p(t),\mathcal{P})=0$
[13, 48]. The second part of the assumption is: for any given constant
$\kappa>0$, it holds that $\inf\left\\{\left\lVert
N(\xi)Ke(\xi)\right\rVert:\mathrm{dist}(\xi,\mathcal{M})\geq\kappa\right\\}>0$.
Similarly, this implies
$\lim\limits_{t\rightarrow\infty}||N\big{(}p(t)\big{)}Ke\big{(}p(t)\big{)}||=0\implies\lim\limits_{t\rightarrow\infty}\mathrm{dist}(p(t),\mathcal{M})=0$
[13, 48].

##  Appendix B Proofs for Section III

###  B-A Proof of Proposition 1

###### Proof.

Since $c\in\mathcal{P}$ is a crossing point, we have $e(c)=0$, and thus the
vector field at the crossing point is simplified to
$\raisebox{2.0pt}{$\chi$}(c)=\nabla_{\times}\phi(c)$ in view of (6). Next we
show that the gradients at the crossing point $\nabla\phi_{i}(c)$ are linearly
dependent, and hence $\raisebox{2.0pt}{$\chi$}(c)=0$. Suppose, on the
contrary, the gradients are not linearly dependent. Then we can use the
implicit function theorem [49] to conclude that there is a unique curve in a
neighborhood of $c$ satisfying $e(\xi)=0$, where $\xi\in\mathbb{R}^{n}$. But
this contradicts the fact that $\mathcal{P}$ is self-intersected. Therefore,
the gradients at the crossing point are indeed linearly dependent. ∎

###  B-B Proof of Proposition 2

Given an autonomous differential equation $\dot{x}(t)=f(x(t))$, where $f$ is
continuously differentiable in $x$, and let $t\mapsto\Psi(t,x_{0})$ be the
solution to the differential equation with the initial condition
$\Psi(0,x_{0})=x_{0}$, then $\Psi$ is a flow [43]. In the literature, the
notation $\Psi^{t}(x_{0})$, which is adopted in the sequel, is often used in
place of $\Psi(t,x_{0})$. To assist the proof of Proposition 2, we state a
more general result in the following lemma regarding any time-invariant
autonomous system that admits a (locally) asymptotically stable limit cycle.
Note that similar to the definition of Lyapunov stability of an equilibrium
point [34, Chapter 4], a limit cycle $\mathcal{L}$ is (locally) asymptotically
stable if for every neighborhood $\mathcal{U}\supseteq\mathcal{L}$ of the
limit cycle $\mathcal{L}$, there exists a _smaller_ neighborhood
$\mathcal{V}\subseteq\mathcal{U}$, such that every trajectory starting from
$\mathcal{V}$ always stays within $\mathcal{U}$ and $\mathcal{L}$ is locally
attractive. The proof of this lemma involves notions about topological and
differential manifolds [24, 31], which are briefly explained here. If
$\mathcal{V}$ and $\mathcal{W}$ are vector spaces, a bijective linear map
$T:\mathcal{V}\to\mathcal{W}$ is called an _isomorphism_. If this isomorphism
exists, then $\mathcal{V}$ and $\mathcal{W}$ are called _isomorphic_ , denoted
by $\mathcal{V}\cong\mathcal{W}$. Suppose $\mathcal{X},\mathcal{Y}$ are
topological spaces. An _open cover_ of $\mathcal{X}$ is a collection of open
subsets of $\mathcal{X}$ whose union is $\mathcal{X}$. A _diffeomorphism_
$f:\mathcal{X}\to\mathcal{Y}$ is a smooth bijection that has a smooth inverse.
If there exists a diffeomorphism between $\mathcal{X}$ and $\mathcal{Y}$, then
$\mathcal{X}$ and $\mathcal{Y}$ are _diffeomorphic_. Let
$f,g:\mathcal{X}\to\mathcal{Y}$ be continuous maps. _A homotopy from $f$ to
$g$_ is a continuous map $H:\mathcal{X}\times[0,1]\to\mathcal{Y}$ such that
$H(x,0)=f(x)$ and $H(x,1)=g(x)$ for all $x\in\mathcal{X}$. If there exists
such a homotopy, then $f$ and $g$ are _homotopic_ , denoted by $f\simeq g$.
Let $\mathcal{A}\subseteq\mathcal{X}$. The _inclusion map of $\mathcal{A}$ in
$\mathcal{X}$_ is $\iota_{\mathcal{A}}:\mathcal{A}\to\mathcal{X}$ defined by
$\iota_{\mathcal{A}}(x)=x$ for $x\in\mathcal{A}$. A continuous map
$r:\mathcal{X}\to\mathcal{A}$ is a _retraction_ if the restriction of $r$ to
$\mathcal{A}$ is the identity map of $\mathcal{A}$, or equivalently if
$r\circ\iota_{\mathcal{A}}=\mathrm{id}_{\mathcal{A}}$, where
$\iota_{\mathcal{A}}$ is the inclusion map of $\mathcal{A}$ in $\mathcal{X}$,
and $\mathrm{id}_{\mathcal{A}}$ is the identity map of $\mathcal{A}$. In this
case, $\mathcal{A}$ is called a _retract of $\mathcal{X}$_. The space
$\mathcal{X}$ is called _contractible_ if the identity map
$\mathrm{id}_{\mathcal{X}}$ of $\mathcal{X}$ is homotopic to a constant map.
Intuitively, this means that the whole space $\mathcal{X}$ can continuously
shrink to a point. For the definition of a tubular neighborhood and an
introduction to fundamental groups and homomorphisms of fundamental groups
induced by continuous maps, see [24, pp. 137-139] and [31, Chapter 7]
respectively.

######  Lemma 5 (Asymptotically stable limit cycles are not GAS).

Consider an autonomous differential equation $\dot{x}(t)=f(x(t))$, where
$f:\mathbb{R}^{n}\to\mathbb{R}^{n}$ is continuously differentiable in $x$.
Suppose there is a (locally) asymptotically stable limit cycle
$\mathcal{L}\subseteq\mathbb{R}^{n}$, then global convergence of trajectories
to the limit cycle is not possible; namely, the domain of attraction of the
limit cycle cannot be $\mathbb{R}^{n}$. In other words, the limit cycle cannot
be globally asymptotically stable (GAS) in $\mathbb{R}^{n}$.

###### Proof.

We prove by contradiction: Suppose that global convergence to the limit cycle
$\mathcal{L}$ holds. Since the limit cycle is compact, it is an embedded
submanifold in $\mathbb{R}^{n}$ [24, Proposition 5.21]. So we can take a
tubular neighborhood $\mathcal{O}\supseteq\mathcal{L}$ of the limit cycle [24,
Theorem 6.24]. Then due to the asymptotic stability of the limit cycle, there
exists a smaller neighborhood $\mathcal{U}\subseteq\mathcal{O}$ of the limit
cycle such that every trajectory starting from $\mathcal{U}$ will remain
within the tubular neighborhood $\mathcal{O}$ perpetually. Since the limit
cycle is compact, we can find a closed ball
$\overline{\mathcal{B}}\subseteq\mathbb{R}^{n}$ centered at
$\bm{0}\in\mathbb{R}^{n}$ sufficiently large such that the limit cycle lies in
its interior (i.e., $\mathcal{L}\subseteq\mathcal{B}$). Due to the global
convergence assumption, for any point $w\in\overline{\mathcal{B}}$, there
exists a time instant $T_{w}>0$ such that $\Psi^{T_{w}}(w)\in\mathcal{U}$,
where $\Psi$ denotes the flow of the differential equation $\dot{x}=f(x)$. Due
to the continuous dependence on initial conditions [34, Theorem 3.5], there
exists an open set $\mathcal{V}_{w}\ni w$, such that
$\Psi^{T_{w}}(\mathcal{V}_{w})\subseteq\mathcal{U}$. Therefore, according to
the uniqueness of solutions to the differential equation [34, Theorem 3.1] and
the asymptotic stability discussed before, we further have
$\Psi^{t}(\mathcal{V}_{w})\subseteq\mathcal{O}$ for all $t\geq T_{w}$. Thus,
for every point $w\in\overline{\mathcal{B}}$, we can associate an open set
$\mathcal{V}_{w}$ and a time instant $T_{w}$ as discussed before. Since
$\mathcal{D}:=\\{\mathcal{V}_{w}\subseteq\mathbb{R}^{n}:w\in\overline{\mathcal{B}}\\}$
is an open cover of the _compact_ ball $\overline{\mathcal{B}}$, there exists
a finite number of points $w_{i}\in\overline{\mathcal{B}},i=1,\dots,k$, and
$\mathcal{V}_{w_{i}}\in\mathcal{D}$, such that
$\bigcup_{i=1}^{k}\mathcal{V}_{w_{i}}\supseteq\overline{\mathcal{B}}$ [44,
Theorem 1.5.8]. Thus, we can take $T>\max_{i=1,\dots,k}\\{T_{w_{i}}\\}$, and
therefore, we have $\Psi^{T}(\mathcal{B})\subseteq\mathcal{O}$.

Let $r:\mathcal{O}\to\mathcal{L}$ be a retraction131313The existence of $r$ is
guaranteed by Proposition 6.25 in [24].  of $\mathcal{O}$ onto $\mathcal{L}$;
i.e., $r\circ i_{\mathcal{L}}=\mathrm{id}_{\mathcal{L}}$, where
$i_{\mathcal{L}}:\mathcal{L}\to\mathcal{O}$ is the inclusion map of
$\mathcal{L}$ in $\mathcal{O}$ and $\mathrm{id}$ is the identity map. Now let
$i_{\mathcal{L}}^{\prime}:\mathcal{L}\to\overline{\mathcal{B}}$ be another
inclusion map, and note that for any $t\in\mathbb{R}$, $\Psi^{t}(\cdot)$ is a
diffeomorphism of $\mathcal{L}$ [43, p. 13]. Then it is easy to check that
$\mathrm{id}_{\mathcal{L}}=\Psi^{-T}\circ r\circ\Psi^{T}\circ
i_{\mathcal{L}}^{\prime}$, where we view $\Psi^{T}$ as a map from
$\overline{\mathcal{B}}$ to $\mathcal{O}$ and $\Psi^{-T}$ a map from
$\mathcal{L}$ to $\mathcal{L}$. It is conventional to use $(\cdot)_{*}$ and
$\pi_{1}(\cdot)$ to denote the homomorphism and the fundamental group
respectively. Then

$\begin{split}(\mathrm{id}_{\mathcal{L}})_{*}&=(\Psi^{-T}\circ
r\circ\Psi^{T}\circ i_{\mathcal{L}}^{\prime})_{*}\\\
&=(\Psi^{-T})_{*}\circ(r)_{*}\circ(\Psi^{T})_{*}\circ(i_{\mathcal{L}}^{\prime})_{*},\end{split}$
(27)

where
$(\mathrm{id}_{\mathcal{L}})_{*}:\pi_{1}(\mathcal{L})\to\pi_{1}(\mathcal{L})$,
$(i_{\mathcal{L}}^{\prime})_{*}:\pi_{1}(\mathcal{L})\to\pi_{1}(\overline{\mathcal{B}})$,
$(\Psi^{T})_{*}:\pi_{1}(\overline{\mathcal{B}})\to\pi_{1}(\mathcal{O})$,
$r_{*}:\pi_{1}(\mathcal{O})\to\pi_{1}(\mathcal{L})$ and
$(\Psi^{-T})_{*}:\pi_{1}(\mathcal{L})\to\pi_{1}(\mathcal{L})$ are the
homomorphisms of fundamental groups induced by the corresponding maps [24,
Proposition A.64, A.65]. Since $\overline{\mathcal{B}}$ is contractible and
$\pi_{1}(\overline{\mathcal{B}})\cong\\{0\\}$, where $\cong$ denotes the
isomorphic relation, both $(i_{\mathcal{L}}^{\prime})_{*}$ and
$(\Psi^{T})_{*}$ are zero morphisms, and so is the composition
$(\Psi^{-T})_{*}\circ(r)_{*}\circ(\Psi^{T})_{*}\circ(i_{\mathcal{L}}^{\prime})_{*}$.
But this contradicts with the left-hand side of (27), where
$(\mathrm{id}_{\mathcal{L}})_{*}$ is the identity map (and an isomorphism) of
$\pi_{1}(\mathcal{L})\cong\mathbb{Z}$. The contradiction implies that global
convergence is not possible. ∎

Based on Lemma 5, we can prove Proposition 2.

###### Proof of Proposition 2.

We consider the autonomous systems (10). Without loss of generality, we assume
that the flow of (10) is complete, since otherwise we can replace the vector
field $\chi$ by
$\raisebox{2.0pt}{$\chi$}/(1+\left\lVert\raisebox{2.0pt}{$\chi$}\right\rVert)$
without causing any difference of the topological properties of the flow that
will be discussed in the sequel [43, Proposition 1.14].

Given $\alpha>0$, we define a neighborhood of the desired path $\mathcal{P}$
by

$\mathcal{E}_{\alpha}=\\{\xi\in\mathbb{R}^{n}:\left\lVert
e(\xi)\right\rVert<\alpha\\}.$ (28)

Therefore, the value of $\left\lVert e(\cdot)\right\rVert$ encodes the
_distance_ to the desired path in view of the definition of $\mathcal{P}$ in
(7). From Lemma 1, we have
$N^{\top}\raisebox{2.0pt}{$\chi$}=N^{\top}(\nabla_{\times}\phi-
NKe)=-N^{\top}NKe$. We define a Lyapunov function candidate

$V(e)=\frac{1}{2}e^{\top}Ke,$ (29)

and take the time derivative of it, obtaining

$\begin{split}\dot{V}(e)&=\frac{1}{2}(\dot{e}^{\top}Ke+e^{\top}K\dot{e})\\\
&=\frac{1}{2}(\raisebox{2.0pt}{$\chi$}^{\top}NKe+e^{\top}KN^{\top}\raisebox{2.0pt}{$\chi$})\\\
&=-e^{\top}Qe=-\left\lVert NKe\right\rVert^{2}\leq 0,\end{split}$ (30)

where the $(n-1)\times(n-1)$ matrix

$Q(\xi)=K^{\top}N^{\top}(\xi)N(\xi)K$ (31)

is positive semidefinite. Based on the LaSalle’s invariance principle [34,
Theorem 4.4], one can show that the desired path $\mathcal{P}$ is the limit
cycle of (10) by construction, and that $\mathcal{P}$ is Lyapunov stable [22].
The claim then easily follows from Lemma 5. ∎

##  Appendix C Proofs for Section IV

###  C-A Proof of Lemma 4

###### Proof.

Due to the twice continuous differentiability of the transformation operator
$G_{l}$, the corresponding Jacobian matrix function
$\bm{\mathrm{D}}G_{l}=\partial G_{l}/\partial
x:\mathbb{R}^{m}\to\mathbb{R}^{m\times m}$ is locally Lipschitz continuous,
where $x$ is the argument of $G_{l}$. Therefore, the product of the vector
field $\chi$ and the Jacobian $\bm{\mathrm{D}}G_{l}$ are also locally
Lipschitz continuous. It follows that
$\big{(}\xi(t),{\xi}^{\mathrm{trs}}(t)\big{)}\in\mathbb{R}^{2m}$, where
${\xi}^{\mathrm{trs}}(t)=G_{l}\big{(}\xi(t)\big{)}$, is the unique solution to
(11) [43]. Recall that $l$ is the (global) Lipschitz constant of $G_{l}$. Fix
$t$, then

$\displaystyle\mathrm{dist}({\xi}^{\mathrm{trs}}(t),{\mathcal{A}}^{\mathrm{trs}})$
$\displaystyle=\inf\\{\left\lVert{\xi}^{\mathrm{trs}}(t)-p\right\rVert:p\in{\mathcal{A}}^{\mathrm{trs}}\\}$
$\displaystyle=\inf\\{\left\lVert
G_{l}\big{(}\xi(t)\big{)}-G_{l}(q)\right\rVert:q\in\mathcal{A}\\}$
$\displaystyle\leq\inf\\{l\left\lVert\xi(t)-q\right\rVert:q\in\mathcal{A}\\}$
$\displaystyle=l\cdot\mathrm{dist}(\xi(t),\mathcal{A}).$

Since $\xi(t)$ asymptotically converges to $\mathcal{A}$, we have
$\mathrm{dist}(\xi(t),\mathcal{A})\to 0$ as $t\to\infty$. In other words, for
any $\epsilon>0$, there exists a $T>0$, such that for all $t\geq T$,
$\mathrm{dist}(\xi(t),\mathcal{A})<\epsilon/l$; hence
$\mathrm{dist}({\xi}^{\mathrm{trs}}(t),{\mathcal{A}}^{\mathrm{trs}})\leq
l\cdot\mathrm{dist}(\xi(t),\mathcal{A})<\epsilon$. Therefore,
$\mathrm{dist}({\xi}^{\mathrm{trs}}(t),{\mathcal{A}}^{\mathrm{trs}})\to 0$ as
$t\to\infty$. Thus the transformed solution ${\xi}^{\mathrm{trs}}(t)$
asymptotically converges to the transformed set
${\mathcal{A}}^{\mathrm{trs}}$. ∎

###  C-B Proof of Corollary 1

###### Proof.

First consider the differential equation
$\dot{\xi}=\raisebox{2.0pt}{$\chi$}(\xi)$. Using the same Lyapunov function
candidate as (29) and the argument in the proof of Proposition 2, we have
$\dot{V}(e)\leq 0$. In addition, the norm of the first term of the scaled
vector field $s\hat{\raisebox{2.0pt}{$\chi$}}(\xi)$ is
$s\left\lVert\nabla_{\times}\phi\right\rVert/\left\lVert\raisebox{2.0pt}{$\chi$}\right\rVert$,
and it is obviously upper bounded in $\mathbb{R}^{m}$. Since the new vector
field $s\hat{\raisebox{2.0pt}{$\chi$}}(\xi)$ differs from the actual vector
field $\raisebox{2.0pt}{$\chi$}(\xi)$ only by the magnitude of each vector,
the two differential equations $\dot{\xi}=\raisebox{2.0pt}{$\chi$}(\xi)$ and
$\dot{\xi}=s\hat{\raisebox{2.0pt}{$\chi$}}(\xi)$ have the same phase portrait
in $\mathbb{R}^{m}\setminus\mathcal{C}$ [43, Proposition 1.14]. Therefore,
from the dichotomy convergence result proved in [22, Theorem 3], the solution
to $\dot{\xi}=s\hat{\raisebox{2.0pt}{$\chi$}}(\xi)$ will converge to either
$\mathcal{P}$ or $\mathcal{C}$ for initial conditions
$\xi(0)\in\mathbb{R}^{m}\setminus\mathcal{C}$.

Note that if the maximal prolonged time is $t^{*}<\infty$, then the solution
to $\dot{\xi}=s\hat{\raisebox{2.0pt}{$\chi$}}(\xi)$ must converge to the
singular set $\mathcal{C}$. This is shown by contradiction. Since
$\big{\lVert}\dot{\xi}\big{\rVert}=s<\infty$, $\xi^{*}:=\lim\limits_{t\to
t^{*}}\xi(t)=\xi(0)+\int_{0}^{t^{*}}\dot{\xi}(t)dt$ exists. Suppose
$\raisebox{2.0pt}{$\chi$}(\xi^{*})\neq 0$, then we can define the solution at
$t=t^{*}$, then it can be further prolonged to $t^{*}+\epsilon$ for some
$\epsilon>0$, contradicting that $t^{*}$ is the maximal prolonged time. This
shows that $\raisebox{2.0pt}{$\chi$}(\xi^{*})=0$ and thus the solution
converges to $\mathcal{C}$.

Finally, suppose $\xi(t)$ is the unique solution to the initial value problem
$\dot{\xi}(t)=s\hat{\raisebox{2.0pt}{$\chi$}}\big{(}\xi(t)\big{)},\;\xi(0)=\xi_{0}$,
then by Lemma 4, $\big{(}\xi(t),{\xi}^{\mathrm{trs}}(t)\big{)}$ is the
solution to (14). Therefore, the transformed trajectory
${\xi}^{\mathrm{trs}}(t)$ asymptotically converges to the transformed desired
path ${\mathcal{P}}^{\mathrm{trs}}$ as $t\rightarrow\infty$ or the transformed
singular set ${\mathcal{C}}^{\mathrm{trs}}$ as $t\rightarrow t^{*}$. ∎

##  Appendix D Theoretical Guarantees for Features in Section V

We use the set defined in (28) in the sequel, which corresponds to a
neighborhood of the higher-dimensional desired path
${\mathcal{P}}^{\mathrm{hgh}}$ (rather than that of the physical desired path
${\mathcal{P}}^{\mathrm{phy}}$). Note that the higher-dimensional desired path
${\mathcal{P}}^{\mathrm{hgh}}$ is always unbounded and non-self-intersecting
(i.e., homeomorphic to the real line), and thus an appealing result about the
_exponential convergence_ of the norm of the path-following error $\left\lVert
e(\cdot)\right\rVert$ can be concluded from [13, Theorem 2] with some
modifications, where one of the conditions of the original theorem is
automatically satisfied. Before we present the main results, we need a lemma
regarding some algebraic calculations as follows.

######  Lemma 6.

It holds that

$\det(N^{\top}N)=\left\lVert\nabla_{\times}\phi\right\rVert^{2},$

where $N$ and $\nabla_{\times}\phi$ are defined in (6).

###### Proof.

If $\nabla\phi_{1},\dots,\nabla\phi_{n-1}$ are linearly dependent, then we
have $\det(N^{\top}N)=\left\lVert\nabla_{\times}\phi\right\rVert^{2}=0$. Now
suppose $\nabla\phi_{1},\dots,\nabla\phi_{n-1}$ are linearly independent. In
this case, it follows that $\left\lVert\nabla_{\times}\phi\right\rVert\neq 0$.
So we can define the normalized vector
$v=\nabla_{\times}\phi/\left\lVert\nabla_{\times}\phi\right\rVert$. Due to
Lemma 1, we have

$\det(N^{\top}N)=\det\left(\begin{bmatrix}v^{\top}\\\
N^{\top}\end{bmatrix}\begin{bmatrix}v&N\end{bmatrix}\right)=\det\left(\begin{bmatrix}v^{\top}\\\
N^{\top}\end{bmatrix}\right)^{2}.$

According to the definition of $\nabla_{\times}\phi$ in (6), we have
$\left\lVert\nabla_{\times}\phi\right\rVert=\det\left(\begin{bmatrix}v^{\top}\\\
N^{\top}\end{bmatrix}\right)$. Then it follows that
$\det(N^{\top}N)=\left\lVert\nabla_{\times}\phi\right\rVert^{2}$. ∎

The result about the local and global exponential vanishing of the norm of the
path-following error now follows.

######  Proposition 4.

Define $\mathcal{E}_{\alpha}$ as in (28) for some $\alpha>0$, and let $\xi(t)$
be the solution to
$\dot{\xi}(t)=\raisebox{2.0pt}{$\chi$}\big{(}\xi(t)\big{)}$, where
$\raisebox{2.0pt}{$\chi$}:\mathbb{R}^{n}\to\mathbb{R}^{n}$ is the constructed
guiding vector field on $\mathbb{R}^{n}$ using the procedure in Section V-A.
If $|f_{i}^{\prime}(w)|$ is upper bounded in $\mathcal{E}_{\alpha}$ for
$i=1,\dots,n-1$, then there exists a positive constant
$\gamma^{\prime}<\alpha$, such that the norm of the path-following error
$\left\lVert e(t)\right\rVert$ locally exponentially converges to $0$ as
$t\to\infty$ when the initial condition
$\xi(0)\in\mathcal{E}_{\gamma^{\prime}}$. If $|f_{i}^{\prime}(w)|$ is upper
bounded globally, then the norm of the path-following error $\left\lVert
e(t)\right\rVert$ globally exponentially converges to $0$ as $t\to\infty$ .

###### Proof.

The proof uses the same idea as in [13, Theorem 2]. Using the construction of
$\phi_{i}$ as in Section V-A, the norm of the “tangential term”
$\left\lVert\nabla_{\times}\phi\right\rVert$ is bounded if and only if
$|f_{i}^{\prime}(w)|$ is bounded. This condition guarantees that the solution
is well defined for $t\to\infty$. In addition, since the last component of the
term $\nabla_{\times}\phi$ is a nonzero constant,
$\left\lVert\nabla_{\times}\phi\right\rVert$ is automatically bounded away
from zero globally in $\mathbb{R}^{n}$. Namely,
$\inf_{\xi\in\mathbb{R}^{n}}\left\lVert\nabla_{\times}\phi(\xi)\right\rVert>0$.
Due to Lemma 6, we have
$\lambda_{1}\dots\lambda_{n-1}=\det\big{(}Q(\xi)\big{)}=(k_{1}\dots
k_{n-1})^{2}\left\lVert\nabla_{\times}\phi\right\rVert^{2}$, where
$\lambda_{i}$ is an eigenvalue of the matrix $Q(\xi)$ defined in (31). The
previous reasoning and the boundedness of $|f_{i}^{\prime}(w)|$ imply that
$\Lambda^{\prime}:=\inf_{\xi\in\mathcal{E}_{\alpha}}\lambda_{\mathrm{min}}\big{(}Q(\xi)\big{)}>0$.
Using the function $V$ as in (29), we have
$\dot{V}\big{(}e(\xi)\big{)}\leq-\Lambda^{\prime}\left\lVert
e(\xi)\right\rVert^{2}\leq
2\Lambda^{\prime}V\big{(}e(\xi)\big{)}/k_{\mathrm{max}}$, further leading to
$||e(\xi)||\leq
c||e_{0}||\exp\left({-\Lambda^{\prime}t/k_{\mathrm{max}}}\right),$ where
$e_{0}=e\big{(}\xi(0)\big{)}$ and
$c=\sqrt{k_{\mathrm{max}}/k_{\mathrm{min}}}$. Therefore, the error
$\left\lVert e(\xi)\right\rVert$ will locally exponentially converges to $0$
as $t$ approaches infinity. If $|f_{i}^{\prime}(w)|$ is upper bounded globally
in $\mathbb{R}^{n}$, then the solution to
$\dot{\xi}(t)=\raisebox{2.0pt}{$\chi$}\big{(}\xi(t)\big{)}$ can be infinitely
prolonged for every initial condition $\xi(0)\in\mathbb{R}^{n}$, and thus the
norm of the path-following error $\left\lVert e(t)\right\rVert$ globally
exponentially converges to $0$ as $t\to\infty$ . ∎

Due to Proposition 4, the dynamical system with the guiding vector field
$\chi$ on a higher-dimensional space is robust against disturbance. To be more
specific, suppose disturbance denoted by $d$ is added to the vector field;
hence the perturbed system below is considered.

$\dot{\xi}(t)=\raisebox{2.0pt}{$\chi$}\big{(}\xi(t)\big{)}+d(t),$ (32)

where $d:\mathbb{R}_{\geq 0}\to\mathbb{R}^{n}$, is a piecewise continuous and
bounded function of time $t$ for all $t\geq 0$. It can be shown that the path-
following error dynamics

$\dot{e}(t)=N^{\top}\big{(}\xi(t)\big{)}\big{(}\raisebox{2.0pt}{$\chi$}(\xi(t))+d(t)\big{)}$

is locally input-to-state stable. This infers the following property of the
path-following error.

######  Proposition 5.

Consider the perturbed system (32). If $|f_{i}^{\prime}(w)|$ is upper bounded
in $\mathcal{E}_{\alpha}$ for $i=1,\dots,n-1$, then the error dynamics is
locally input-to-state stable. Namely, there exists positive constants
$\delta$ and $r$, such that if the initial condition
$\xi(0)\in\mathcal{E}_{\delta}$ and the disturbance is uniformly bounded by
$r$ (i.e., $\left\lVert d(t)\right\rVert\leq r$ for $t\geq 0$), then the norm
of the path-following error $\left\lVert e(t)\right\rVert$ is uniformly
ultimately bounded. In addition, if the disturbance is vanishing; that is,
$\lim_{t\to\infty}\left\lVert d(t)\right\rVert=0$, then the norm of the path-
following error $\left\lVert e(t)\right\rVert$ is also vanishing as
$t\to\infty$.

###### Proof.

The proof is similar to [13, Theorem 3], except that we consider a guiding
vector field on $\mathbb{R}^{n}$. ∎

## References

  * [1] B. Siciliano and O. Khatib, _Springer Handbook of Robotics_. Secaucus, NJ, USA: Springer-Verlag New York, Inc., 2007.
  * [2] S. Lacroix, G. Roberts, E. Benard, M. Bronz, F. Burnet, E. Bouhoubeiny, J.-P. Condomines, C. Doll, G. Hattenberger, F. Lamraoui _et al._ , “Fleets of enduring drones to probe atmospheric phenomena with clouds,” in _EGU General Assembly Conference Abstracts_ , vol. 18, 2016.
  * [3] V. M. Goncalves, L. C. A. Pimenta, C. A. Maia, B. C. O. Dutra, and G. A. S. Pereira, “Vector fields for robot navigation along time-varying curves in $n$-dimensions,” _IEEE Transactions on Robotics_ , vol. 26, no. 4, pp. 647–659, Aug 2010.
  * [4] A. M. Rezende, V. M. Gonçalves, G. V. Raffo, and L. C. Pimenta, “Robust fixed-wing UAV guidance with circulating artificial vector fields,” in _2018 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS)_. IEEE, 2018, pp. 5892–5899.
  * [5] M. M. Michałek and T. Gawron, “VFO path following control with guarantees of positionally constrained transients for unicycle-like robots with constrained control input,” _Journal of Intelligent & Robotic Systems_, vol. 89, no. 1-2, pp. 191–210, 2018.
  * [6] H. G. de Marina, Z. Sun, M. Bronz, and G. Hattenberger, “Circular formation control of fixed-wing UAVs with constant speeds,” in _2017 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS)_ , Sept 2017, pp. 5298–5303.
  * [7] Y. Liang and Y. Jia, “Combined vector field approach for 2D and 3D arbitrary twice differentiable curved path following with constrained UAVs,” _Journal of Intelligent & Robotic Systems_, vol. 83, no. 1, pp. 133–160, 2016.
  * [8] Y. A. Kapitanyuk, S. A. Chepinskiy, and A. A. Kapitonov, “Geometric path following control of a rigid body based on the stabilization of sets,” _IFAC Proceedings Volumes_ , vol. 47, no. 3, pp. 7342 – 7347, 2014, 19th IFAC World Congress.
  * [9] K. Łakomy and M. M. Michałek, “The VFO path-following kinematic controller for robotic vehicles moving in a 3D space,” in _Robot Motion and Control (RoMoCo), 2017 11th International Workshop on_. IEEE, 2017, pp. 263–268.
  * [10] M. Michalek and K. Kozlowski, “Vector-field-orientation feedback control method for a differentially driven vehicle,” _IEEE Transactions on Control Systems Technology_ , vol. 18, no. 1, pp. 45–65, 2010.
  * [11] W. Caharija, K. Y. Pettersen, M. Bibuli, P. Calado, E. Zereik, J. Braga, J. T. Gravdahl, A. J. Sørensen, M. Milovanović, and G. Bruzzone, “Integral line-of-sight guidance and control of underactuated marine vehicles: Theory, simulations, and experiments,” _IEEE Transactions on Control Systems Technology_ , vol. 24, no. 5, pp. 1623–1642, 2016.
  * [12] E. Borhaug, A. Pavlov, and K. Y. Pettersen, “Integral LOS control for path following of underactuated marine surface vessels in the presence of constant ocean currents,” in _2008 47th IEEE Conference on Decision and Control_. IEEE, 2008, pp. 4984–4991.
  * [13] W. Yao and M. Cao, “Path following control in 3D using a vector field,” _Automatica_ , vol. 117, p. 108957, 2020.
  * [14] A. P. Aguiar, J. P. Hespanha, and P. V. Kokotović, “Performance limitations in reference tracking and path following for nonlinear systems,” _Automatica_ , vol. 44, no. 3, pp. 598–610, 2008.
  * [15] M. M. Seron, J. H. Braslavsky, P. V. Kokotovic, and D. Q. Mayne, “Feedback limitations in nonlinear systems: From bode integrals to cheap control,” _IEEE Transactions on Automatic Control_ , vol. 44, no. 4, pp. 829–833, 1999\.
  * [16] Y. A. Kapitanyuk, A. V. Proskurnikov, and M. Cao, “A guiding vector-field algorithm for path-following control of nonholonomic mobile robots,” _IEEE Transactions on Control Systems Technology_ , vol. 26, no. 4, pp. 1372–1385, 2017.
  * [17] D. A. Lawrence, E. W. Frew, and W. J. Pisano, “Lyapunov vector fields for autonomous unmanned aircraft flight control,” _Journal of Guidance, Control, and Dynamics_ , vol. 31, no. 5, pp. 1220–1229, 2008.
  * [18] P. B. Sujit, S. Saripalli, and J. B. Sousa, “Unmanned aerial vehicle path following: A survey and analysis of algorithms for fixed-wing unmanned aerial vehicless,” _IEEE Control Systems_ , vol. 34, no. 1, pp. 42–59, Feb 2014\.
  * [19] W. F. Phillips, _Mechanics of flight_. John Wiley & Sons, 2004.
  * [20] B. Rubí, R. Pérez, and B. Morcego, “A survey of path following control strategies for UAVs focused on quadrotors,” _Journal of Intelligent & Robotic Systems_, pp. 1–25, 2019.
  * [21] W. Caharija, K. Y. Pettersen, P. Calado, and J. Braga, “A comparison between the ILOS guidance and the vector field guidance,” _IFAC-PapersOnLine_ , vol. 48, no. 16, pp. 89–94, 2015.
  * [22] W. Yao, Y. A. Kapitanyuk, and M. Cao, “Robotic path following in 3D using a guiding vector field,” in _IEEE Conference on Decision and Control_ , 2018, pp. 4475–4480.
  * [23] H. G. De Marina, Y. A. Kapitanyuk, M. Bronz, G. Hattenberger, and M. Cao, “Guidance algorithm for smooth trajectory tracking of a fixed wing UAV flying in wind flows,” in _2017 IEEE International Conference on Robotics and Automation (ICRA)_. IEEE, 2017, pp. 5740–5745.