Coordinated trajectory tracking of multiple vertical take-off
and landing UAVs ⋆

arXiv:1711.03261v1 [cs.SY] 9 Nov 2017

Yao Zou a, Ziyang Meng a,
a

Department of Precision Instrument, Tsinghua University, Beijing 100084, P. R. China.

Abstract
This paper investigates the coordinated trajectory tracking problem of multiple vertical takeooff and landing (VTOL) unmanned
aerial vehicles (UAVs). The case of unidirectional information flow is considered and the objective is to drive all the follower
VTOL UAVs to accurately track the trajectory of the leader. Firstly, a novel distributed estimator is developed for each
VTOL UAV to obtain the leader’s desired information asymptotically. With the outputs of the estimators, the solution to the
coordinated trajectory tracking problem of multiple VTOL UAVs is transformed to individually solving the tracking problem
of each VTOL UAV. Due to the under-actuated nature of the VTOL UAV, a hierarchical framework is introduced for each
VTOL UAV such that a command force and an applied torque are exploited in sequence, then the position tracking to the
estimated desired position and the attitude tracking to the command attitude are achieved. Moreover, an auxiliary system
with proper parameters is implemented to guarantee the singularity-free command attitude extraction and to obviate the
use of the unavailable desired information. The stability analysis and simulations effectively validate the achievement of the
coordinated trajectory tracking of multiple VTOL UAVs with the proposed control approach.
Key words: Unmanned aerial vehicle (UAV); Coordinated trajectory tracking; Directed graph; Distributed estimators

1

Introduction

The past few decades have witnessed a rapid development in the formation control of unmanned aerial vehicles (UAVs). Replacing a single monolithic UAV with a
formation of multiple micro ones can effectively improve
efficiency without costly expense [1]. More recently, the
vertical takeoff and landing (VTOL) UAV, as a representative UAV, has received increasing interest, due
to its capacities of hovering and low-speed/low-altitude
flight [2]. Additionally, the VTOL UAV is a canonical nonlinear system with under-actuation property [3],
which raises a lot of technical problems for control theory research. Therefore, the formation control of multiple VTOL UAVs deserves intensive studies.
⋆ This paper was not presented at any IFAC meeting. This
work has been supported in part by National Natural Science
Foundation of China under Grant 61503249 and Beijing Municipal Natural Science Foundation under Grant 4173075,
and the National Key Research and Development Program
of China under Grant 2016YFB0500900/2.
Corresponding author: Ziyang Meng.
Email addresses: zouyao@mail.tsinghua.edu.cn
(Yao Zou), ziyangmeng@mail.tsinghua.edu.cn
(Ziyang Meng).

Preprint submitted to Automatica

Generally, the study of formation control problem is categorized into leaderless and leader-follower formation
control problem. The leaderless formation requires its
members to simply reach a prescribed pattern [4]. For
example, a distributed control algorithm is proposed in
[5] such that the formation of VTOL UAVs with an identical velocity was achieved. The special case with communication delay was also studied [6] for the leaderless
formation objective and corresponding control solutions
were proposed. Another formation protocol was developed in [7] to realize a time-varying formation of VTOL
UAVs without a leader, and the obtained theoretical results were verified with practical experiments.
Compared with the leaderless formation, the objective
of the leader-follower formation is that followers reach
an agreement with the desired information associated
with a leader while forming the prescribed pattern [8].
This may lead the formation to complete some intricate
missions, where the leader is responsible for performing the desired trajectory of the formation and it is
delivered via the communication network between the
leader and followers. Although [9,10,11] proposed control approaches to achieve the leader-follower formation
of VTOL UAVs, the leader’s desired information was
required to be available to all the followers. In prac-

10 November 2017

control approach project and the stability analysis. Section 4 performs some simulations to verify the theoretical results. And section 5 draws final conclusions.

tice, due to limited information exchange and communication constraints, the leader’s desired information
is only accessible to a portion of the followers. To
achieve the leader-follower formation under restricted
communication networks, distributed algorithms were
implemented with a local information exchange mechanism [12,13,14,15,16,17]. Using backstepping and filtering strategies, a distributed control algorithm was
developed in [18] to realize the asymptotically stable
leader-follower formation of VTOL UAVs. A distributed
formation and reconfiguration control approach is designed in [19] to accomplish the leader-follower formation without inter-vehicle collisions. With feedback
linearization technique, [20] proposed a leader-follower
formation protocol for VTOL UAVs, which ensured
their heading synchronization as well. [21] presented
a distributed control strategy over a switched graph
and derived necessary and sufficient conditions on the
time-varying leader-follower formation of VTOL UAVs.
However, the network graphs among the followers used
in [18,19,20,21] are undirected, which means that each
pair of the followers interacts bidirectionally. This undirected graph condition is quite restrictive, which, due
to communication constraints, can hardly be met in
practical applications.

Notations. Rm×n denotes the m × n Euclidean space.
Given a vector x = [x1 , x2 , · · · , xn ]T , define sgn(x)
=
Pn
T
[sgn(x1 ), sgn(x
),
·
·
·
,
sgn(x
)]
,
and
kxk
=
|x
2
n
1
i|
i=1
√
and kxk = xT x are its 1-norm and 2-norm. Given
a square matrix A = [aij ] ∈ Rn×n , define λmin (A) and
λmax (A)qas its minimum and maximum eigenvalues, and
Pn Pn
2
kAk =
i=1
j=1 aij is its F-norm. In is an n × n
identity matrix and 1n is an n-dimensional vector with
all entries being one. Furthermore, given a vector x =
[x1 , x2 , x3 ]T , superscript × represents a transformation
from x to a skew-symmetric matrix:


0

−x3 x2






x× = 
 x3 0 −x1  .
−x2 x1 0
2

This paper proposes a coordinated trajectory tracking
control approach for multiple VTOL UAVs with local
information exchange, where the desired trajectory information is described by a leader. The network graph
among the follower VTOL UAVs is assumed to be directed. This effectively relaxes the restrictive assumption that the graph is symmetric. By applying a novel
distributed estimator, the leader’s desired information
is accurately estimated for each follower VTOL UAV.
Based on the hierarchial framework, a command force
and an applied torque are synthesized for each VTOL
UAV such that the coordinated trajectory tracking is
achieved for a group of VTOL UAVs. Compared with
the aforementioned work, the main contributions of this
paper are three-fold. First, in contrast to the work in
[5,6,7], where only a prescribed pattern is formed with
a constant velocity, the leader-follower tracking of multiple VTOL UAVs is achieved by introducing a novel
distributed estimator. Second, the coordinated tracking
is achieved with weak connectivity, where the network
graph among the followers is directed, rather than the
limited undirected one used in [18,19,20,21]. Third, instead of solely discussing the position loop [7,21], a complete VTOL UAV system is studied based on a hierarchical framework, where an auxiliary system is proposed
to ensure the non-singular command attitude extraction
and to avoid the use of the unavailable desired information.

2.1

Background
Problem statement

Suppose that there are n follower VTOL UAVs in a team,
which are labeled by V = {1, 2, . . . , n}. Each UAV is
a six-dof (degree of freedom) rigid body and operates
in two reference frames: inertia frame I = {OI xI yI zI }
which is attached to the earth and body frame B =
{OB xB yB zB } which is fixed to the fuselage. To establish
the model of the UAVs, rotation matrix Ri ∈ SO(3) ,
{R ∈ R3×3 | det R = 1, RT R = RRT = I3 } and unit
quaternion Qi = [σi , qiT ]T ∈ Q , {Q ∈ R × R3 | σ 2 +
q T q = 1} are applied to represent the attitude of each
UAV. In terms of Euler formula [22], an explicit relation
between these two attitude representations is derived as
Ri (Qi ) = (σi2 − qiT qi )I3 + 2qi qiT − 2σi qi× .

(1)

Based on Euler-Newton formulae, the kinematics and
dynamics of the i-th VTOL UAV are given by
ṗi = vi ,
Ti
v̇i = −gê3 +
Ri (Qi )ê3 ,
mi
1
Q̇i = Gi (Qi )ωi ,
2
Ji ω̇i = −ωi× Ji ωi + τi ,

(2)
(3)
(4)
(5)

where pi = [pxi , pyi , pzi ]T and vi = [vix , viy , viz ]T denote
the position and velocity of the center of gravity of the
UAV in frame I, respectively, mi is the total mass, g
is the local gravitational acceleration, ê3 , [0, 0, 1]T ,

The remaining sections are arranged as follows. Section
2 describes the problem to be solved and provides some
useful preliminaries. Section 3 states the main results
in detail, including the distributed estimator design, the

2

Ti denotes the applied thrust along ê3 , Qi = [σi , qiT ]T
and Ri (Qi ) are the unit quaternion and rotation matrix, Gi (Qi ) = [−qi , σi I3 − qi× ]T , ωi = [ωix , ωiy , ωiz ]T
denotes the angular velocity of the UAV in frame B,
Ji = diag{Jix , Jiy , Jiz } is the inertial matrix with respect
to frame B, and τi denotes the applied torque in frame B.

di0 > 0 if node i is accessible
to #the leader and di0 = 0
"
0 01×n
otherwise; and L̄ ,
, where M = [mij ] ,
−d0 M
L + diag{d10 , d20 , · · · , dn0 }.
Assumption 2.2 The leader-follower graph Gn+1 has a
directed spanning tree with the leader being the root.

In addition to n followers, there is a leader, labeled by
0, to represent the global desired information including
the desired position pr and its derivatives. The control
objective is to design applied thrust Ti and torque τi
for each follower VTOL UAV described by (2)-(5) such
that all the followers track the leader while maintaining
a prescribed formation pattern. More specifically, given
a desired position offset δi , the objective is to guarantee
that

Some important properties associated with matrix M
are given in Lemma 2.1 [23].
Lemma 2.1 Under Assumption 2.2, M is a nonsingular M-matrix with the properies that all its eigenvalues have positive real parts, and there exists a positive definite diagonal matrix Θ = diag{θ1 , θ2 , · · · , θn } such that
Ξ = MT Θ+ΘM is strictly diagonally dominant and positive definite, where [1/θ1 , 1/θ2 , · · · , 1/θn ]T = M−1 1n .

lim (pi (t)−pr (t)) = δi , lim (vi (t)− ṗr (t)) = 0, ∀i ∈ V.

t→∞

t→∞

(6)
Due to communication constraints, the leader’s desired
information is only available to a subset of the followers and the followers only have access to their neighbors’ information. To solve such a coordinated tracking
problem via local information exchange, distributed algorithms are implemented. Moreover, it follows from (6)
that limt→∞ (pi (t) − pj (t)) = δij , where δij = δi − δj ,
∀i, j ∈ V. This means that the followers form a pattern determined by δi while tracking the leader. Therefore, a proper position offset δi is essential such that the
proposed algorithm ensures followers’ convergence to a
well-defined and unique formation.

2.3

Consider the vector differential equation
ẋ = f (x, t),

(7)

where f : Rn × R → Rn is measurable and essentially
locally bounded.
In what follows, the definitions of Filippov solution, generalized gradient and regular function are given according to [24,25,26].

Assumption 2.1 The desired position pr and its
(3)
derivatives ṗr , p̈r and pr are bounded.
2.2

Filippov solution and non-smooth analysis

Definition 2.1 (Filippov solution) A vector function x(t) is called a solution of (7) on [t0 , t1 ], if x(t)
is absolutely continuous on [t0 , t1 ] and for almost all
t ∈ [t0 , t1 ], ẋ ∈ K[f ](x, t), where

Graph theory

K[f ](x, t) =

Communication topology among UAVs is described by a
graph Gn , (V, E), which is composed of a node set V ,
{1, 2, · · · , n} and an edge set E ⊆ V × V. For a directed
graph, (i, j) ∈ E means that the information of node i is
accessible to node j, but not conversely. All neighbours
of node i are included in set Ni = {j ∈ V | (j, i) ∈ E}.
A path from node i to node j is a sequence of edges.

\ \

ρ>0 µN =0

cof (B(x, ρ) − N, t),

T

µN =0 denotes the intersection over all sets N of
Lebesgue measure zero, co(·) denotes the vector convex
closure, and B(x, ρ) denotes the open ball of radius ρ
centered at x.

For a follower graph Gn , its adjacent matrix D = [dij ] ∈
Rn×n is defined such that dij > 0 if (j, i) ∈ E and dij =
0 otherwise, and the associated nonsymmetric Laplan×n
cian
is defined such that lii =
Pn matrix L = [lij ] ∈ R
j=1,j6=i dij and lij = −dij for j 6= i. For a leaderfollower graph Gn+1 , {V̄, Ē} (leader is labeled as 0)
with V̄ = {0, 1, · · · , n} and Ē ⊆ V̄ × V̄, we define D̄ ∈
R(n+1)×(n+1) and L̄ ∈ R(n+1)×(n+1) as its adjacent matrix and
Laplacian matrix. Specifically,
" nonsymmetric
#
0 01×n
D̄ ,
, where d0 = [d10 , d20 , · · · , dn0 ]T and
d0 D

Definition 2.2 (Generalized gradient) For a locally
Lipschitz function V : Rn × R → R, its generalized gradient at (x, t) is defined as
∂V (x, t) = co{lim ▽V (x, t) | (xi , ti ) → (x, t), (xi , ti ) ∈
/ ΩV },
where ΩV is the set of measure zero where the gradient of V is not defined. Furthermore, the generalized derivative of V along system (7) is defined as
"
#
T
K[f
](x,
t)
˙
.
Ṽ , φ∈∂V φT
1
3

Definition 2.3 (Regular) f (x, t) : Rn × R is called
regular if
(1) for all ν ≥ 0, the usual one-sided directional derivative f ′ (x; ν) exists;
(2) for all ν ≥ 0, f ′ (x; ν) = f o (x; ν), where the generalized directional derivative f o (x; ν) is defined as
f o (x; ν) = lim sup
y→x t↓0

command attitude extraction. For i ∈ V, a distributed
estimator is proposed as follows:
p̂˙i =v̂i − kp
v̂˙ i =âi − kv

f (y + tν) − f (y)
.
t

Lemma 2.2 Let system (7) be essentially locally
bounded and 0 ∈ K[f ](x, t) in a region Rn × [0, ∞).
Suppose that f (0, t) is uniformly bounded for all t ≥ 0.
Let V : Rn × [0, ∞) → R be locally Lipschitz in t
and regular such that W1 (x) ≤ V (t, x) ≤ W2 (x) and
Ṽ˙ (x, t) ≤ −W (x), where W (x) and W (x) are continu-

(8a)

dij (v̂i − v̂j ),

(8b)

j=0

j=0

(8c)

where p̂0 = pr , v̂0 = ṗr and â0 = p̈r are specified,
dij is the (i, j)-th entry of the adjacent matrix D associated with the follower graph Gn , kp , kv , ka and la
are positive parameters, and Γi = diag{µxi , µyi , µzi } with
µki = tanh(γ̇ik )γ̇ik and Γ̄i = diag{µ̄xi , µ̄yi , µ̄zi } with µ̄ki =
1 − tanh2 (γik ) for k = x, y, z. Next, define the estimation
errors p̄i = p̂i −pr , v̄i = v̂i − ṗr and āi = âi − p̈r for i ∈ V.
It then follows from (8) that their dynamics satisfy

2

that ϕ ≤ −W , ∀ϕ ∈ Ṽ˙ . Then, all Filippov solutions of
system (7) are bounded and satisfy limt→∞ W (x(t)) = 0.
Main results

p̄˙ i =v̄i − kp

Due to the under-actuated nature of the VTOL UAV, a
hierarchical strategy is applied to solve the coordinated
trajectory tracking problem of multiple VTOL UAV systems. First, a distributed estimator using local information interaction is designed for each follower UAV to estimate the leader’s desired information. Then, the coordinated trajectory tracking problem of multiple VTOL
UAVs is transformed into the asymptotic stability problem of each individual error system. Next, a command
force and an applied torque are exploited for each UAV
to asymptotically stabilize the position and attitude error systems, respectively. Finally, the stability of each
error system is analyzed.
3.1

dij (p̂i − p̂j ),

j=0

ous positive definite functions, W (x) is a continuous positive semi-definite function, and Ṽ˙ (x, t) ≤ −W (x) means

3

j=0
n
X

γ̈i = − la γ̇i + 2Γi γ̇i


n
n
X
ka −1 X
− Γ̄i
dij (â˙ i − â˙ j )
dij (âi −âj )+
kγ
j=0
j=0


n
n
X
X
dij (â˙ i − â˙ j ) ,
+sgn 
dij (âi − âj ) +

The Lyapunov stability criterion for non-smooth systems is given in Lemma 2.2 [27].

1

n
X

v̄˙ i =āi − kv

n
X

j=1
n
X

mij p̄j ,

(9a)

mij v̄j ,

(9b)

j=1

¨i =kγ Γ̄i γ̈i − 2kγ Γ̄i Γi γ̇i − p(4)
ā
r
!
n
n
X
X
= − la kγ Γ̄i γ̇i − ka
mij āj +
mij ā˙ j
j=1

− ka sgn

n
X

mij āj +

mij ā˙ j

j=1

j=1

= − la ā˙ i − ka

j=1

n
X

n
X

mij āj +

− ka sgn

Distributed estimator design

j=1

Since the leader’s desired information including the desired position pr and its derivatives is not available to all
the followers, a distributed estimator is firstly designed
for each VTOL UAV to estimate them.

mij āj +

− p(4)
r

mij ā˙ j

j=1

j=1

n
X

n
X

!

n
X
j=1

mij ā˙ j

!

!

(9c)

+ Np ,

where mij denotes the (i, j)-th entry of M defined in
(4)
(3)
Section 2, and Np = la pr − pr is bounded according
to Assumption 2.1. Equivalently, the error dynamics (9)
can be rewritten as

For i ∈ V, we define p̂i , v̂i and âi = kγ tanh(γi ) as the
estimations of pr , ṗr and p̈r , respectively, where γi is an
auxiliary variable and parameter kγ ≥ supt≥0 p̈r (t). As
will be shown subsequently, the definition of âi using the
hyperbolic tangent function enables the control parameters to be chosen explicitly in case of singularity in the

p̄˙ =v̄ − kp (M ⊗ I3 )p̄,
v̄˙ =ā − kv (M ⊗ I3 )v̄,
˙
¨ = − la ā˙ − ka (M ⊗ I3 )(ā + ā)
ā
˙ + 1n ⊗ Np ,
− ka sgn ((M ⊗ I3 )(ā + ā))

4

(10a)
(10b)
(10c)

According to Definition 2.2, the generalized derivative
of Ls along (11) satisfies

where p̄, v̄ and ā are the column stack vectors of p̄i , v̄i
and āi , respectively, and operator ⊗ denotes the kronecker
define a sliding surface si =
Pn product. Moreover,
Pn
la j=1 mij āj + j=1 mij ā˙ j for i ∈ V, and correspond˙ It
ingly, its column stack vector s = (M ⊗ I3 )(la ā + ā).
follows from (10c) that the dynamics of s satisfies
ṡ = (M ⊗ I3 ) (−ka s − ka sgn(s) + 1n ⊗ Np ) .

\
L̃˙ s = (φ + s)T (Θ ⊗ I3 )K[(M ⊗ I3 )(−ka s − ka sgn(s)
φ∈∂Ls

+ 1n ⊗ Np )],
\
= (φ + s)T (ΘM⊗ I3 )(−ka s− ka ∂ksk1 + 1n ⊗ Np ),

(11)

φ∈∂Ls

(15)

Theorem 3.1 indicates that the developed distributed
estimator (8) with appropriate parameters enables the
estimation errors p̄i , v̄i and āi for each VTOL UAV to
converge to zero asymptotically.



sji ∈ R+
{1}
j
where ∂|si | = {−1} sji ∈ R− , ∀i ∈ V, j = x, y, z,

[−1, 1] sj = 0
i
and the calculation of K is applied using the same argument given in [24].

Theorem 3.1 Under Assumptions 2.1 and 2.2, if the
estimator parameters kp , kv , la and ka are chosen based
on
kΘk2
,
λmin (Ξ)2
kp λmin (Ξ)kΘk2
la >
,
λmin (MT ΘM)(kp kv λmin (Ξ)2 − kΘk2)
√
2 nkΘMkN̄p
ka >
,
λmin (Ξ)
kp kv >

If L̃˙ s 6= ∅, suppose ϕ ∈
/ L̃˙ s , then we know that ∀φ ∈
T
∂ksk1, ϕ = (φ + s) (ΘM ⊗ I3 )(−ka s − ka δ + 1n ⊗ Np )
for some δ ∈ ∂ksk1 . Choose φ = arg minδ∈∂ksk1 (δ +
T
Θ
s)T ( ΘM+M
⊗ I3 )(δ + s). According to [24], for all
2
δ ∈ ∂ksk1 , we have that

(12)
(13)
(14)

(φ + s)T (ΘM ⊗ I3 )(δ + s)


ΘM + MT Θ
T
⊗ I3 (φ + s)
≥(φ + s)
2
1
= (φ + s)T (Ξ ⊗ I3 )(φ + s).
2

where N̄p = supt≥0 kNp (t)k, and Θ and Ξ are given in
Lemma 2.1, the distributed estimator (8) ensures that
limt→∞ p̄i (t) = 0, limt→∞ v̄i (t) = 0 and limt→∞ āi (t) = 0,
∀i ∈ V.
Proof: The proof is divided into two parts: first, the
sliding surface si , i ∈ V is proven in Proposition 3.1 to
converge to zero asymptotically; then, the final result is
shown in Proposition 3.2.

It then follows that L̃˙ s further satisfies
ka
L̃˙ s ≤ − (φ + s)T (Ξ ⊗I3 )(φ + s) +(φ + s)T (1n ⊗Np )
2
√
ka λmin (Ξ)
≤−
kφ + sk2 + nkφ + skkΘMkN̄p
2


√
ka λmin (Ξ)
kφ + sk − nkΘMkN̄p kφ + sk.
=−
2
(16)

Proposition 3.1 Under Assumptions 2.1 and 2.2, if the
estimator parameter ka satisfies (14), the distributed estimator (8) guarantees that limt→∞ si (t) = 0, ∀i ∈ V.
Proof: Obviously, system (11) is non-smooth; thereafter,
the solution of (11) is studied in the sense of Filippov
and the non-smooth framework is applied. The stability
of system (11) is to be proven based on Lemma 2.2.

Note that if φ = 0, then s = 0 and kφ+sk = 0, and if φ 6=
0, then kφ + sk ≥ 1. Hence, if the estimator parameter
ka√ satisfies (14), there exists a constant k̄a satisfying
2 nkΘMkN̄p
(Ξ)
≤ k̄a < ka such that −( k̄a λmin
kφ + sk −
2
√ λmin (Ξ)
nkΘMkN̄p)kφ + sk ≤ 0. Therefore, it follows that

Pn
We first propose a Lyapunov function Ls = i=1 θi ksi k1+
sT (Θ ⊗ I3 )s, where θi is the i-th diagonal entry of Θ.
Note that Ls is non-smooth but regular [24]. It can be
derived that Ls is bounded by

(ka − k̄a )λmin (Ξ)
L̃˙ s ≤ −
kφ + sk2 .
2

s

W1 (s) ≤ L ≤ W2 (s),

(17)

In addition, each entry of s has the same sign as its
counterpart in sgn(s), and thus, it follows that kφ+sk ≥
ksk. We finally have that

√
where W1 (s) = λmin (Θ)ksk2 and W2 (s) = nλmax (Θ)
· (ksk1 + ksk2 ). In terms of Lemma 2.2, the stable result
can be deduced if only the generalized derivative of Ls
satisfies L̃˙ s ≤ W (s), where W (s) is a continuous positive
semi-definite function.

(ka − k̄a )λmin (Ξ)
ksk2 = −W (s).
L̃˙ s ≤ −
2

5

(18)

this case, L̇e further satisfies

Since W (s) ≥ 0 has been ensured, it follows that
Rt
W (s(τ ))dτ is bounded for t ≥ 0, which implies that s
0
is bounded. Hence, it follows from (11) that ṡ is bounded,
which implies that s is uniformly continuous in t. This
means that W (s) is uniformly continuous in t. Based on
Lemma 2.2, we can conclude that limt→∞ W (s(t)) = 0,
which further implies that limt→∞ s(t) = 0.

L̇e ≤ − λmin (Ω)kzk2 + kΘMkkzkksk
√
≤ − ϑ1 Le + ϑ2 ksk Le ,
where ϑ1 =

(19)

1
Le = p̄T (Θ ⊗ I3 )p̄ + v̄ T (Θ⊗I3 )v̄ + āT (MT ΘM ⊗ I3 )ā.
2

λ1 kzk2 ≤ Le ≤ λ2 kzk2,
where λ1 = min{λmin (Θ), 21 λmin (MT ΘM)} and λ2 =
max{λmax (Θ), 12 λmax (MT ΘM)}. The derivative of Le
along (10) and (19) satisfies
L̇e = − p̄T ((MT Θ + ΘM) ⊗ I3 )p̄ + 2p̄T (Θ ⊗ I3 )v̄

3.2

Since the leader’s desired information has been estimated via the distributed estimator (8) for each follower
VTOL UAV, the remaining problem is to transform the
coordinated trajectory tracking problem into the simultaneous tracking problem for the decoupled VTOL UAV
group. This is verified as follows.

+ āT (MT Θ ⊗ I3 )(−la (M ⊗ I3 )ā + s)

= − kp p̄T (Ξ ⊗ I3 )p̄ + 2p̄T (Θ ⊗ I3 )v̄ − kv v̄ T (Ξ ⊗ I3 )v̄
+ 2v̄ T (Θ ⊗ I3 )ā − la āT (MT ΘM ⊗ I3 )ā

+ āT (MT Θ ⊗ I3 )s
≤ − kp λmin (Ξ)kp̄k2 + 2kΘkkp̄kkv̄k − kv λmin (Ξ)kv̄k2
+ 2kΘkkv̄kkāk − la λmin (MT ΘM)kāk2
+ kΘMkkākksk

≤ − z1T Ωz1 + kΘMkkakksk,

Define the position error pei = pi − pr − δi and the velocity error vie = vi − ṗr for i ∈ V. Using the estimations p̂i and v̂i obtained from the distributed estimator
(8), we rewrite pei and vie as pei = pi − p̂i − δi + p̄i and
vie = vi − v̂i + v̄i . Since Theorem 3.1 has shown that
limt→∞ p̄i (t) = 0 and limt→∞ v̄i (t) = 0, ∀i ∈ V, the
coordinated tracking control objective (6) can be transformed into the following simultaneous tracking objective:

(20)

where z1 = [kp̄k, kv̄k, kāk]T and


Ω=


−kΘk
0

Problem transformation

T

− v̄ ((M Θ + ΘM) ⊗ I3 )v̄ + 2v̄ (Θ ⊗ I3 )ā

kp λmin (Ξ)

(22)

Remark 3.1 According to (8), singularity may occur in
the distributed estimator when some diagonal entry of
Γ̄i equals to zero, and this corresponds to the case where
some entry of the auxiliary variable γi tends to infinity. Theorem 3.1 has shown that, with a bounded initial
value, the estimation error āi for each UAV is bounded
all the time. This implies that Γ̄i is always positive definite. Consequently, no singularity is introduced in the
developed distributed estimator (8).

It is bounded by



and kzk = kz1 k has

When V e = 0, it can be shown that D+ V e ≤ ϑ2 ksk.
Thus, it follows that D+ V e satisfies (22) all the time
[28]. For system ẏ = −ϑ1 y + ϑ2 ksk with respect to
y ∈ [0, ∞), it can be proven in terms of input-to-state
stability theory [28] that limt→∞ y(t) = 0 given the fact
that limt→∞ s(t) = 0. According to Comparison Principal [28], it follows that limt→∞ V e (t) = 0, which further implies that limt→∞ p̄(t) = 0, limt→∞ v̄(t) = 0 and
limt→∞ ā(t) = 0.

Define z = [p̄T , v̄ T , āT ]T and assign a Lyapunov function

T

kΘMk
√
,
λ1√
e

V̇ e ≤ −ϑ1 V e + ϑ2 ksk.

Proof: Consider the definition of the sliding surface s,
then the dynamics of the estimator error ā satisfies

T

ϑ2 =

been applied. Next, take V = 2Le . When Le 6= 0, it
follows from (21) that the derivative of V e satisfies

Proposition 3.2 Under Assumptions 2.1 and 2.2, if
the estimator parameters kp , kv and la satisfy (12) and
(13), then limt→∞ si (t) = 0 is sufficient to ensure that
limt→∞ p̄i (t) = 0, limt→∞ v̄i (t) = 0 and limt→∞ āi (t) =
0, ∀i ∈ V.

(M ⊗ I3 )ā˙ = −la (M ⊗ I3 )ā + s.

λmin (Ω)
,
λ2

(21)

−kΘk

kv λmin (Ξ)
−kΘk

0
−kΘk

−la λmin (MT ΘM)




.


lim (pi (t) − p̂i (t)) = δi , lim (vi (t) − v̂i (t)) = 0, ∀i ∈ V.

t→∞

t→∞

We next define p̄ei = pi − p̂i − δi and v̄ie = vi − v̂i for
i ∈ V. It follows from (2), (3) and (8) that their dynamics

If the estimator parameters kp , kv and la are chosen
based on (12) and (13), then Ω is positive definite. In

6

Based on the above discussions, by introducing the distributed estimator (8), the coordinated trajectory tracking problem for multiple VTOL UAV systems (2)-(5) can
be transformed into the simultaneous asymptotic stability problem for each error system (23), (27) and (28).
Lemma 3.2 summarizes this point.

satisfy
p̄˙ ei =v̄ie +

n
X

mij p̄j ,

(23a)

j=1

v̄˙ ie =ui − gê3 − âi +

n
X

mij v̄j +

j=1

Ti
(Ri (Qi )− Ri (Qci ))ê3 ,
mi

Lemma 3.2 Consider the i-th error system (23),
(27) and (28). If a command force ui and an applied
torque τi can be developed such that limt→∞ p̄ei (t) = 0,
limt→∞ v̄ie (t) = 0, limt→∞ qie (t) = 0 and limt→∞ ωie (t) =
0, the coordinated trajectory tracking of multiple VTOL
UAV systems (2)-(5) is achieved in the sense of (6).

(23b)

Ti
Ri (Qci )ê3 is the command force with Qci
where ui = m
i
being the command unit quaternion. Moreover, once
the command force ui can be determined, in view of
kRi (Qci )ê3 k = 1, the applied thrust Ti is derived as

Ti = mi kui k, ∀i ∈ V.

3.3

(24)

In this subsection, a command force ui for each VTOL
UAV will be synthesized. The main difficulties here are
that the command force ui should comply with the nonsingular condition (25) and the desired position pr and
its derivatives are not available in the command force
ui and the subsequent applied torque τi due to limited
communication.

Now that the control strategy is based on a hierarchical
framework, the command unit quaternion Qci , as the attitude tracking objective for each VTOL UAV, should
be extracted from the command force ui . Based on minimal rotation principle, a viable extraction algorithm is
proposed in Lemma 3.1 [5].

To address the above dilemmas, we introduce the virtual
position error p̃ei = p̄ei − ηi and the virtual velocity error
ṽie = v̄ie − η̇i for i ∈ V, where ηi is an auxiliary variable. It
follows from (23) that the dynamics of p̃ei and ṽie satisfy

Lemma 3.1 For i ∈ V, if the command force ui =
[uxi , uyi , uzi ]T satisfies the non-singular condition:
ui ∈
/ U , {u ∈ R3 | u = [0, 0, uz ]T , uz ≤ 0},
the command unit quaternion
tracted as

σic =

s

1 g − uzi
+
, qc =
2
2kui k i

Qci

=

[σic , (qic )T ]T


uyi

1
 −ux
i
2kui kσic 
0




.


(25)

p̃˙ei = ṽie + ~1 ,
ṽ˙ ie = ui − η̈i − gê3 − âi + ~2 ,

is ex-

1
Gi (Qei )ωie ,
2

(29a)
(29b)

Pn
Pn
where ~1 =
j=1 mij p̄j and ~2 =
j=1 mij v̄j +
Ti
c
mi (Ri (Qi ) − Ri (Qi ))ê3 .

(26)

Lemma 3.3 Consider the i-th virtual position error
system (29). If a command force ui can be synthesized such that limt→∞ p̃ei (t) = 0, limt→∞ ṽie (t) =
0, limt→∞ ηi (t) = 0 and limt→∞ η̇i (t) = 0, then
limt→∞ p̄ei (t) = 0 and limt→∞ v̄ie (t) = 0 are achieved.

Next, define the attitude error Qei = [σie , (qie )T ]T =
(Qci )−1 ⊙ Qi for i ∈ V, where operator ⊙ is the
unit quaternion product. According to [29], Qei =
[±1, 0, 0, 0]T corresponds to the extract attitude tracking. The dynamics of Qei satisfies
Q̇ei =

Command force development

To guarantee the condition in Lemma 3.3, for i ∈ V, we
propose the following command force:
ui = gê3 + âi − kη (tanh(ηi + η̇i ) + tanh(η̇i )),

(30)

(27)

and introduce a dynamic system with respect to the auxiliary variable ηi :

where ωie = ωi − Ri (Qei )T ωic is the angular velocity error
with ωic being the command angular velocity. Please refer
to [29] for the derivations of ωic and its derivative ω̇ic . In
addition, it follows from (5) and Ṙi (Qei ) = Ri (Qei )(ωie )×
that the dynamics of ωie satisfies

η̈i = −kη (tanh(ηi + η̇i ) + tanh(η̇i )) + lp p̃ei + lv ṽie ,
(31)

Ji ω̇ie = − ωi× Ji ωi + τi + Ji [(ωie )× Ri (Qei )T ωic
−

Ri (Qei )T ω̇ic ].

where kη , lp and lv are positive control parameters. Substituting (30) and (31) into (29) yields
p̃˙ei = ṽie + ~1 ,
ṽ˙ e = − lp p̃ei − lv ṽie + ~2 .

(28)

7

(32a)
(32b)

A proper control parameter kη should be chosen such
that the non-singular condition (25) is met. Specifically,
kη <

1
(g − kγ ).
2

functions of the derivatives of the command force ui .
Their expressions are presented as follows:
u̇i =kγ Γ̄i γ̇i − kη [Di η̇i + (Di + Si )η̈i ],
üi = − 2kγ Γ̄i Γi γ̇i + kγ Γ̄i γ̈i − kη [D̄i η̇i + (D̄i + Di + S̄i )η̈i

(33)

(3)

In such a case, the third row of the command force ui
satisfies

+ (Di + Si )ηi ],

where Γi and Γ̄i have been specified below (8),
Di = diag{ǫxi , ǫyi , ǫzi } with ǫji = 1 − tanh2 (ηij + η̇ij ),
Si = diag{νix , νiy , νiz } with νij = 1 − tanh2 (η̇ij ),
D̄i = {ǭxi , ǭyi , ǭzi } with ǭji = −2 tanh(ηij + η̇ij )(1 −
tanh2 (ηij + η̇ij ))(ηij + η̇ij ) and S̄i = {ν̄ix , ν̄iy , ν̄iz } with
ν̄ij = −2 tanh(η̇ij )(1 − tanh2 (η̇ij ))η̇ij , for j = x, y, z, and

uzi =g + âzi − kη (tanh(ηiz + η̇iz ) + tanh(η̇iz ))
≥g − kγ − 2kη > 0,
where âzi = kγ tanh(γiz ) and the property that
| tanh(·)| < 1 have been applied. To this end, kη satisfying (33) is sufficient to guarantee that the developed
command force ui in (30) for each UAV strictly satisfies
the non-singular condition (25).

(3)

ηi

Remark 3.2 By defining âi = kγ tanh(γi ) in the distributed estimator (8) and introducing the auxiliary dynamics (31), the developed command force ui for i ∈ V is
equipped with a saturation property. Based on this property, the choice of the control parameter kη is independent on any estimator state.

From the above derivations, it it trivial to see that the
desired information is not used in the developed applied
torque τi for the UAV without accessibility to the leader.
3.5

Remark 3.3 It follows from (24), (30) and kâi k ≤
√
3kγ that the resulting applied thrust Ti is bounded by
√
√
Ti ≤ mi (g + 2 3kη + 3kγ ), ∀i ∈ V,

(34)

Theorem 3.2 Consider n follower VTOL UAV systems
(2)-(5) with Assumptions 2.1 and 2.2. The synthesized
command force ui in (30) and applied torque τi in (36)
guarantee the coordinated trajectory tracking of multiple
VTOL UAVs in the sense of (6).

Applied torque development

Proof: Theorem 3.1 has shown that, for i ∈ V, the
distributed estimator developed in (8) enables the estimation errors p̄i and v̄i to converge to zero asymptotically. Based on this, it follows from Lemmas 3.2 and
3.3 that the coordinated trajectory tracking objective is
achieved, if the following results are guaranteed by the
synthesized command force ui and applied torque τi :
Th2.i) limt→∞ qie (t) = 0 and limt→∞ ωie (t) = 0, ∀i ∈ V,
Th2.ii) limt→∞ p̃ei (t) = 0 and limt→∞ ṽie (t) = 0, ∀i ∈ V,
Th2.iii) limt→∞ ηi (t) = 0 and limt→∞ η̇i (t) = 0, ∀i ∈ V.
They will be proven in Propositions 3.3-3.5, respectively.

Define a sliding surface ri = lq qie + ωie for i ∈ V, where
lq > 0. From (27) and (28), the dynamics of ri satisfies
lq
Ji ṙi = Ji (σie I3 + (qie )× )ωie − ωi× Ji ωi + τi
2
+ Ji [(ωie )× Ri (Qei )T ωic − Ri (Qei )T ω̇ic ].

(35)

We propose an applied torque τi for each UAV as follows:
lq
Ji (σie I3 + (qie )× )ωie + ωi× Ji ωi
2
− Ji [(ωie )× Ri (Qei )T ωic − Ri (Qei )T ω̇ic ],

τi = − kq ri −

(36)
Proposition 3.3 Consider the attitude error system
(27) and (28). The developed applied torque τi in (36)
guarantees that limt→∞ qie (t) = 0 and limt→∞ ωie (t) = 0,
∀i ∈ V.

where kq > 0. Substituting (36) into (35) yields
Ji ṙi = −kq ri .

Stability analysis

Theorem 3.2 summarizes the final stability result associated with the coordinated trajectory tracking of multiple VTOL UAV systems (2)-(5) controlled by the developed applied thrust and torque.

which means that each Ti is upper bounded by a constant
associated with the individual mass mi and the specified
parameters kη and kγ .
3.4

= − kη [Di η̇i + (Di + Si )η¨i ] + lp p̃˙ ei + lv ṽ˙ ie .

(37)

Proof: It follows from (37) that, for i ∈ V, the developed applied torque τi enables the sliding surface ri
to converge to zero asymptotically. Then, assign a nonnegative function yi = 21 [kqie k2 + (1 − σie )2 ] = 1 − σie ≤ 2

It follows from (36) that, the command angular velocity
ωic and its derivative ω̇ic are necessary to determine each
applied torque τi . According to [29], ωic and ω̇ic are the

8

for i ∈ V. With the definition of ri , the derivative of yi
along (27) satisfies

It is trivial to verify that
1
1
(k tanh(ηi + η̇i )k2 + k tanh(η̇i )k2 + kη̇i k2 ).
2
kη
(39)
The derivative of Lηi along (31) satisfies
Lηi ≥

1
lq
1
ẏi = (qie )T ωie = − kqie k2 + (qie )T ri
2
2
2
lq e 2 1
2
≤ − kqi k + kri k
4
4
lq
1
= (−2yi + yi2 ) + kri k2 .
4
4

L̇ηi = tanh(ηi + η̇i )T η̇i + [tanh(ηi + η̇i ) + tanh(η̇i )
1
+ η̇i ]T [−kη (tanh(ηi + η̇i ) + tanh(η̇i )) + εi ]
kη
= − kη k tanh(ηi + η̇i ) + tanh(η̇i )k2 − η̇iT tanh(η̇i )
1
+ [tanh(ηi + η̇i ) + tanh(η̇i ) + η̇i ]T εi
(40)
kη
≤ − kη k tanh(ηi + η̇i ) + tanh(η̇i )k2 − η̇iT tanh(η̇i )
q
+ 2 Lηi kεi k
q
≤2 Lηi ρ(kεi (0)k, 0).
(41)

l

It can be shown that system ẏ = 4q (−2y + y 2 ) with
respect to y ∈ [0, 2] is asymptotically stable. For system
l
ẏ = 4q (−2y + y 2 ) + 41 kri k2 with respect to y ∈ [0, 2], by
using input-to-state stability theory [28], it follows that
limt→∞ y(t) = 0 given the fact limt→∞ ri (t) = 0. According to Comparison Principal [28], limt→∞ yi (t) = 0
is obtained, that is, limt→∞ qie (t) = 0, which, together with limt→∞ ri (t) = 0, further implies that
limt→∞ ωie (t) = 0, ∀i ∈ V.
Proposition 3.4 Consider the virtual position error system (29) with Assumptions 2.1 and 2.2. If
limt→∞ qie (t) = 0, limt→∞ p̄i (t) = 0 and limt→∞ v̄i (t) =
0 are achieved, the developed command force ui in (30)
guarantees that limt→∞ p̃ei (t) = 0 and limt→∞ ṽie (t) = 0,
∀i ∈ V.

Integrating both sides of (41), we obtain that
q
q
Lηi (t) − Lηi (0) ≤ ρ(kεi (0)k, 0)t, ∀t ≥ 0,

which indicates that Lηi cannot escape to infinity in finite
time. In addition, it follows from (40) that L̇ηi satisfies

e
Proof: It follows from (1) that (Ri (Qei ) − I3 )ê3 = ϕ×
i qi ,
ey ex
e T
c
where ϕi = [−qi , qi , σi ] . In terms of kRi (Qi )k = 1,
kQei k = 1 and (34), it follows from limt→∞ qie (t) = 0
Ti
Ti
that each m
(Ri (Qi )−Ri (Qci ))ê3 = m
Ri (Qci )(Ri (Qei )−
i
i
I3 )ê3 converges to zero asymptotically. This, together
with limt→∞ p̄i (t) = 0 and limt→∞ v̄i (t) = 0, guarantees
that the perturbation items ~1 and ~2 in the virtual position error system (32) converge to zero asymptotically.
Furthermore, it can be shown that system

L̇ηi ≤ − tanh(η̄i )T Λ tanh(η̄i )


1
+ kηi + η̇i k + (1 + )kη̇i k ρ(kεi (0)k, t)
kη
≤ − min{1, kη }kD tanh(η̄i )k2 + c2 kη̄i kρ(kεi (0)k, t)
≤ − c1 k tanh(η̄i )k2 + c2 kη̄i kρ(kεi (0)k, t),
(43)
where η̄i =

p̃˙ei = ṽie ,
ṽ˙ e = − lp p̃ei − lv ṽie ,

0

tanh(χ)T dχ+

Z

0

η̇i

tanh(χ)T dχ+

#

, c1 =

Λ = I3 ⊗

"

kη

kη

kη kη + 1

#

,D=

√
3− 5
2

c (χ )2 ∆

tanh(∆ )

L̇ηi ≤ −c1 (χη )2 kη̄i k2 + c2 kη̄i kρ(kεi (t1 )k, t − t1 ). (44)

Proof: Denote εi = lp p̃ei + lv ṽie for i ∈ V. It follows
from limt→∞ p̃ei (t) = 0 and limt→∞ ṽie (t) = 0 that
limt→∞ εi (t) = 0. To this end, there exists a KL-class
function ρ(kεi (0)k, t) such that kεi (t)k ≤ ρ(kεi (0)k, t).
For i ∈ V, the following Lyapunov function is proposed:
ηi +η̇i

11

+ η̇iT , η̇iT ]T ,

η
t1 ) ≤ 1 cη2 η for t ∈ [t1 , ∞), where χη <
is a
∆η
constant. In particular, k tanh(η̄i )k ≤ χη kη̄i k. Thus, for
t ≥ t1 , L̇ηi satisfies

Proposition 3.5 Consider the auxiliary system (31). If
limt→∞ p̃ei (t) = 0 and limt→∞ ṽie (t) = 0 are achieved,
then limt→∞ ηi (t) = 0 and limt→∞ η̇i (t) = 0, ∀i ∈ V.

Z

"

[ηiT

min{1, kη } and c2 = 1+ k1η . Since
12
Lηi cannot escape in finite time, there exist t1 and ∆η
such that kη̄i (t)k ≤ ∆η for t ∈ [0, t1 ] and ρ(kεi (t1 )k, t −

I3 ⊗

is asymptotically stable. Thus, it follows from inputto-state stability theory [28] that limt→∞ p̃ei (t) = 0
and limt→∞ ṽie (t) = 0, ∀i ∈ V given the fact that
limt→∞ ~(t) = 0, where ~ = [~T1 , ~T2 ]T .

Lηi =

(42)

This implies that L̇ηi is negative outside the set

Zi = η̄i ∈ R6 | kη̄i k ≤

1
kη̇i k2 .
2kη


c2
ρ(kεi (t1 )k, t − t1 ) .
c1 (χη )2

Thus, η̄i is bounded and ultimately converges to the
set Zi . In view of limt→∞ ρ(kεi (t1 )k, t − t1 ) = 0, it
9

100
80

←80 sec

pz (m)

60

←60 sec

p1
p2
p3
p4
pr

←40 sec

40
20

←20 sec
0

←0 sec

−20
5
0

py (m)

0

−5

10

5

px (m)

Fig. 2. Snapshots of coordinated trajectory tracking of
VTOL UAVs.

Fig. 1. Leader-follower graph Gn+1 .

follows that limt→∞ η̄i (t) = 0, which implies that
limt→∞ ηi (t) = 0 and limt→∞ η̇i (t) = 0, ∀i ∈ V.

pex (m)

2

Since Th2.i)-Th2.iii) have been proven, it can be concluded that the coordinated trajectory tracking of multiple VTOL UAVs is achieved in the sense of (6).

0
−2

pey (m)

−4
20

Simulations

In this section, simulations are performed to verify the
proposed distributed control approach on a formation
of four VTOL UAVs described by (2)-(5). The inertial
parameters are assumed to be identical: mi = 0.85(kg)
and Ji = diag{4.856, 4.856, 8.801} × 10−2 (kgm2 ), i =
1, 2, 3, 4. The leader-follower graph Gn+1 is illustrated
in Fig. 1, where each arrow denotes the corresponding
information flow. Furthermore, define dij = 1 if follower j is accessible to follower i, and dij = 0, otherwise, for i, j = 0, 1, · · · , 4. The desired trajectory is
described as pr (t) = [5 cos(0.2t), 5 sin(0.2t), t]T (m), and
the desired position offsets of the followers relative to
the leader are δ1 = [2, 2, 0]T (m), δ2 = [2, −2, 0]T (m),
δ3 = [−2, −2, 0]T (m) and δ4 = [−2, 2, 0]T (m), respectively. This indicates that the desired formation pattern is a square. The distributed estimator states of
each follower UAV are initialized as zeros. The follower
UAVs are initially situated at p1 (0) = [5, 3, −1]T (m),
p2 (0) = [9, −4, 1]T (m), p3 (0) = [4, −2, −3]T (m) and
p4 (0) = [−1, 4, −2]T (m) and their initial attitudes are
Qi (0) = [1, 0, 0, 0]T , i = 1, 2, 3, 4. The estimator and
control parameters are chosen as follows: kγ = 0.5, kp =
kv = 8 based on (12), la = 12 based on (13), ka = 4
based on (14), lp = lv = kη = 4 and lq = kq = 16. The
simulation results are illustrated in Figs. 2-4.

pe1
5

10

15

20

25

30

35

40

pe2
0

pe3
−2
10

pez (m)

4

−5

5

10

15

5

10

15

20

25

30

35

40

20

25

30

35

40

0

pe4

−1
−2
−3
0

time (sec)

vex (m/sec)

Fig. 3. Position error of follower VTOL UAVs.
4
2
0

vey (m/sec)

0

vez (m/sec)

v1e

−2
5

10

15

20

25

30

35

40

2

v2e

0

v3e

−2
0
2

5

10

15

20

25

30

35

40

v4e

1
0
−1
0

5

10

15

20

time (sec)

25

30

35

40

Fig. 4. Velocity error of follower VTOL UAVs.

that each error component converges to zero asymptotically. Consequently, the simulation results validate
that the proposed distributed control approach effectively guarantees the coordinated trajectory tracking of
multiple VTOL UAVs in the sense of (6).

Fig. 2 exhibits the evolution of the VTOL UAV formation with respect to the leader in the three-dimensional
space, where the formation is depicted every twenty seconds. It can be seen that the follower UAVs reach the prescribed square pattern while tracking the leader. Figs. 3
and 4 describe the position and velocity errors of the follower UAVs with respect to the leader. It can be observed

5

Conclusion

A distributed control strategy is proposed in this paper
to achieve the coordinated trajectory tracking of multiple VTOL UAVs with local information exchange. The

10

connectivity of the network graph is weak in the sense
that we only require the graph to contain a directed
spanning tree. A novel distributed estimator is firstly
designed for each VTOL UAV to obtain the leader’s desired information asymptotically. Then, under the hierarchical framework, a command force and an applied
torque are exploited for each VTOL UAV to fulfill the
accurate tracking to the desired information asymptotically. In addition, an auxiliary system is introduced in
the control development to avoid the non-singular command attitude extraction and the use of the unavailable
desired information. Simulations are carried out to validate the theoretical results.

[14] Z. Li, X. Liu, W. Ren, and L. Xie, “Distributed tracking
control for linear multiagent systems with a leader of bounded
unknown input,” IEEE Transactions on Automatic Control,
vol. 58, no. 2, pp. 518–523, 2013.

References

[18] J. Ghommam, L. F. Luque-Vega, B. Castillo-Toledo, and
M. Saad, “Three-dimensional distributed tracking control
for multiple quadrotor helicopters,” Journal of The Franklin
Institute, vol. 353, no. 10, pp. 2344–2372, 2016.

[15] W. Yu, G. Chen, and M. Cao, “Distributed leader-follower
flocking control for multi-agent dynamical systems with timevarying velocities,” Systems & Control Letters, vol. 59, no. 9,
pp. 543–552, 2010.
[16] J. Qin, C. Yu, and B. Anderson, “On leaderless and leaderfollowing consensus for interacting clusters of second-order
multi-agent systems,” Automatica, vol. 74, pp. 214–221, 2016.
[17] T. Yang, Z. Meng, G. Shi, Y. Hong, and K. H. Johansson,
“Network synchronization with nonlinear dynamics and
switching interactions,” IEEE Transactions on Automatic
Control, vol. 61, no. 10, pp. 3103–3108, 2016.

[1] F. Giulietti, L. Pollini, and M. Innocenti, “Autonomous
formation flight,” IEEE Control Systems Magazine, vol. 20,
no. 6, pp. 34–44, 2000.

[19] F. Liao, R. Teo, J. Wang, X. Dong, F. Lin, and K. Peng,
“Distributed formation and reconfiguration control of vtol
uavs,” IEEE Transactions on Control Systems Technology,
vol. 25, no. 1, pp. 270–277, 2017.

[2] M. Hua, T. Hamel, P. Morin, and C. Samson, “Introduction
to feedback control of underactuated vtol vehicles: A review
of basic control design ideas and principles,” IEEE Control
Systems Magazine, vol. 22, no. 1, pp. 61–75, 2013.

[20] A. Mahnmood and Y. Kim, “Leader-following formation
control of quadcopters with heading synchronization,”
Aerospace Science and Technology, vol. 47, pp. 68–74, 2015.

[3] Z. Zuo, “Trajectory tracking control design with commandfiltered compensation for a quadrotor,” IET Control Theory
and Applications, vol. 4, no. 11, pp. 2343–2355, 2010.

[21] X. Dong, Y. Zhou, Z. Ren, and Y. Zhong, “Timevarying formation tracking for second-order multi-agent
systems subjected to switching topologies with application
to quadrotor formation flying,” IEEE Transactions on
Industrial Electronics, vol. DOI: 10.1109/TIE.2016.2593656,
2016.

[4] Y. Zhang and Y. Tian, “Consentability and protocol design
of multi-agent systems with stochastic switching topology,”
Automatica, vol. 45, no. 5, pp. 1195–1201, 2009.
[5] A. Abdessameud and A. Tayebi, “Formation control of vtol
uavs,” in 48th IEEE Conference on Decision and Control,
Shanghai, P. R. China, 2009, pp. 3454–3459.

[22] M. Shuster, “A survey of attitude representations,” The
Journal of the Astronautical Sciences, vol. 41, no. 4, pp. 439–
517, 1993.

[6] A. Abdessameud, I. G. Polushin, and A. Tayebi, “Motion
coordination of thrust-propelled underactuated vehicles with
intermittent and delayed communications,” Systems &
Control Letters, vol. 79, pp. 15–22, 2015.

[23] Z. Qu, Cooperative control of dynamical systems: applications
to autonomous vehicles. Berlin: Springer, 2009.

[7] X. Dong, B. Yu, Z. Shi, and Y. Zhong, “Time-varying
formation control for unmanned aerial vehicles: Theories
and experiments,” IEEE Transactions on Control Systems
Technology, vol. 23, no. 1, pp. 340–348, 2015.

[24] B. Paden and S. Sastry, “A calculus for computing
filipovs differential inclusion with application to the variable
structure control of robot manipulators,” IEEE Transactions
on Circuits and Systems, vol. CAS-34, no. 1, pp. 73–82, 1987.

[8] Y. Hong, G. Chen, and L. Bushnell, “Distributed observers
design for leader-following control of multi-agent networks,”
Automatica, vol. 44, no. 3, pp. 846–850, 2008.

[25] D. Shevitz and B. Paden, “Lyapunov stability theory of
nonsmooth systems,” IEEE Transactions on Automatic
Control, vol. 39, no. 9, pp. 1910–1914, 1994.

[9] B. Yun, B. Chen, K. Lum, and T. Lee, “Design and
implementation of a leader-follower cooperative control
system for unmanned helicopters,” Journal of Control Theory
and Applications, vol. 8, no. 1, pp. 61–68, 2010.

[26] F. Clarke, Optimization and nonsmooth analysis. New York:
Wiley, 1983.
[27] N. Fischer, R. Kamalapurkar, and W. Dixon, “Lasalleyoshizawa corollaries for nonsmooth systems,” IEEE
Transactions on Automatic Control, vol. 58, no. 9, pp. 2333–
2338, 2013.

[10] D. Mercado, R. Castro, and R. Lozano, “Quadrotors flight
formation control using a leader-follower,” in European
Control Conference, Zurich, Switzerland, 2013, pp. 3858–
3863.

[28] H. Khalil, Nonlinear Systems, Third Edition.
Prentice-Hall, 2002.

[11] D. Lee, “Distributed backstepping control of multiple thrustpropelled vehicles on a balanced graph,” Automatica, vol. 48,
no. 11, pp. 2971–2977, 2012.

New Jersey:

[29] Y. Zou, “Nonlinear robust adaptive hierarchical sliding mode
control approach for quadrotors,” International Journal of
Robust and Nonlinear Control, vol. DOI: 10.1002/rnc.3607,
2016.

[12] A. Loria, J. Dasdemir, and N. A. Jarquin, “Leader-follower
formation and tracking control of mobile robots along straight
paths,” IEEE Transactions on Control Systems Technology,
vol. 24, no. 2, pp. 727–732, 2016.
[13] G. Wen, Y. Zhao, Z. Duan, and W. Yu, “Containment of
higher-order mmlti-leader multi-agent systems: A dynamic
output approach,” IEEE Transactions on Automatic Control,
vol. 61, no. 4, pp. 1135–1140, 2016.

11

