# Distributed coordinated path following using guiding vector fields

Weijia Yao, Hector Garcia de Marina, Zhiyong Sun, Ming Cao Weijia Yao and Ming
Cao are with the University of Groningen, the Netherlands.
<EMAIL_ADDRESS>Hector Garcia de Marina is with Universidad Complutense
de Madrid, Spain<EMAIL_ADDRESS>Zhiyong Sun is with the Eindhoven University
of Technology, the Netherlands<EMAIL_ADDRESS>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). The work of Yao is
supported in part by the China Scholarship Council. The work of Hector Garcia
de Marina is supported by the grant Atraccion de Talento 2019-T2/TIC-13503
from the Government of the Autonomous Community of Madrid.

###### Abstract

It is essential in many applications to impose a scalable coordinated motion
control on a large group of mobile robots, which is efficient in tasks
requiring repetitive execution, such as environmental monitoring. In this
paper, we design a guiding vector field to guide multiple robots to follow
possibly different desired paths while coordinating their motions. The vector
field uses a path parameter as a virtual coordinate that is communicated among
neighboring robots. Then, the virtual coordinate is utilized to control the
relative parametric displacement between robots along the paths. This enables
us to design a saturated control algorithm for a Dubins-car-like model. The
algorithm is distributed, scalable, and applicable for any smooth paths in an
$n$-dimensional configuration space, and global convergence is guaranteed.
Simulations with up to fifty robots and outdoor experiments with fixed-wing
aircraft validate the theoretical results.

##  I Introduction

One of the grand challenges in robotics is the reliable and systematic
coordination of a large number of robots without a central nexus [1]. Tasks
involving vast geographic areas or volumes, such as search and rescue,
patrolling, and environmental monitoring, are practically impossible for a
single robot, while multiple coordinated robots can deal with such scenarios
efficiently. In these applications, it is fundamental for a group of robots to
accurately follow possibly different paths, while coordinating their motions
to form a formation which satisfies some geometric or parametric constraints.
In this paper, we design a _coordinating vector field_ to guide multiple
robots to converge to their respective desired paths and realize the
aforementioned coordinated behavior in a distributed way through local
information exchange.

Related work: Many path-following algorithms for a single robot have been
studied in the literature, such as Line-of-sight (LOS) [2] among others [3].
Among several tested algorithms, those using a guiding vector field are shown
to achieve the highest path-following accuracy while requiring the least
control efforts [3]. However, to the best of our knowledge, there are only a
few studies on extending the guiding vector field approach to realize multi-
robot coordinated path following tasks. Utilizing the guiding vector field in
[4], the work [5] employs a distributed algorithm to change the radii of the
circular paths for each fixed-wing aircraft flying at a constant speed such
that the aircraft asymptotically follow a circular path of the same radius
while keeping user-specified distances. Another work [6] designs different
potential functions to derive a guiding vector field for multiple robots, but
their motions cannot be explicitly coordinated (e.g., the multi-robot distance
cannot be controlled). These two studies only consider some common desired
paths, such as a circle or a straight line in 2D. The paper [7] provides a
decentralized control law such that multiple robots can circulate a closed
curve in 3D, but the curve must have a specific form, and some strong
assumptions on the integral curves are required.

There are many studies on coordinated path-following _without_ employing a
guiding vector field. The robots’ motions in [8] are coordinated by using a
virtual structure when following predefined paths. However, each robot needs
to broadcast its states and reference trajectories to the rest of the robots,
so the communication burden increases rapidly as the size of the multi-robot
systems increases. Recently, [9] proposes a distributed hybrid control law
based on [10] to coordinate the robot motions such that robots keep a constant
parametric separation. A similar task is achieved by another recent work in
[11]. However, these two studies, as well as many others (e.g., [12, 13, 14]),
are restricted to planar or simple desired paths, such as a straight line.
Another related work [15] only considers planar simple closed curves and unit-
speed particles, and the convergence result is valid locally. In [16], an
output-regulation controller is applied to certain robots while other robots
are controlled through local interactions. This approach is effective, but the
desired motion must be periodic and generated by a linear exosystem.

Multi-robot coordinated path following control has been extensively studied in
the specialized cases of circular formation control and circumnavigation
control. In circular formation control, multiple robots either need to move to
or initially start on a specified circle, and then move along the circle until
some desired arc distances are achieved via local communication [17, 18, 19,
20, 21, 22, 23]. Circumnavigation control additionally requires all robots to
keep moving along the circle to fulfill some tasks, such as escorting an
important target [24, 25, 26, 27, 28, 29, 30, 31, 32]. However, these studies
are all restricted to 2D circular paths, and many of them impose additional
requirements on the robots’ initial positions (e.g., to maintain the order
preservation [22, 21]). The paper [33] presents a decentralized controller for
robots characterized by the double-integrator model such that they can
generate 2D desired geometric patterns composed of simple closed curves. This
work is extended in [34] for robots described by the single-integrator model
to circulate a closed curve. However, these two studies are both restricted to
2D closed curves.

Contributions: We propose a coordinating vector field to enable distributed
coordinated path following for multiple robots, and a corresponding saturated
control algorithm is designed for a Dubins-car-like model. In particular, we
use the path parameter as an additional virtual coordinate for the guiding
vector field, and this virtual coordinate is further communicated among
neighboring robots to achieve motion coordination. There are many appealing
features of our approach: 1) It is distributed and scalable, and thus
applicable for a team of an arbitrary number of robots; 2) It can deal with
complex paths such as self-intersecting, non-convex, or non-closed ones in an
$n$-dimensional configuration space. 3) Given a fixed communication frequency,
the communication burden is low since every two neighboring robots only need
to transmit and receive at most two scalars (i.e., the virtual coordinate and
its time derivative); 4) It guarantees that the coordinated path-following
control objectives are achieved globally regardless of robots’ initial
positions. 5) The approach has been verified by simulations with a large
multi-robot system and by experiments with fixed-wing aircraft in a windy
environment. We also show that our approach is promising in area or volume
coverage tasks by exploiting a 2D or 3D Lissajous curve with _irrational
coefficients_ , as it is guaranteed that robots can eventually visit every
point of an area or volume while keeping prescribed distances.

Notations: The notation $\mathbb{Z}_{i}^{j}$ represents the set
$\\{m\in\mathbb{Z}:i\leq m\leq j\\}$. We use boldface for a vector
$\bm{v}\in\mathbb{R}^{n}$, and its $i$-th entry is denoted by $v_{i}$ for
$i\in\mathbb{Z}_{1}^{n}$. Consider a system consisting of $N$ robots, and a
vector $\bm{u_{i}}\in\mathbb{R}^{n}$ for $i\in\mathbb{Z}_{1}^{N}$, then for
$j\in\mathbb{Z}_{1}^{n}$, the notation $u_{ij}$ refers to the $j$-th entry of
$\bm{u_{i}}$. Graphs: We call $\mathcal{V}:=\\{1,\dots,N\\}$ the _node set_
representing the robots, and
$\mathcal{E}\subseteq(\mathcal{V}\times\mathcal{V})$ is the _edge set_
encoding the communication links. We call
$\mathcal{N}_{i}:=\\{j\in\mathcal{V}:(i,j)\in\mathcal{E}\\}$ the set of
neighbors of robot $i$. In this paper, we only consider undirected graphs,
which implies that if $(i,j)\in\mathcal{E}$, then Robot $i$ and Robot $j$ can
share information bidirectionally. Please see [35] for an introduction to
graph theory.

##  II Preliminaries

Guiding vector fields with virtual coordinates for robot navigation enjoy many
desirable features, such as guaranteeing global convergence to the desired
path and enabling self-intersecting desired path following [36]. We present
here a brief introduction. Suppose the desired path
${\mathcal{P}}^{\mathrm{phy}}$ is parameterized by
$x_{1}=f_{1}(w),\dots,x_{n}=f_{n}(w),$ where $x_{i}$ is the $i$-th coordinate,
$w\in\mathbb{R}$ is the parameter of the path and $f_{i}$ is twice
continuously differentiable (i.e., $f_{i}\in C^{2}$), for
$i\in\mathbb{Z}_{1}^{n}$. To derive a corresponding guiding vector field, one
needs to describe the desired path as the intersection of $n$ hyper-surfaces
[37, 38, 39, 40, 41]. To this end, taking $w$ as an additional argument, one
defines $n$ _surface functions_ $\phi_{i}:\mathbb{R}^{n+1}\to\mathbb{R}$ as
follows: $\phi_{i}(\bm{\xi})=x_{i}-f_{i}(w)$ for $i\in\mathbb{Z}_{1}^{n}$,
where $\bm{\xi}=(x_{1},\dots,x_{n},w)\in\mathbb{R}^{n+1}$ is the generalized
coordinate with an additional (virtual) coordinate $w$. Therefore, the desired
path with an additional coordinate is the intersection of the hyper-surfaces
described by the zero-level set of these functions; i.e.,
${\mathcal{P}}^{\mathrm{hgh}}:=\\{\bm{\xi}\in\mathbb{R}^{n+1}:\phi_{i}(\bm{\xi})=0,i\in\mathbb{Z}_{1}^{n}\\}$.
The projection of ${\mathcal{P}}^{\mathrm{hgh}}$ on the plane spanned by the
first $n$ coordinates is the original desired path
${\mathcal{P}}^{\mathrm{phy}}$, so we can use the higher-dimensional guiding
vector field $\raisebox{2.0pt}{$\chi$}:\mathbb{R}^{n+1}\to\mathbb{R}^{n+1}$
corresponding to ${\mathcal{P}}^{\mathrm{hgh}}$ to make a robot follow the
original desired path ${\mathcal{P}}^{\mathrm{phy}}$ by projecting to the
$n$-dimensional space [36], where

$\displaystyle\raisebox{2.0pt}{$\chi$}(\bm{\xi})=\times\big{(}\nabla\phi_{i},\dots,\nabla\phi_{n}\big{)}-\sum\nolimits_{i=1}^{n}k_{i}\phi_{i}\nabla\phi_{i},$

(1)

with $k_{i}>0$, $\nabla\phi_{i}\in\mathbb{R}^{n+1}$ being the gradient of
$\phi_{i}$ with respect to the generalized coordinate $\bm{\xi}$, and
$\times(\cdot)$ being the wedge product [42, Chapter 7.2]. The physical
interpretation of the vector field $\chi$ is clear: the second term
$-\sum_{i=1}^{n}k_{i}\phi_{i}\nabla\phi_{i}$ is a weighted sum of all the
gradients, which guides the robot trajectory towards the intersection of the
hyper-surfaces (i.e., the desired path), while the first term
$\times\big{(}\nabla\phi_{i},\dots,\nabla\phi_{n}\big{)}$, being orthogonal to
all the gradients [42, Proposition 7.2.1], provides a propagation direction
along the desired path.

In our previous work [36, 43], this higher-dimensional guiding vector field in
(1) is proved to possess no singular points, where the vector field becomes
zero, thanks to the additional dimension (i.e., the additional virtual
coordinate $w$). However, like many other current studies (e.g., [38, 39,
44]), the vector field is used for the guidance of only one single robot. We
aim to extend the vector field to include a coordination component and achieve
motion coordination among multiple robots. Since the additional virtual
coordinate $w$ not only helps eliminate singular points, but also acts as the
path parameter, one idea is to utilize the virtual coordinate $w$ to
coordinate the robots’ motions. Thus, the question is:

Question: Suppose there are $N$ robots. Based on the higher-dimensional vector
field in (1), how can one design an extra _coordination mechanism_ using the
virtual coordinates $w_{i}$, $i\in\mathbb{Z}_{1}^{N}$, such that 1) Each robot
can follow their desired paths; 2) All robots can coordinate their motions by
controlling the virtual coordinates $w_{i},i\in\mathbb{Z}_{1}^{N}$ via local
information exchange among neighbors?

The extra coordination mechanism and the precise meaning of motion
coordination will become clear in the next section, where a mathematical
problem formulation based on dynamical systems theory will be presented.

##  III Design of guiding vector fields for multi-robot systems

###  III-A Mathematical problem formulation

The $i$-th robot is required to follow a path in $\mathbb{R}^{n}$,
parameterized by $n$ parametric equations:

$x_{i1}=f_{i1}(w_{i})\quad\dots\quad x_{in}=f_{in}(w_{i}),$ (2)

where $x_{ij}$ is the $j$-th coordinate, $f_{ij}\in C^{2}$ is the $j$-th
parametric function for the $i$-th robot, for
$i\in\mathbb{Z}_{1}^{N},j\in\mathbb{Z}_{1}^{n}$, and $w_{i}\in\mathbb{R}$ is
the parameter of the desired path. To derive the path-following guiding vector
field, we use the parameter $w_{i}$ as an additional virtual coordinate, and
the higher-dimensional desired path is described by

$\mathcal{P}_{i}:=\\{\bm{\xi_{i}}\in\mathbb{R}^{n+1}:\phi_{i1}(\bm{\xi_{i}})=0,\dots,\phi_{in}(\bm{\xi_{i}})=0\\},$

where $\bm{\xi_{i}}:=(x_{i1},\dots,x_{in},w_{i})\in\mathbb{R}^{n+1}$ denotes
the _generalized coordinate_ of the $i$-th robot, and the surface functions
are $\phi_{ij}=x_{ij}-f_{ij}(w_{i})$ for $i\in\mathbb{Z}_{1}^{N}$ and
$j\in\mathbb{Z}_{1}^{n}$. Define
$\bm{\Phi_{i}}(\bm{\xi_{i}}):=(\phi_{i1}(\bm{\xi_{i}}),\dots,\phi_{in}(\bm{\xi_{i}}))^{\top}\in\mathbb{R}^{n}$.
Observe that $\bm{\xi}_{i}\in\mathcal{P}_{i}$ if and only if
$\left\lVert\bm{\Phi_{i}}(\bm{\xi}_{i})\right\rVert=0$. Therefore, one can use
$\bm{\Phi_{i}}(\bm{\xi_{i}})$ to quantify the distance to the desired path
$\mathcal{P}_{i}$; we call $\bm{\Phi_{i}}(\bm{\xi_{i}})$ the _path-following
error_ to $\mathcal{P}_{i}$. The aim is to design guiding controllers such
that the norm $\left\lVert\bm{\Phi_{i}}(\bm{\xi_{i}})\right\rVert$ converge to
zero eventually. By combining (1) and (2) , we obtain the analytic expression
of the _path-following guiding vector field_ $\bm{\prescript{\rm
pf}{}{\raisebox{2.0pt}{$\chi$}}_{i}}:\mathbb{R}^{n+1}\to\mathbb{R}^{n+1}$ for
the $i$-th robot:

(3)

for $i\in\mathbb{Z}_{1}^{N}$, where $k_{ij}>0$ are constant gains, and
$f^{\prime}_{ij}$ are the derivatives of $f_{ij}$ with respect to the argument
$w_{i}$.

To achieve coordination in $w_{i}$, thus indirectly coordinate the positions
of robots, we introduce the _coordination component_ $\bm{\prescript{\rm
co}{}{\raisebox{2.0pt}{$\chi$}}_{i}}:\mathbb{R}^{N}\to\mathbb{R}^{n+1}$ for
the $i$-th robot, $i\in\mathbb{Z}_{1}^{N}$, which is defined to be

$\displaystyle\bm{\prescript{\rm
co}{}{\raisebox{2.0pt}{$\chi$}}_{i}}(\bm{w})=\big{(}0,\cdots,0,c_{i}(\bm{w})\big{)}^{\top},$
(4)

where $\bm{w}=(w_{1},\dots,w_{N})^{\top}$ and
$c_{i}:\mathbb{R}^{N}\to\mathbb{R}$ is called the _coordination function_ to
be designed later to enable coordination among robots through local
interactions via the neighboring virtual coordinates $w_{j}$ for
$j\in\mathcal{N}_{i}$. Specifically, the virtual coordinates of multiple
robots $w_{i}(t)-w_{j}(t)$ are to be steered to converge to the _desired
relative states_ $\Delta_{ij}\in\mathbb{R}$ for $(i,j)\in\mathcal{E}$. We will
see how to design $\Delta_{ij}$ in the beginning of Section III-B.

We design the $i$-th _coordinating vector field_
$\bm{\mathfrak{X}_{i}}:\mathbb{R}^{n+N}\to\mathbb{R}^{n+1}$ to be the weighted
sum of the path-following vector field $\bm{\prescript{\rm
pf}{}{\raisebox{2.0pt}{$\chi$}}_{i}}$ and the coordination component
$\bm{\prescript{\rm co}{}{\raisebox{2.0pt}{$\chi$}}_{i}}$ as below:

$\displaystyle\bm{\mathfrak{X}_{i}}(\bm{\xi_{i}},\bm{w})=\bm{\prescript{\rm
pf}{}{\raisebox{2.0pt}{$\chi$}}_{i}}(\bm{\xi_{i}})+k_{c}\bm{\prescript{\rm
co}{}{\raisebox{2.0pt}{$\chi$}}_{i}}(\bm{w}),$ (5)

where $k_{c}>0$ is a parameter to adjust the weighted contribution of
$\bm{\prescript{\rm pf}{}{\raisebox{2.0pt}{$\chi$}}_{i}}$ and
$\bm{\prescript{\rm co}{}{\raisebox{2.0pt}{$\chi$}}_{i}}$ to
$\bm{\mathfrak{X}_{i}}$. With a larger value of $k_{c}$, the motion
coordination is achieved faster. The coordinating vector field
$\bm{\mathfrak{X}_{i}}$ represents the desired moving direction for a robot,
so the robot is guided. Thus it is key to study the guidance result, or
precisely, the convergence results of the integral curves of the vector field
$\bm{\mathfrak{X}_{i}}$ for $i\in\mathbb{Z}_{1}^{N}$. Precisely, we define
$\bm{\xi}:=(\bm{\xi_{i}}^{\top},\dots,\bm{\xi_{N}}^{\top})^{\top}\in\mathbb{R}^{(n+1)N}$
and
$\bm{\mathfrak{X}}(\bm{\xi}):=(\bm{\mathfrak{X}_{1}},\dots,\bm{\mathfrak{X}_{N}})^{\top}\in\mathbb{R}^{(n+1)N}$,
and study the integral curves of $\bm{\mathfrak{X}}(\bm{\xi})$; i.e., the
trajectories to the differential equation

$\dot{\bm{\xi}}=\bm{\mathfrak{X}}(\bm{\xi}),$ (6)

given an initial condition $\bm{\xi_{0}}\in\mathbb{R}^{(n+1)N}$ at $t=0$. Now
we can formally formulate the problem as follows:

######  Problem 1.

Design the coordinating vector field $\bm{\mathfrak{X}_{i}}$ in (5) for
$i\in\mathbb{Z}_{1}^{N}$, such that the trajectories of (6), given an initial
condition $\bm{\xi_{0}}\in\mathbb{R}^{(n+1)N}$, fulfill the two control
objectives:

  1. 1.

(Path Following) Robot $i$’s path following error to its desired path
$\mathcal{P}_{i}$ converges to zero asymptotically for all
$i\in\mathbb{Z}_{1}^{N}$. That is,
$\left\lVert\bm{\Phi_{i}}(\bm{\xi_{i}}(t)\right\rVert\to 0$ as $t\to\infty$,
$\forall i\in\mathbb{Z}_{1}^{N}$.

  2. 2.

(Motion Coordination) Each robot’s motion is coordinated distributedly subject
to the communication graph $\mathcal{G}=(\mathcal{V},\mathcal{E})$ such that
their additional virtual coordinates satisfy $w_{i}(t)-w_{j}(t)-\Delta_{ij}\to
0$ as $t\to\infty$ for $(i,j)\in\mathcal{E}$.

We require two mild standing assumptions that are reasonable in practical
applications:

######  Assumption 1.

The communication graph $\mathcal{G}=(\mathcal{V},\mathcal{E})$ is connected.

######  Assumption 2.

The first and second derivatives of $f_{ij}(\cdot)$ are bounded for
$i\in\mathbb{Z}_{1}^{N},j\in\mathbb{Z}_{1}^{n}$.

###  III-B Design of the coordination component

Given the desired path $\mathcal{P}_{i}$, we can design the desired relative
states $\Delta_{ij}$ starting from a particular reference configuration
$\bm{w}^{*}:=(w^{*}_{1},\cdots,w^{*}_{N})^{\top}$. Hence,
$\bm{\Delta}^{*}=D^{\top}\bm{w}^{*}$ is the stacked vector of
$\Delta_{ij},(i,j)\in\mathcal{E}$ [45], where
$D\in\mathbb{R}^{N\times|\mathcal{E}|}$ is the incidence matrix [35, p. 23].
Now we propose to employ the following consensus control algorithm [46, p.
25]:

$c_{i}=-\sum_{j\in\mathcal{N}_{i}}\big{(}w_{i}-w_{j}-\Delta_{ij}\big{)},\;\forall
i\in\mathbb{Z}_{1}^{N}.$ (7)

Equation (7) can be written in a compact form as

$\bm{c}(\bm{w})=-L(\bm{w}-\bm{w}^{*})=-L\bm{\tilde{w}},$ (8)

where $\bm{c}(\bm{w})=(c_{1}(\bm{w}),\dots,c_{N}(\bm{w}))^{\top}$,
$L=L(\mathcal{G})$ is the Laplacian matrix and
$\bm{\tilde{w}}:=\bm{w}-\bm{w}^{*}$. Combining (3), (4), (5) and (7), we have
attained the coordinating vector field $\bm{\mathfrak{X}_{i}}$ for
$i\in\mathbb{Z}_{1}^{N}$.

######  Remark 1.

From (4), (5) and (7), one observes that neighboring information exchange only
happens in the coordination component $c_{i}$. Notably, the communication
burden is low: Robot $i$ transmits only a scalar $w_{i}$ to the neighboring
Robot $j\in\mathcal{N}_{i}$. ∎

###  III-C Convergence analysis

The convergence analysis of trajectories to (6) is nontrivial given that the
right-hand side of (6) is _not_ a gradient of any potential function, since
the path-following vector field in (1) contains a wedge product of all the
gradients. In this subsection, we show that the coordinating vector field (5)
enables multiple robots to follow their desired paths while being coordinated
by the virtual coordinates such that $w_{i}(t)-w_{j}(t)$ converge to
$\Delta_{ij}$ for $(i,j)\in\mathcal{E}$ as $t\to\infty$.

For simplicity, we first consider Robot $i$, and most of the function
arguments are ignored henceforth unless ambiguity arises. Define $K_{i}:={\rm
diag}\\{k_{i1},\dots,k_{in}\\}$ and
$\bm{f^{\prime}_{i}}(w_{i}):=(f_{i1}^{\prime}(w_{i}),\dots,f_{in}^{\prime}(w_{i}))^{\top}$.
Then one can calculate that
$\displaystyle(\nabla\phi_{i1}^{\top}\bm{\prescript{\rm
pf}{}{\raisebox{2.0pt}{$\chi$}}_{i}},\cdots,\nabla\phi_{in}^{\top}\bm{\prescript{\rm
pf}{}{\raisebox{2.0pt}{$\chi$}}_{i}})^{\top}=-K_{i}\bm{\Phi_{i}}-\bm{f^{\prime}_{i}}(w_{i})\bm{f^{\prime}_{i}}(w_{i})^{\top}K_{i}\bm{\Phi_{i}}$
and  $\displaystyle(\nabla\phi_{i1}^{\top}\bm{\prescript{\rm
co}{}{\raisebox{2.0pt}{$\chi$}}_{i}}(\bm{w}),\cdots,\nabla\phi_{in}^{\top}\bm{\prescript{\rm
co}{}{\raisebox{2.0pt}{$\chi$}}_{i}}(\bm{w}))^{\top}=-c_{i}(\bm{w})\bm{f^{\prime}_{i}}(w_{i})$
. Now we consider all robots. Define $\mathfrak{F}:={\rm
diag}\\{\bm{f^{\prime}_{1}},\dots,\bm{f^{\prime}_{N}}\\}\in\mathbb{R}^{nN\times
N}$, $K:={\rm diag}\\{K_{1},\dots,K_{N}\\}\in\mathbb{R}^{nN\times nN}$ and
$\bm{\Phi}:=(\bm{\Phi_{1}}^{\top},\cdots,\bm{\Phi_{N}}^{\top})^{\top}$. Then
we have

$\displaystyle\dot{\bm{\Phi}}=-K\bm{\Phi}-\mathfrak{F}\mathfrak{F}^{\top}K\bm{\Phi}-k_{c}\mathfrak{F}\bm{c(\bm{w})}.$
(9)

One can also calculate that
$\dot{w}_{i}=\begin{bmatrix}0&\cdots&0&1\end{bmatrix}(\bm{\prescript{\rm
pf}{}{\raisebox{2.0pt}{$\chi$}}_{i}}(\bm{\xi_{i}})+k_{c}\bm{\prescript{\rm
co}{}{\raisebox{2.0pt}{$\chi$}}_{i}}(\bm{\xi_{i}}))=(-1)^{n}+\bm{f^{\prime}_{i}}(w_{i})^{\top}K_{i}\bm{\Phi_{i}}(\bm{\xi_{i}})+k_{c}c_{i}(\bm{w})$.
Therefore, we have that

$\displaystyle\dot{\bm{\tilde{w}}}=(-1)^{n}\bm{1}+\mathfrak{F}^{\top}K\bm{\Phi}-k_{c}L\bm{\tilde{w}},$
(10)

where $\bm{1}\in\mathbb{R}^{N}$ is a vector consisting of all ones. The
Laplacian matrix can be re-written as $L=DD^{\top}$, where
$D\in\mathbb{R}^{N\times|\mathcal{E}|}$ is the incidence matrix. We define the
_composite error vector_
$\bm{e}:=(\bm{\Phi}^{\top},(D^{\top}\bm{\tilde{w}})^{\top})^{\top}\in\mathbb{R}^{nN+|\mathcal{E}|}$
and the composite gain matrix to be $\mathfrak{K}:={\rm
diag}\\{K,k_{c}I_{|\mathcal{E}|}\\}\in\mathbb{R}^{(nN+|\mathcal{E}|)\times(nN+|\mathcal{E}|)}$,
where $I_{m}$ is the $|\mathcal{E}|$-by-$|\mathcal{E}|$ identity matrix.
Therefore, from (8), (9) and (10), and noting that $D^{\top}\bm{1}=\bm{0}$, we
have the following composite error dynamics:

(11)

######  Theorem 1.

The coordinating vector fields $\bm{\mathfrak{X}_{i}}$ for
$i\in\mathbb{Z}_{1}^{N}$ designed by combining (3), (4), (5) and (7) solve
Problem 1 globally in the sense that the initial states
$\bm{\xi_{0}}\in\mathbb{R}^{(n+1)\times N}$ can be chosen arbitrarily.

###### Proof.

The proof is in the extended version [47]. ∎

######  Remark 2.

Our work’s focus is the distributed control algorithm based on guiding vector
fields, so collision avoidance is not elaborated here. However, our proposed
approach can incorporate existing collision avoidance algorithms [48, 49]. For
example, we can modify the _nominal_ guiding vector field (5) in a minimally
invasive way using _safety barrier certificates_ [50, 49]. Specifically, we
add an extra term $\bm{u_{i}^{\mathrm{col}}}\in\mathbb{R}^{n+1}$ to the vector
field in (5), and this term is calculated by a quadratic program
$\min\left\lVert\bm{u_{i}^{\mathrm{col}}}(t)\right\rVert^{2}$ subject to the
constraints $\dot{B}_{ij}(\bm{\xi_{i}},\bm{\xi_{j}})\leq
1/B_{ij}(\bm{\xi_{i}},\bm{\xi_{j}})$ for $j\neq i$ and
$i,j\in\mathbb{Z}_{1}^{N}$, where $B_{ij}(\bm{\xi_{i}},\bm{\xi_{j}})$ is a
control barrier function [51] (e.g.,
$B_{ij}(\bm{\xi_{i}},\bm{\xi_{j}})=1/(\left\lVert\bm{p_{i}}-\bm{p_{j}}\right\rVert^{2}-R^{2})$,
where $\bm{p}_{i}$ is the physical position and $R$ is the safe distance
between robots). The collision avoidance behavior is shown in the
supplementary video, but the theoretical analysis is left for future work. ∎

##  IV Controller for A Dubins-car-like model

If a robot’s motion can be approximately modeled by the single-integrator
model, then the coordinating vector field in (5) can be used directly as the
velocity input to the robot. For the unicycle model, one can use feedback
linearization to transform it into the single-integrator model [52] to utilize
the guiding vector field directly. However, we will design a controller for
unicycles traveling at a constant speed (i.e., the Dubins-car model) without
using the feedback linearization technique. Note that the control algorithm
design idea in this section is applicable to robot models of which the motions
are characterized by the robot’s orientations, such as the car-like model and
the underwater glider model [53]. These models (approximately) represent many
different robotic systems in reality; thus, the design methodology is widely
applicable.

Different from the unicycle model, which allows backward or stationary motion,
we use the following Dubins-car-like 3D model that describes fixed-wing
aircraft dynamics:

$\displaystyle\dot{p}_{i1}=v\cos\theta_{i},\;\dot{p}_{i2}=v\sin\theta_{i},\;\dot{p}_{i3}=u^{z}_{i},\;\dot{\theta_{i}}=u^{\theta}_{i},$
(12)

where $v$ is a constant airspeed,
$(p_{i1},p_{i2},p_{i3})^{\top}\in\mathbb{R}^{3}$ is the position of the $i$-th
aircraft’s center of mass, $\theta_{i}$ is the yaw angle, and $u^{z}_{i}$ and
$u^{\theta}_{i}$ are two control inputs to be designed. Since the essential
role of a guiding vector field is to provide the desired yaw angle to guide
the flight of a fixed-wing aircraft, the core idea behind the control
algorithm design is to align the aircraft’s flying direction with that given
by the guiding vector field. To this end, as the guiding vector field has an
additional coordinate, one needs to add an additional virtual coordinate
$p_{i4}$ such that the aircraft’s generalized position is
$(p_{i1},p_{i2},p_{i3},p_{i4})\in\mathbb{R}^{4}$. Correspondingly, its
generalized velocity is
$(\dot{p}_{i1},\dot{p}_{i2},\dot{p}_{i3},\dot{p}_{i4})=(v\cos\theta_{i},v\sin\theta_{i},u^{z}_{i},\dot{p}_{i4})$,
where $\dot{p}_{i4}$ is virtual and extra design freedom. To align the
aircraft heading $v(\cos\theta_{i},\sin\theta_{i})$ with the counterpart of
the coordinating vector field $\bm{\mathfrak{X}_{i}}$, we need to “partially
normalize” the vector field $\bm{\mathfrak{X}_{i}}$ such that its first two
entries form a vector of the same length as $v$; that is,
$v\bm{\underline{\mathfrak{X}_{i}}}:=v\bm{\mathfrak{X}_{i}}/\sqrt{\mathfrak{X}_{i1}^{2}+\mathfrak{X}_{i2}^{2}}.$
Subsequently, we design the yaw angular control input $u^{\theta}_{i}$ such
that the aircraft heading $v(\cos\theta_{i},\sin\theta_{i})$ aligns gradually
with the vector formed by the first two entries of
$v\bm{\underline{\mathfrak{X}_{i}}}$ (i.e.,
$v({\underline{\mathfrak{X}_{i}}}{}_{1},{\underline{\mathfrak{X}_{i}}}{}_{2})$).
For the last two entries of the generalized velocity
$(\dot{p}_{i1},\dot{p}_{i2},\dot{p}_{i3},\dot{p}_{i4})$, one can simply equate
them with those of $v\bm{\underline{\mathfrak{X}_{i}}}$, respectively:

(13)

The control algorithm design method [36, 54] is extended here to handle the
issue with the actuator saturation in the yaw control input $u^{\theta}_{i}$
as described below. First, we define the saturation function
$\mathrm{Sat}_{a}^{b}:\mathbb{R}\to\mathbb{R}$ by $\mathrm{Sat}_{a}^{b}(x)=x$
for $x\in[a,b]$, $\mathrm{Sat}_{a}^{b}(x)=a$ for $x\in(-\infty,a)$ and
$\mathrm{Sat}_{a}^{b}(x)=b$ for $x\in(b,\infty)$, where $a,b\in\mathbb{R},a<b$
are some constants. The saturation function $\mathrm{Sat}_{a}^{b}$ is
Lipschitz continuous. For convenience, we call the time interval when
$\mathrm{Sat}_{a}^{b}(x(t))=b$ the _upper saturation period_ , and the time
interval when $\mathrm{Sat}_{a}^{b}(x(t))=a$ the _lower saturation period_. We
use $\overline{\bm{v}}$ to denote the normalization of a vector $\bm{v}$
(i.e., $\overline{\bm{v}}=\bm{v}/\left\lVert\bm{v}\right\rVert$). We also
define
$\bm{\mathfrak{X}^{p}_{i}}=({\overline{\mathfrak{X}_{i}}}_{1},{\overline{\mathfrak{X}_{i}}}_{2})^{\top}$,
which is the vector formed by the first two entries of the normalized vector
field $\overline{\bm{\mathfrak{X}_{i}}}$. One can easily calculate that
$\overline{\bm{\mathfrak{X}_{i}}^{p}}=({\underline{\mathfrak{X}_{i}}}{}_{1},{\underline{\mathfrak{X}_{i}}}{}_{2})$.
Therefore, $\overline{\bm{\mathfrak{X}_{i}}^{p}}$ represents the orientation
given by the vector field $\bm{\mathfrak{X}_{i}}$.

Suppose we are given 3D physical desired paths
${\mathcal{P}}^{\mathrm{phy}}_{i}\subseteq\mathbb{R}^{3}$ parameterized by
(2), and the coordinating vector field
$\bm{\mathfrak{X}_{i}}:\mathbb{R}^{3+N}\to\mathbb{R}^{3+1}$ in (5). We denote
$\bm{h_{i}}=\overline{\bm{h_{i}}}=(\cos\theta_{i},\sin\theta_{i})^{\top}$ as
the orientation of the aircraft, the (signed) angle difference directed from
$\overline{\bm{\mathfrak{X}^{p}_{i}}}$ to $\overline{\bm{h_{i}}}$ by
$\sigma_{i}\in(-\pi,\pi]$ and define the rotation matrix
$E=\left[\begin{smallmatrix}0&-1\\\ 1&0\end{smallmatrix}\right]$. Then after
introducing another natural assumption, one can reach Theorem 2 showing that
the aircraft’s orientation $\overline{\bm{h_{i}}}$ will converge to that of
the vector field $\overline{\bm{\mathfrak{X}_{i}}^{p}}$ asymptotically (i.e.,
$\sigma_{i}$ converges to zero).

######  Assumption 3.

The desired paths can be traveled by robots with saturated angular velocity,
and the eventual coordination is feasible [55].

This assumption is reasonable. It restricts the design of the desired path and
the reference state vector $\bm{\omega}^{*}$ (hence $\bm{\Delta}^{*}$) to
guarantee the feasibility of the robot motions, and avoid the collision of
robots in the steady-state, respectively. For example, if the parametric
desired path is periodic and closed, then one can simply choose
$\Delta_{ij}=T/N$, where $T$ is the period, for any two neighboring indices
$i,j\in\mathbb{Z}_{1}^{N}$.

Figure 1: The first simulation results. The path is parameterized by
$x_{i1}=15\sin(2w_{i})$,
$x_{i2}=30\sin(w_{i})\sqrt{0.5(1-0.5\sin^{2}(w_{i}))}$ and
$x_{i3}=5+5\cos(2w_{i})-2$ for all $i\in\mathbb{Z}_{1}^{N}$. The period of
this closed curve is $T=2\pi$ and the desired differences between two adjacent
robots’ virtual coordinates are $T/(2N)$. The desired relative states
$\Delta_{ij}$ are constructed from the reference $w_{i}^{*}=(i-1)T/(2N)$ for
$i\in\mathbb{Z}_{1}^{N}$. The control gains are:
$k_{i1}=k_{i2}=k_{i3}=1,k_{c}=300$ for $i\in\mathbb{Z}_{1}^{N}$, where $k_{c}$
is large to accelerate the motion coordination. On top, the trajectories of
robots, where squares and circles symbolize the trajectories’ initial and
final positions, respectively. On the bottom left, the path-following errors
$e_{ij}$ for $i\in\mathbb{Z}_{1}^{N}$ and $j\in\mathbb{Z}_{1}^{n}$. On the
bottom right, the coordination error $w_{i}-w_{j}-\Delta_{ij}$ for
$i,j\in\mathbb{Z}_{1}^{N},i<j$.

######  Theorem 2.

Assume that the vector field satisfies
$\mathfrak{X}_{i1}(\bm{\xi_{i}})^{2}+\mathfrak{X}_{i2}(\bm{\xi_{i}})^{2}>\gamma>0$
for $\bm{\xi_{i}}\in\mathbb{R}^{3+1}$, $i\in\mathbb{Z}_{1}^{N}$, where
$\gamma$ is a positive constant. Let the angular velocity control input
$u^{\theta}_{i}$ in model (12) be
$\dot{\theta_{i}}=u^{\theta}_{i}=\mathrm{Sat}_{a}^{b}(\dot{\theta}^{d}_{i}-k_{\theta}\overline{\bm{h_{i}}}^{\top}E\overline{\bm{\mathfrak{X}^{p}_{i}}}),$
where
$\dot{\theta}^{d}_{i}=-\overline{\bm{\mathfrak{X}^{p}_{i}}}^{\top}E\bm{\dot{\mathfrak{X}}^{p}_{i}}/\left\lVert\bm{\mathfrak{X}^{p}_{i}}\right\rVert$,
$k_{\theta}>0$ is a constant, and $a<0$, $b>0$ are constants for the
saturation function $\mathrm{Sat}_{a}^{b}$. If the angle difference
$\sigma_{i}$ satisfies the following conditions:

  1. 1.

The initial angle difference $\sigma_{i}(t=0)\neq\pi$;

  2. 2.

$\sigma_{i}(t)\in[0,\pi)$ during the upper saturation period, and
$\sigma_{i}(t)\in(-\pi,0]$ during the lower saturation period,

then $\sigma_{i}$ will vanish asymptotically (i.e., $\sigma_{i}(t)\to 0$).

###### Proof.

The proof is in the extended version [47]. ∎

Although Condition 2 in Theorem 2 might be difficult to check in practice, _it
conveys the core message that the saturation, albeit allowed, should not last
for a long time_. However, this condition is only sufficient, while in
practice, violating this condition does not immediately entail instability of
the algorithm. As the aircraft is guided by the vector field, the aircraft can
re-orient its heading towards the desired path even if it temporarily deviates
from the desired path due to saturation or other constraints (such as path
curvature), as long as the _vector field orientation change rate_ 111As
$\dot{\overline{\bm{\mathfrak{X}^{p}_{i}}}}=\dot{\theta}^{d}_{i}E\overline{\bm{\mathfrak{X}^{p}_{i}}}$,
and $E\overline{\bm{\mathfrak{X}^{p}_{i}}}$ is orthogonal to
$\overline{\bm{\mathfrak{X}^{p}_{i}}}$, the term $\dot{\theta}^{d}_{i}$ in (2)
encodes how fast the vector field changes its orientation along the trajectory
of the aircraft (i.e., $\dot{\theta}^{d}_{i}$ is the _vector field orientation
change rate_). $\dot{\theta}^{d}_{i}$ does not saturate the control input
persistently. The subsequent fixed-wing aircraft experiment verifies the
effectiveness of the control law. Nevertheless, we can remove Condition 2 in
Theorem 2 by imposing an upper bound on the magnitude of
$\dot{\theta}^{d}_{i}$. Due to the page limit, the detail is in the extended
version [47]. To reduce $|\dot{\theta}^{d}_{i}|$ and avoid possible
saturation, one can scale down the path parameter in the parametric functions
in (2) (e.g. by changing $f_{ij}(w_{i})$ to $f_{ij}(\beta w_{i})$, where
$0<\beta<1$), or choose another desired path with a smaller curvature. Another
solution is to impose an additional constraint
$|u^{\theta}_{i}|\leq\min\\{|a|,|b|\\}$ on the quadratic program mentioned in
Remark 2.

######  Remark 3.

Theorem 2 guarantees that each robot’s moving direction will align with the
vector field eventually. Therefore, it can be shown that the convergence
results of the robots’ trajectories are the same as the integral curves of the
guiding vector field in (6). The rigorous theoretical analysis relies on the
input-to-state stability (ISS) property [56, Chapter 4.9], and is left for
future work. ∎

##  V Simulations and Experiments

###  V-A Simulations

We use _cycle graphs_ as the communication topology, and thus each robot is
only allowed to communicate with its two adjacent neighbors. This simple
communication graph results in very low communication and computation load. In
the first simulation, we let $N=50$ robots follow a 3D _self-intersecting_
bent “$\infty$”-shaped curve. As shown in Fig. 1, all robots follow the
“$\infty$”-shaped path successfully and keep desired positions (in terms of
$w_{i}$) between each other. In the second simulation, we employ $N=3$ robots
to show that our algorithm is applicable for complicated and _open_ curves
(i.e., aperiodic curve), and demonstrate its potential application to volume
coverage in 3D (see Fig. 2). In the third simulation, we show how different
robots can follow different desired paths while still coordinating their
motions to form some (time-varying) formation shapes (see Fig. 3).

Figure 2: The second simulation results. We choose a 3D Lissajous curve with
_irrational_ coefficients: $x_{i1}=\cos(n_{x}w_{i})+m_{x}$,
$x_{i2}=\cos(n_{y}w_{i})+m_{y}$ and $x_{i3}=\cos(n_{z}w_{i})+m_{z}$ for
$i\in\mathbb{Z}_{1}^{N}$, where
$n_{x}=\sqrt{2},n_{y}=4.1,n_{z}=7.1,m_{x}=0.1,m_{y}=0.7,m_{z}=0$. This is an
open curve bounded in a cube (as $n_{x}$ is irrational). Therefore, it is
ideal for a volume coverage task. The control gains for the coordinating
vector field are: $k_{i1}=k_{i2}=k_{i3}=k_{c}=1$ for $i\in\mathbb{Z}_{1}^{N}$.
Squares and circles symbolize trajectories’ initial and final positions, and
the solid blue lines are the trajectories during the last $30$ time steps. The
first sub-figure shows the trajectories of three robots (red curves). The
magenta curve represents part of the Lissajous curve. The last two sub-figures
correspond to the Y-Z side views of the trajectories for $9.4$ and $76.8$
seconds, respectively.  Figure 3: The third simulation results. We let $N=21$
robots follow three different paths, where the first seven robots follow a
large circle of radius $a=10$, the last seven robots follow a small circle of
radius $b=5$, and the remaining seven robots follow an ellipse with the
semimajor axis $a$ and the semiminor axis $b$. The robots are distributed
evenly in the path through the parameter $w_{i}$; i.e., we construct the
desired relative states $\Delta_{ij}$ from the references
$w_{i}^{*}(t)=(i-1)2\pi/N$, for $i\in\mathbb{Z}_{1}^{N}$. The control gains
are: $k_{i1}=k_{i2}=1,k_{c}=100$ for all $i\in\mathbb{Z}_{1}^{N}$, where the
large value is chosen for $k_{c}$ to accelerate the motion coordination. The
squares and $*$ symbolize the trajectories’ initial and final positions,
respectively. The dashed black line shows the communication links between
robots and the resulting formation shape. During the steady-state, these
robots do not collide into one another, since the distributed coordination
guarantees that the adjacent neighbors satisfy $|w_{i}-w_{j}|=2\pi/N$, where
$i$ and $j$ are any neighboring indices.

###  V-B Flights with fixed-wing aircraft

In this experiment, two autonomous fixed-wing aircraft (i.e., Autonomous
Opterra 1.2m) similar to [5, 57] are employed to validate Theorem 2. The
aircraft are equipped with the open-source software/hardware project
_Paparazzi_[58]. The codes related to the proposed algorithm are in [59]. We
choose the following 3D Lissajous curve
$f_{1}(w)=225\cos(w),f_{2}(w)=225\cos(2w+\pi/2),f_{3}(w)=-20\cos(2w)$, which
is a bent “$\infty$”-shaped path. The mission requires both aircraft to have
$\Delta_{12}=0$; i.e., to fly in a tight formation. The collision is avoided
by biasing the GPS measurement of one aircraft by a constant distance of $1m$
in the horizontal plane; i.e., when the aircraft achieve $\Delta_{12}=0$, they
are displaced physically. We choose $k_{1}=k_{2}=0.002,k_{3}=0.0025$,
$k_{c}=0.01,k_{\theta}=1$, and the communication frequency is $10$ Hz. In the
experiment, the weather forecast reported 14 degrees Celsius and a South wind
of 10 km/h. In Figure 4, the telemetry shows that both airplanes converge to
fly together while following the path. The experiment shows that once an
aircraft goes faster than its partner, the algorithm guides the airplane to
deviate from the curve to travel more distances to “wait” for its partner.
Nevertheless, these deviations from the desired path are within the order of
one or two meters (see Figure 5).

Note that the employed aircraft do _not_ control their ground speeds. In fact,
they have a reference signal in their throttle to keep a safe airspeed, and
the aircraft increase/decrease such a reference to ascend/descend while
following the 3D path. A traditional trajectory tracking algorithm would
_force_ the aircraft to track an _open-loop_ $f_{1}(t),f_{2}(t)$ and
$f_{3}(t)$; i.e., it requires controlling the airspeed/ground speed of the
aircraft. Such a requirement is demanding (or unsatisfiable) if the aircraft
are not equipped with the required sensors or actuators (e.g., spoilers/flaps
while descending), and the wind always affects the speed of the airplane. By
contrast, our guidance algorithm is free from such a requirement since the
parameter $w$ is in a _closed-loop_ with the aircraft state, and thereby
adapting automatically to the position and velocity of the aircraft.

Figure 4: Two aircraft (blue and red trajectories) flying together tracking a
3D bent “$\infty$”-shaped path. Although we show that the positions are
overlapped, the same position corresponds to two different ones in reality
since their GPS receptors are biased with respect to each other by around 1
meter in the XY plane. Figure 5: Distance errors (in XYZ) of the two aircraft
(blue and red) to the desired 3D path, and errors with respect to the desired
$\Delta w=0$. The horizontal axes denote time in seconds.

##  VI Conclusions and future work

We design a coordinating vector field for multiple robots to follow (possibly
different) desired paths while coordinating their motions distributedly, and a
control algorithm is designed for the Dubins-car-like model. In particular, we
use the path parameter $w_{i}$ as an additional virtual coordinate for the
guiding vector field, which is further utilized as local information
communicated among neighboring robots to achieve motion coordination. The
vector field is rigorously proved to enable robots to follow their respective
desired paths while coordinating their motions. Our algorithm is a significant
extension of the circumnavigation algorithms and circular formation control
algorithms since it enables multiple robots to cooperatively move along
arbitrary paths and form different (possibly varying) desired shapes.
Moreover, our algorithm requires a very low amount of local communication.
Simulations and experiments support the theoretical results. Future work
includes detailed analysis of the collision avoidance using control barrier
functions, and performance evaluation in the presence of localization errors.
We will compare our algorithm with the approaches mentioned in the
Introduction in our future work, and will also consider the problem on smooth
manifolds.

## References

  * [1] G.-Z. Yang, J. Bellingham, P. E. Dupont, P. Fischer, L. Floridi, R. Full, N. Jacobstein, V. Kumar, M. McNutt, R. Merrifield _et al._ , “The grand challenges of science robotics,” _Science robotics_ , vol. 3, no. 14, p. 7650, 2018.
  * [2] T. I. Fossen, M. Breivik, and R. Skjetne, “Line-of-sight path following of underactuated marine craft,” _IFAC Proceedings Volumes_ , vol. 36, no. 21, pp. 211–216, 2003.
  * [3] 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\.
  * [4] 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, July 2018.
  * [5] 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)_. IEEE, 2017, pp. 5298–5303.
  * [6] K. Nakai and K. Uchiyama, “Vector fields for UAV guidance using potential function method for formation flight,” in _AIAA Guidance, Navigation, and Control (GNC) Conference_ , 2013, p. 4626.
  * [7] L. C. Pimenta, G. A. Pereira, M. M. Gonçalves, N. Michael, M. Turpin, and V. Kumar, “Decentralized controllers for perimeter surveillance with teams of aerial robots,” _Advanced Robotics_ , vol. 27, no. 9, pp. 697–709, 2013\.
  * [8] J. Ghommam, H. Mehrjerdi, M. Saad, and F. Mnif, “Formation path following control of unicycle-type mobile robots,” _Robotics and Autonomous Systems_ , vol. 58, no. 5, pp. 727–736, 2010.
  * [9] H. Chen, Y. Cong, X. Wang, X. Xu, and L. Shen, “Coordinated path-following control of fixed-wing unmanned aerial vehicles,” _IEEE Transactions on Systems, Man, and Cybernetics: Systems_ , 2021.
  * [10] Y. Lan, G. Yan, and Z. Lin, “Synthesis of distributed control of coordinated path following based on hybrid approach,” _IEEE Transactions on Automatic Control_ , vol. 56, no. 5, pp. 1170–1175, 2011.
  * [11] Y. Wang, D. Wang, and S. Zhu, “Cooperative moving path following for multiple fixed-wing unmanned aerial vehicles with speed constraints,” _Automatica_ , vol. 100, pp. 82–89, 2019.
  * [12] M. Burger, A. Pavlov, E. Borhaug, and K. Y. Pettersen, “Straight line path following for formations of underactuated surface vessels under influence of constant ocean currents,” in _2009 American Control Conference_. IEEE, 2009, pp. 3065–3070.
  * [13] A. Doosthoseini and C. Nielsen, “Coordinated path following for unicycles: A nested invariant sets approach,” _Automatica_ , vol. 60, pp. 17–29, 2015\.
  * [14] L. A. V. Reyes and H. G. Tanner, “Flocking, formation control, and path following for a group of mobile robots,” _IEEE Transactions on Control Systems Technology_ , vol. 23, no. 4, pp. 1268–1282, 2015.
  * [15] F. Zhang and N. E. Leonard, “Coordinated patterns of unit speed particles on a closed curve,” _Systems & control letters_, vol. 56, no. 6, pp. 397–407, 2007.
  * [16] L. Sabattini, C. Secchi, M. Cocetti, A. Levratti, and C. Fantuzzi, “Implementation of coordinated complex dynamic behaviors in multirobot systems,” _IEEE Transactions on Robotics_ , vol. 31, no. 4, pp. 1018–1032, 2015.
  * [17] M. Iqbal, S.-i. Azuma, J. Leth, and T. D. Ngo, “Circular and concentric formation of kinematic unicycles,” in _2019 International Conference on Robotics and Automation (ICRA)_. IEEE, 2019, pp. 1–6.
  * [18] D. St-Onge, C. Pinciroli, and G. Beltrame, “Circle formation with computation-free robots shows emergent behavioural structure,” in _2018 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS)_. IEEE, 2018, pp. 5344–5349.
  * [19] S. Zhu, D. Wang, and C. B. Low, “Cooperative control of multiple UAVs for moving source seeking,” _Journal of Intelligent & Robotic Systems_, vol. 74, no. 1-2, pp. 333–346, 2014.
  * [20] Z. Sun, H. G. de Marina, G. S. Seyboth, B. D. O. Anderson, and C. Yu, “Circular formation control of multiple unicycle-type agents with nonidentical constant speeds,” _IEEE Transactions on Control Systems Technology_ , vol. 27, no. 1, pp. 192–205, 2018.
  * [21] C. Wang, G. Xie, and M. Cao, “Controlling anonymous mobile agents with unidirectional locomotion to form formations on a circle,” _Automatica_ , vol. 50, no. 4, pp. 1100–1108, 2014.
  * [22] ——, “Forming circle formations of anonymous mobile agents with order preservation,” _IEEE Transactions on Automatic Control_ , vol. 58, no. 12, pp. 3248–3254, 2013.
  * [23] R. Sepulchre, D. A. Paley, and N. E. Leonard, “Stabilization of planar collective motion: All-to-all communication,” _IEEE Transactions on Automatic Control_ , vol. 52, no. 5, pp. 811–824, 2007.
  * [24] J. Fonseca, J. Wei, K. H. Johansson, and T. A. Johansen, “Cooperative decentralized circumnavigation with application to algal bloom tracking,” _arXiv preprint arXiv:1903.05993_ , 2019.
  * [25] I. Shames, S. Dasgupta, B. Fidan, and B. D. O. Anderson, “Circumnavigation using distance measurements under slow drift,” _IEEE Transactions on Automatic Control_ , vol. 57, no. 4, pp. 889–903, 2011.
  * [26] M. Deghat, I. Shames, B. D. O. Anderson, and C. Yu, “Localization and circumnavigation of a slowly moving target using bearing measurements,” _IEEE Transactions on Automatic Control_ , vol. 59, no. 8, pp. 2182–2188, 2014.
  * [27] R. Zheng and D. Sun, “Circumnavigation by a mobile robot using bearing measurements,” in _2014 IEEE/RSJ International Conference on Intelligent Robots and Systems_. IEEE, 2014, pp. 4643–4648.
  * [28] W. Yao, H. Lu, Z. Zeng, J. Xiao, and Z. Zheng, “Distributed static and dynamic circumnavigation control with arbitrary spacings for a heterogeneous multi-robot system,” _Journal of Intelligent & Robotic Systems_, vol. 94, no. 3-4, pp. 883–905, 2019.
  * [29] W. Yao, Z. Zeng, X. Wang, H. Lu, and Z. Zheng, “Distributed encirclement control with arbitrary spacing for multiple anonymous mobile robots,” in _2017 36th Chinese Control Conference (CCC)_. IEEE, 2017, pp. 8800–8805.
  * [30] W. Yao, S. Luo, H. Lu, and J. Xiao, “Distributed circumnavigation control with dynamic spacing for a heterogeneous multi-robot system,” in _Robot World Cup_. Springer, 2018, pp. 374–386.
  * [31] J. Ma, W. Yao, W. Dai, H. Lu, J. Xiao, and Z. Zheng, “Cooperative encirclement control for a group of targets by decentralized robots with collision avoidance,” in _2018 37th Chinese Control Conference (CCC)_. IEEE, 2018, pp. 6848–6853.
  * [32] H. Lu, W. Yao, and L. Chen, “Distributed multi-robot circumnavigation with dynamic spacing and time delay,” _Journal of Intelligent & Robotic Systems_, pp. 1–18, 2019.
  * [33] A. H. Mong-ying and V. Kumar, “Pattern generation with multiple robots,” in _Proceedings 2006 IEEE International Conference on Robotics and Automation (ICRA)_. IEEE, 2006, pp. 2442–2447.
  * [34] A. H. Mong-ying, S. Loizou, and V. Kumar, “Stabilization of multiple robots on stable orbits via local sensing,” in _Proceedings 2007 IEEE International Conference on Robotics and Automation (ICRA)_. IEEE, 2007, pp. 2312–2317.
  * [35] M. Mesbahi and M. Egerstedt, _Graph theoretic methods in multiagent networks_. Princeton University Press, 2010, vol. 33.
  * [36] W. Yao, H. G. de Marina, B. Lin, and M. Cao, “Singularity-free guiding vector field for robot navigation,” _IEEE Transactions on Robotics_ , vol. 37, no. 4, 2021.
  * [37] 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.
  * [38] V. M. Gonçalves, 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.
  * [39] 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.
  * [40] W. Yao and M. Cao, “Path following control in 3D using a vector field,” _Automatica_ , vol. 117, p. 108957, 2020.
  * [41] W. Yao, B. Lin, B. D. O. Anderson, and M. Cao, “Topological analysis of vector-field guided path following on manifolds,” 2020, submitted. [Online]. Available: http://tiny.cc/yao_topo
  * [42] A. Galbis and M. Maestre, _Vector analysis versus vector calculus_. Springer Science & Business Media, 2012.
  * [43] W. Yao, H. G. de Marina, and M. Cao, “Vector field guided path following control: Singularity elimination and global convergence,” in _2020 59th IEEE Conference on Decision and Control (CDC)_. IEEE, 2020, pp. 1543–1549.
  * [44] D. R. Nelson, D. B. Barber, T. W. McLain, and R. W. Beard, “Vector field path following for miniature air vehicles,” _IEEE Transactions on Robotics_ , vol. 23, no. 3, pp. 519–529, 2007.
  * [45] H. G. de Marina, “Maneuvering and robustness issues in undirected displacement-consensus-based formation control,” _IEEE Transactions on Automatic Control_ , 2020.
  * [46] W. Ren and R. W. Beard, _Distributed consensus in multi-vehicle cooperative control_. Springer, 2008, vol. 27, no. 2.
  * [47] W. Yao, H. G. de Marina, Z. Sun, and M. Cao, “Guiding vector fields for multi-robot coordinated navigation,” 2021, in preparation.
  * [48] W. Yao, B. Lin, and M. Cao, “Integrated path following and collision avoidance using a composite vector field,” in _2019 IEEE 58th Conference on Decision and Control (CDC)_. IEEE, 2019, pp. 250–255.
  * [49] L. Wang, A. D. Ames, and M. Egerstedt, “Safety barrier certificates for collisions-free multirobot systems,” _IEEE Transactions on Robotics_ , vol. 33, no. 3, pp. 661–674, 2017.
  * [50] ——, “Safe certificate-based maneuvers for teams of quadrotors using differential flatness,” in _2017 IEEE International Conference on Robotics and Automation (ICRA)_. IEEE, 2017, pp. 3293–3298.
  * [51] A. D. Ames, S. Coogan, M. Egerstedt, G. Notomista, K. Sreenath, and P. Tabuada, “Control barrier functions: Theory and applications,” in _2019 18th European Control Conference (ECC)_. IEEE, 2019, pp. 3420–3431.
  * [52] X. Yun and Y. Yamamoto, “On feedback linearization of mobile robots,” Department of Computer and Information Science, University of Pennsylvania, technical report, 1992.
  * [53] B. Siciliano and O. Khatib, _Springer Handbook of Robotics_. Secaucus, NJ, USA: Springer-Verlag New York, Inc., 2007.
  * [54] 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.
  * [55] L. J. Colombo and D. V. Dimarogonas, “Motion feasibility conditions for multiagent control systems on Lie groups,” _IEEE Transactions on Control of Network Systems_ , vol. 7, no. 1, pp. 493–502, 2019.
  * [56] H. Khalil, _Nonlinear Systems_ , 3rd ed. Prentice Hall, 2002.
  * [57] Y. A. Kapitanyuk, H. G. de Marina, A. V. Proskurnikov, and M. Cao, “Guiding vector field algorithm for a moving path following problem,” _IFAC-PapersOnLine_ , vol. 50, no. 1, pp. 6983–6988, 2017.
  * [58] B. Gati, “Open source autopilot for academic research-the paparazzi system,” in _2013 American Control Conference_. IEEE, 2013, pp. 1478–1481.
  * [59] “Code repository for this paper,” 2020. [Online]. Available: https://github.com/noether/paparazzi/tree/gvf_param_multi/sw/airborne/modules/guidance/gvf_parametric