# Modeling and Flocking Consensus Analysis for Large-Scale UAV Swarms

Authors: Li Bing; Li Jie; Huang KeWei
Journal: Mathematical Problems in Engineering
 (2013)

---

## Abstract

Recently, distributed coordination control of the unmanned aerial vehicle (UAV) swarms has been a particularly active topic in intelligent system field. In this paper, through understanding the emergent mechanism of the complex system, further research on the flocking and the dynamic characteristic of UAV swarms will be given. Firstly, this paper analyzes the current researches and existent problems of UAV swarms. Afterwards, by the theory of stochastic process and supplemented variables, a differential-integral model is established, converting the system model into Volterra integral equation. The existence and uniqueness of the solution of the system are discussed. Then the flocking control law is given based on artificial potential with system consensus. At last, we analyze the stability of the proposed flocking control algorithm based on the Lyapunov approach and prove that the system in a limited time can converge to the consensus direction of the velocity. Simulation results are provided to verify the conclusion.


---

## Body

## 1. Introduction

UAV is an advanced system with high autonomy for intelligent combat [1]. In the future, UAVs will be used for complex tasks, such as surveillance, reconnaissance, and precision strike missions. Many organizations have foreseen that in the near future, swarms of UAVs will replace single ones for more complicated missions in more uncertain and possibly hostile environments [2]. Therefore, many researchers are studying groups of cooperative UAVs.(A) Related Work on the UAV Swarms Problem. The new challenges imposed by UAV swarms have attracted many researchers. New control mechanisms, application domains, simulation models, and simulation tools have been developed to tackle issues in different aspects of the swarm. Currently, a completely new topic is opening up in the area of UAV swarms performing different missions cooperatively. [3], the method of evolutionary pinning control is applied to UAV swarms successfully. Path planning and routing are investigated in [4, 5], using multiobjective evolutionary algorithm. The path planning problem in three-dimensional environment without any obstacles is addressed in [6, 7] and with only static obstacles in [8]. [9], cooperative searching problem is discussed for the purpose of detecting moving and evading targets in a hazardous environment. A similar cooperative searching problem is also discussed in [10, 11]. Reference [12] investigates the automatic target recognition (ATR) problem in UAV control and proposes a distributed strategy for UAV swarms. Task allocation problem is discussed in [13–18] using different methodologies. Some applications of using a UAV swarm to search and destroy targets could be found in [19–22]. [23], a swarm simulator for target searching is implemented with Java. Garcia introduces a multi-UAVs simulator implemented with X-Plane—a commercial flight simulator [24]. Russell et al. present a parallel swarm simulation environment which utilizes an existing parallel emulation and simulation tool called SPEEDS [25]. MASON [26] is a general purpose multiagents simulation library utilized in our previous work, along with MATLAB based UAV simulator [27].(B) Related Work on the Consensus Flocking Problem. Most research on formation of agent swarms uses distributed techniques by Reynolds’ seminal work on the mobility of flocks [28], which prescribes three fundamental operations for each robot to realize distributed flocking—separation, alignment, and cohesion [29, 30]. One of the earliest attempts to realize flocking through a set of basic behaviors including safe wandering, aggregation, dispersion, and homing to implement flocking is by Mataric [31]. Kelley and Keating realize flocking with robots based on leader-following behavior [32]. Hayes and Dormiani-Tabatabaei [33] propose a flocking algorithm based on two behaviors: collision avoidance and velocity-matching flock centering. Holland et al. [34] propose a flocking algorithm for UAV similar to Reynolds’. A host is used as an intermediate station for receiving each UAV’s range, bearing, and velocity and sending them to other UAVs to simulate the sensing process of one UAV for perceiving range, bearing, and heading of its neighbors. Ferrante et al. [35] introduce a new communication strategy called the information aware communication for alignment behavior. Recently, Stranieri et al. [36] perform flocking with a swarm of behaviorally heterogeneous mobile robots.In this paper, we consider models for flocking swarms. Firstly, a mathematical model of cooperative system is established by using Markov stochastic process and calculus analysis. Then, the control law for UAV swarm is established based on artificial potential field. At last, we analyze the stability of the proposed flocking control algorithm based on the Lyapunov approach and prove the conclusion that the system in a limited time can converge to the consensus direction of the velocity. Simulation results are provided to verify the conclusion.


## 2. The Model of the UAVs Swarms

### 2.1. Differential Integral Model

LetC(t) denote the state of the UAV swarms at time t;C(t)=0 identifies the state that UAV swarms are stable at time t. The state of UAV i at time t is denoted by ci(t)=(pi(t),oi(t)), in which the first element pi(t)=(xi(t),yi(t),zi(t)) is the UAV’s position in the environment at time t, and the second element oi(t)  is the UAV’s orientation. The UAV’s dynamics is subject to its physical curvature radius constraints, leading to the fact that it can only change its orientation by at most one step, which is described as go straight, go up, go down, turn left, turn upper left, turn lower left, turn right, turn upper right, and turn lower right.In order to obtain Markov random process, the new state of process is derived by supplement of variable [37, 38], which is described as follows:
(1)Pi(t)=P(C(t)=i),i=0,1,2,…,N,Pjk(x,t)=P(C(t)=j,x<yj,k(t)<x+dx),cj=N+1,…,M,k=0,1,2,…,N,
where yi(t) is the dwell time after state i. So it is easy to verify that {U(t),0≤t≤T0}=(C(t),y(t)) is a broad Markov random process.The probability of state transition afterΔt can be obtained using total probability theorem:
(2)P0(t+Δt)=(∑i=0MP(C(t)=i,C(t+Δt)=i))=P0(t)(1-∑i=1Nλ0iΔt+o(Δt))+∑i=1Nλi0Pk(t)Δt+∑j=N+1M∫0TPj(x,t)rj(x)Δtdx+o(Δt).
According to (2) we can get the all probability:
(3)Pk(t+Δt)=Pk(t)(1-∑i=1,i≠kNλkiΔt+o(Δt))+∑i=1,i≠kNλikPk(t)Δt+o(Δt),+∑i=1,i≠kNλliikPkk=1,2,…,N,Pj(x,t)=P(C(t)=j;=P(lthetimeofself-organisedCsystemisx,=P(lC(t+Δt)=j)=Pj(x,t)(1-rj(x)Δt+o(Δt)),=Pj(x,t)1hhhj=N+1,…,M,
where λij is the average sustained rate of each state and rj(x) is the average repair rate at state j. Similarly, the expression of state transition rate for Pj(x+Δx,t+Δt) can be derivated.Differentiate the expression for state transition probability to derive its limit. Then the mathematical model can be described using integral-differential equations as follows:(4)dP0(t)dt+∑i=0,i≠kNλ0kP0(t)=∑i=0,i≠kNλi0Pi(t)+∑j=N+1M∫0Trj(x)pj(x,t)dxdPk(t)dt+∑i=0,i≠kNλkiPk(t)=∑i=0,i≠kNλikPi(t),dPkdtk=1,2,…,N,∂pj(x,t)∂x+∂pj(x,t)∂t+rj(x)pj(x,t)=0,cccj=N+1,…,M.
The boundary and initial conditions are
(5)P0(0)=1,P1(0)=0,…,Pi(0)=0,…,PN(0)=0,Pj(0,t)=∑i=0NλijPi(t),j=N+1,…,M.Theorem 1.
The reliability of coordination system has uniqueness and nonnegative solution onC[0,T].Proof.
According to the initial conditions we can get the analytic solution of the partial differential equation [39, 40].
Set(6)a0=∑i=1Nλ0i,…,aj=∑i=1Nλji,…,aN=∑i=1NλNi,Pj(x,t)=Pi(0,t-x)exp(-∫0xri(μ)dμ),Pj(x,t)=Pi(0,t-x)expi=N+1,…,MdP0(t)dt=-a0P0(t)+∑i=0,i≠kNλi0Pi(t)+∑j=N+1M∫0TPj(0,t-x)exp(-∫0xri(μ)dμ)hhhhhhh×ri(x)dx,P0(t)=exp(-a0t)+∑i=1Nλi0∫0T0Pi(s)exp(-a0(t-s))ds+∫0T0[∑i=N+1M∫0sPi(0,τ)exp(-∫0s-τri(μ)dμ)hhhhhhhhhhhhhh×ri(s-τ)dτ∑i=N+1M∫0sPi(0,τ)exp(-∫0s-τri(μ)dμ)]×exp(-a0(t-s))ds=exp(-a0t)+∑i=1Nλi0∫0T0Pi(s)exp(-a0(t-s))×ds+∑i=N+1M∫0T0Pi(0,τ)dτ×∫0T0-τexp(-a0(T0-τ)∫0vri(μ)dμ+a0v200×∫0T0-τexp-∫0vri(μ)dμ)ri(v)dv=exp(-a0t)+∑i=1Nλi0∫0T0Pi(s)k0ds+∑i=N+1M∫0T0Pi(0,τ)Ki(t-τ)dτ,
where
(7)k0=exp(-a0(T-s))ki(T-τ)=∫0T-τk0*exp(a0v-∫0vri(μ)dμ)ri(v)dv,hhhhhhhhhhhhhhhhhhhhhhhi=N+1,…,M,Pj(t)=∑i=0,i≠jN∫0Texp(-aj(T-s))λijPi(s)ds=∑i=0,i≠jN∫0TkjλijPi(s)ds,j=1,2,…,N,
where kj=exp(-aj(T-s)),j=1,2,…,N.
So we can get the following equation:(8)Pj(0,t)=∑i=0Nλij∑l=0,l≠jN∫0TkjλljPl(s)ds=∫0Tkj(∑i=0Nλij)(∑l=0,l≠jNλljPl(s))ds,hhhhhhhhhhhhhhhhj=N+1,…,M.
Assuming
(9)P(t)=(P0(t),P1(t),…,PN(t),PN+1(0,t),…,PM(0,t)),f(t)=(f0(t),f1(t),…,fN(t),fN+1(0,t),…,fM(0,t))=(exp(-a0t),0,…,0,…,0),
then, the solution of the system can be converted into vectors format as follows:(10)
P
(
t
)
=
f
(
t
)
+
∫
0
T
k
(
T
-
s
)
P
(
s
)
d
s
,
k
(
T
-
τ
)
=
(
0
λ
10
k
0
⋯
λ
N
0
k
0
k
N
+
1
(
T
-
s
)
⋯
k
M
(
T
-
s
)
λ
01
k
1
0
⋯
λ
N
1
k
1
0
⋯
0
⋮
⋮
⋯
⋮
⋮
⋯
⋮
λ
0
N
k
N
λ
2
N
k
N
⋯
0
0
⋯
0
(
∑
l
=
0
,
l
≠
j
N
λ
l
j
)
λ
0
N
+
1
(
∑
l
=
0
,
l
≠
j
N
λ
l
j
)
λ
1
N
+
1
⋯
(
∑
l
=
0
,
l
≠
j
N
λ
l
j
)
λ
N
N
+
1
0
⋯
0
⋮
⋮
⋯
⋮
⋮
⋯
⋮
(
∑
l
=
0
,
l
≠
j
N
λ
l
j
)
λ
0
M
(
∑
l
=
0
,
l
≠
j
N
λ
l
j
)
λ
1
M
⋯
(
∑
l
=
0
,
l
≠
j
N
λ
l
j
)
λ
N
M
0
⋯
0
)
.Any component ofF(t) and G(t-h) vector is nonnegative. The functions Fj(t) and Gj(t-h),j=1,2,…,Mare limitary on the domain 0<T<+∞. The solution of integral equation is unique and nonnegative on C[0,T]. So the reliability of coordination system has unique and nonnegative solution on C[0,T].


### 2.2. Probabilistic Analysis Based on State Transformation

The behavior evolution of the UAV swarm system is a limited Markov decision process. Suppose that the probability distribution of the system state isP(N,t)at time t. Then at time t+τ, the probability distribution is P(N,t+τ). According to the relationship of the probability density at different time, the marginal probability density P(N,t+τ) is P(N,t+τ)=∫P(η,t)P(N,t+τ∣η,t)dη.And the time derivative of theP(N,t+τ) is
(11)dP(N,t)dt=limτ→0P(N,t+τ)-P(N,t)τ.
Define Wt(η,N) as the transition probability density from state η to state N in unit time during time interval [t,t+τ]. So the transition probability from state η to state N during time interval [t,t+τ] is τWt(η,N). Then the probability by which the transition does not happen is
(12)1-τ∫Wt(η,N)δ(N-η)dn,
where δ(N-η)=P(N,t+τ∣η,t), when τ=0. Thus,
(13)dP(N,t)dt=∫{W(η,n)P(η,t)-W(n,η)P(n,t)}dη.
Equation (13) describes the evolution of the system states over time, which is the primary equation model of the UAV swarms behavior.


## 2.1. Differential Integral Model

LetC(t) denote the state of the UAV swarms at time t;C(t)=0 identifies the state that UAV swarms are stable at time t. The state of UAV i at time t is denoted by ci(t)=(pi(t),oi(t)), in which the first element pi(t)=(xi(t),yi(t),zi(t)) is the UAV’s position in the environment at time t, and the second element oi(t)  is the UAV’s orientation. The UAV’s dynamics is subject to its physical curvature radius constraints, leading to the fact that it can only change its orientation by at most one step, which is described as go straight, go up, go down, turn left, turn upper left, turn lower left, turn right, turn upper right, and turn lower right.In order to obtain Markov random process, the new state of process is derived by supplement of variable [37, 38], which is described as follows:
(1)Pi(t)=P(C(t)=i),i=0,1,2,…,N,Pjk(x,t)=P(C(t)=j,x<yj,k(t)<x+dx),cj=N+1,…,M,k=0,1,2,…,N,
where yi(t) is the dwell time after state i. So it is easy to verify that {U(t),0≤t≤T0}=(C(t),y(t)) is a broad Markov random process.The probability of state transition afterΔt can be obtained using total probability theorem:
(2)P0(t+Δt)=(∑i=0MP(C(t)=i,C(t+Δt)=i))=P0(t)(1-∑i=1Nλ0iΔt+o(Δt))+∑i=1Nλi0Pk(t)Δt+∑j=N+1M∫0TPj(x,t)rj(x)Δtdx+o(Δt).
According to (2) we can get the all probability:
(3)Pk(t+Δt)=Pk(t)(1-∑i=1,i≠kNλkiΔt+o(Δt))+∑i=1,i≠kNλikPk(t)Δt+o(Δt),+∑i=1,i≠kNλliikPkk=1,2,…,N,Pj(x,t)=P(C(t)=j;=P(lthetimeofself-organisedCsystemisx,=P(lC(t+Δt)=j)=Pj(x,t)(1-rj(x)Δt+o(Δt)),=Pj(x,t)1hhhj=N+1,…,M,
where λij is the average sustained rate of each state and rj(x) is the average repair rate at state j. Similarly, the expression of state transition rate for Pj(x+Δx,t+Δt) can be derivated.Differentiate the expression for state transition probability to derive its limit. Then the mathematical model can be described using integral-differential equations as follows:(4)dP0(t)dt+∑i=0,i≠kNλ0kP0(t)=∑i=0,i≠kNλi0Pi(t)+∑j=N+1M∫0Trj(x)pj(x,t)dxdPk(t)dt+∑i=0,i≠kNλkiPk(t)=∑i=0,i≠kNλikPi(t),dPkdtk=1,2,…,N,∂pj(x,t)∂x+∂pj(x,t)∂t+rj(x)pj(x,t)=0,cccj=N+1,…,M.
The boundary and initial conditions are
(5)P0(0)=1,P1(0)=0,…,Pi(0)=0,…,PN(0)=0,Pj(0,t)=∑i=0NλijPi(t),j=N+1,…,M.Theorem 1.
The reliability of coordination system has uniqueness and nonnegative solution onC[0,T].Proof.
According to the initial conditions we can get the analytic solution of the partial differential equation [39, 40].
Set(6)a0=∑i=1Nλ0i,…,aj=∑i=1Nλji,…,aN=∑i=1NλNi,Pj(x,t)=Pi(0,t-x)exp(-∫0xri(μ)dμ),Pj(x,t)=Pi(0,t-x)expi=N+1,…,MdP0(t)dt=-a0P0(t)+∑i=0,i≠kNλi0Pi(t)+∑j=N+1M∫0TPj(0,t-x)exp(-∫0xri(μ)dμ)hhhhhhh×ri(x)dx,P0(t)=exp(-a0t)+∑i=1Nλi0∫0T0Pi(s)exp(-a0(t-s))ds+∫0T0[∑i=N+1M∫0sPi(0,τ)exp(-∫0s-τri(μ)dμ)hhhhhhhhhhhhhh×ri(s-τ)dτ∑i=N+1M∫0sPi(0,τ)exp(-∫0s-τri(μ)dμ)]×exp(-a0(t-s))ds=exp(-a0t)+∑i=1Nλi0∫0T0Pi(s)exp(-a0(t-s))×ds+∑i=N+1M∫0T0Pi(0,τ)dτ×∫0T0-τexp(-a0(T0-τ)∫0vri(μ)dμ+a0v200×∫0T0-τexp-∫0vri(μ)dμ)ri(v)dv=exp(-a0t)+∑i=1Nλi0∫0T0Pi(s)k0ds+∑i=N+1M∫0T0Pi(0,τ)Ki(t-τ)dτ,
where
(7)k0=exp(-a0(T-s))ki(T-τ)=∫0T-τk0*exp(a0v-∫0vri(μ)dμ)ri(v)dv,hhhhhhhhhhhhhhhhhhhhhhhi=N+1,…,M,Pj(t)=∑i=0,i≠jN∫0Texp(-aj(T-s))λijPi(s)ds=∑i=0,i≠jN∫0TkjλijPi(s)ds,j=1,2,…,N,
where kj=exp(-aj(T-s)),j=1,2,…,N.
So we can get the following equation:(8)Pj(0,t)=∑i=0Nλij∑l=0,l≠jN∫0TkjλljPl(s)ds=∫0Tkj(∑i=0Nλij)(∑l=0,l≠jNλljPl(s))ds,hhhhhhhhhhhhhhhhj=N+1,…,M.
Assuming
(9)P(t)=(P0(t),P1(t),…,PN(t),PN+1(0,t),…,PM(0,t)),f(t)=(f0(t),f1(t),…,fN(t),fN+1(0,t),…,fM(0,t))=(exp(-a0t),0,…,0,…,0),
then, the solution of the system can be converted into vectors format as follows:(10)
P
(
t
)
=
f
(
t
)
+
∫
0
T
k
(
T
-
s
)
P
(
s
)
d
s
,
k
(
T
-
τ
)
=
(
0
λ
10
k
0
⋯
λ
N
0
k
0
k
N
+
1
(
T
-
s
)
⋯
k
M
(
T
-
s
)
λ
01
k
1
0
⋯
λ
N
1
k
1
0
⋯
0
⋮
⋮
⋯
⋮
⋮
⋯
⋮
λ
0
N
k
N
λ
2
N
k
N
⋯
0
0
⋯
0
(
∑
l
=
0
,
l
≠
j
N
λ
l
j
)
λ
0
N
+
1
(
∑
l
=
0
,
l
≠
j
N
λ
l
j
)
λ
1
N
+
1
⋯
(
∑
l
=
0
,
l
≠
j
N
λ
l
j
)
λ
N
N
+
1
0
⋯
0
⋮
⋮
⋯
⋮
⋮
⋯
⋮
(
∑
l
=
0
,
l
≠
j
N
λ
l
j
)
λ
0
M
(
∑
l
=
0
,
l
≠
j
N
λ
l
j
)
λ
1
M
⋯
(
∑
l
=
0
,
l
≠
j
N
λ
l
j
)
λ
N
M
0
⋯
0
)
.Any component ofF(t) and G(t-h) vector is nonnegative. The functions Fj(t) and Gj(t-h),j=1,2,…,Mare limitary on the domain 0<T<+∞. The solution of integral equation is unique and nonnegative on C[0,T]. So the reliability of coordination system has unique and nonnegative solution on C[0,T].


## 2.2. Probabilistic Analysis Based on State Transformation

The behavior evolution of the UAV swarm system is a limited Markov decision process. Suppose that the probability distribution of the system state isP(N,t)at time t. Then at time t+τ, the probability distribution is P(N,t+τ). According to the relationship of the probability density at different time, the marginal probability density P(N,t+τ) is P(N,t+τ)=∫P(η,t)P(N,t+τ∣η,t)dη.And the time derivative of theP(N,t+τ) is
(11)dP(N,t)dt=limτ→0P(N,t+τ)-P(N,t)τ.
Define Wt(η,N) as the transition probability density from state η to state N in unit time during time interval [t,t+τ]. So the transition probability from state η to state N during time interval [t,t+τ] is τWt(η,N). Then the probability by which the transition does not happen is
(12)1-τ∫Wt(η,N)δ(N-η)dn,
where δ(N-η)=P(N,t+τ∣η,t), when τ=0. Thus,
(13)dP(N,t)dt=∫{W(η,n)P(η,t)-W(n,η)P(n,t)}dη.
Equation (13) describes the evolution of the system states over time, which is the primary equation model of the UAV swarms behavior.


## 3. Flocking Control of UAV Swarms

### 3.1. Flocking Control Law

In this section, first we design a distributed flocking control law. Assuming that each UAV senses its own position and velocity and is able to obtain its neighbors’ position and velocity, the UAV swarms form flocking behaviour model structure control law as follows:(14)Ui=Uniformity(v)∑j=1,j≠iNf(pi,j)(pi-pj)+β(pi-pgoal)+∑j=1,j≠iNai,j(pi,j)(vj-vi),
where Uniformity(v)=(α∑i≠j∥vi-vj∥2). Note that alignment at a common velocity is equivalent to Uniformity(v)=0.Pi,j is the distance between the individual i and j.  Uij is potential function and satisfies the condition [29, 30]:(i)
UiR(∥pi-pj∥)→∞,∥pi-pj∥<Rmin,(ii)
Rmin≤∥pi-pj∥≤Rmax,∃!minUiR.


### 3.2. Stability Analysis

Consider the following positive semidefinite function:(15)E=12(∑i=1N(∑j=1,j≠iNUi,j+KiUi,goal)+viTvi).
In order to facilitate writing, we simplify the certification process variable substitution as follows:
(16)U-i=(∑j=1,j≠iNUi,j+KiUi,goal),E˙=12∑i=1NU-˙i-∑i=1NviTUi=12∑i=1NU-˙i+∑i=1NviT(-∑j=1,j≠iN∇piUij(pi,j)-Ki∇piUiA(pi,j)hhhhhhhhhhhhhhhh-∑j=1,j≠iNαi,j(vi-vj))=∑i=1NviT∇piU-i+∑i=1NviT(-∑j=1,j≠iN∇piUij(pi,j)hhhhhhhhhhhhhhhhhhh-Ki∇piUi,goal(pi,goal)hhhhhhhhhhhhhhhhhhh-∑j=1,j≠iNαi,j(vi-vj))=∑i=1NviT∇pi(∑j=1,j≠iNUi,j+KiUi,goal)+∑i=1NviT(-∑j=1,j≠iN∇piUij-Ki∇piUi,goal+∑i=1NviThh-∑j=1,j≠iNαi,j(vi-vj))=∑i=1NviT(-∑j=1,j≠iNαi,j(vi-vj))=-vT(LC⊗I3)v,
where Lc is UAV swarms system satisfying the Laplacian matrix of the communication conditions. Therefore, the quadratic form is explicitly described as follows:
(17)E˙=-vxTLcvx-vyTLcvy-vzTLcvz≤0.
Consider the following collections: {vi,Pi,j∣E≤C} is a closed set. The following is to verify that it is a compact set, and there is a clear conclusion that Pi,j≤C. Similarly viTvi≤C,∥vi∥≤C1/2, and according to the definition of the potential field we obtain ∥Pi,j∥≤Ui,j-1(C(N-1)). According to the LaSalle invariance principle, the system will converge to the largest invariant set in the area and meet E˙=0. According to E˙=0, when the system enters the steady state, the speed of each individual is equal, and all individuals move to the target position Pgoal, making the overall potential energy minimum.Theorem 2.
Consider the UAV swarms consisting ofN UAVs. The position of individual i is pi. All individuals in the swarms will eventually build up to the spherical region:
(18)∥pi-pc∥≤2λmax(G)(N-1)lCλmin(ATG+GA).Proof.
Consider(19)Pc˙=VcVc˙=1N∑i=1Nv˙i=1N∑i=1N(-∑j=1,j≠iN∇piUij(∥pi-pj∥))-βi∇piUiA(∥pi-pc∥)-∑j=1,j≠iNαi,j(vi-vj)=0,
where Pc(t)=(1)/N∑i=1Npi(t),Vc=(1)/N∑i=1Nαi,jvi.
By making the variable replacementεip=Pi-Pc,εip=vi-vc, we get
(20)ε˙ip=εiv,ε˙ip=(-∑j=1,j≠iN∇piUij(∥pi-pj∥)-βi∇piUiA(∥pi-pc∥)∑j=1,j≠iNαi,j(vi-vj))=-βiεip-N(∑j=1,j≠iNαi,j)εiv-∑j=1,j≠iN∇piUij(∥pi-pj∥).
Then
(21)ζi˙=(ε˙ipε˙iv)=(01-βi-N(∑j=1,j≠iNαi,j))ζi-(0∑j=1,j≠iN∇piUij(∥pi-pj∥)),A=(01-βi-N(∑j=1,j≠iNαi,j)),B=(0∑j=1,j≠iN∇piUij(∥pi-pj∥)).
Since Ki,N>0,
(22)λ(A)=-N(∑j=1,j≠iNαi,j)2±(N(∑j=1,j≠iNαi,j))2-4ki2,-(∑j=1,j≠iNαi,j)Re(λ(A))<0.
The largest and the smallest eigenvalues of symmetric positive definite matrix S are λmax(S) and λmin(S), respectively. The symmetric positive definite matrixSwith appropriate dimensions satisfies the following conclusion [41, 42].
(23)λmin(S)∥ζi∥2≤ζiTSζi≤λmax(S)∥ζi∥2.
Finally, select Lyapunov function(24)Ei=ζiTGζi,GT=G,λ(G)>0.
Time derivative can be obtained:(25)Ei˙=ζiT˙Gζi+ζiTGζi˙=(Aζi-B)TGζi+ζiTG(Aζi-B)=-ζiT(ATG+GA)ζi-2ζiTGB,Uij(∥pi-pj∥)∥pi-pj∥≤lC.
Therefore, according to the above formula we obtain
(26)Ei˙≤-λmin(ATG+GA)∥ζi∥2+2λmax(G)∥ζi∥(N-1)lC.
When
(27)∥ζi∥>2λmax(G)(N-1)lCλmin(ATG+GA),E˙i<0. The system continues to move closer to the population centre. Therefore, eventually the system stabilizes at a known system of
(28)∥ζi∥≤2λmax(G)(N-1)lCλmin(ATG+GA).


## 3.1. Flocking Control Law

In this section, first we design a distributed flocking control law. Assuming that each UAV senses its own position and velocity and is able to obtain its neighbors’ position and velocity, the UAV swarms form flocking behaviour model structure control law as follows:(14)Ui=Uniformity(v)∑j=1,j≠iNf(pi,j)(pi-pj)+β(pi-pgoal)+∑j=1,j≠iNai,j(pi,j)(vj-vi),
where Uniformity(v)=(α∑i≠j∥vi-vj∥2). Note that alignment at a common velocity is equivalent to Uniformity(v)=0.Pi,j is the distance between the individual i and j.  Uij is potential function and satisfies the condition [29, 30]:(i)
UiR(∥pi-pj∥)→∞,∥pi-pj∥<Rmin,(ii)
Rmin≤∥pi-pj∥≤Rmax,∃!minUiR.


## 3.2. Stability Analysis

Consider the following positive semidefinite function:(15)E=12(∑i=1N(∑j=1,j≠iNUi,j+KiUi,goal)+viTvi).
In order to facilitate writing, we simplify the certification process variable substitution as follows:
(16)U-i=(∑j=1,j≠iNUi,j+KiUi,goal),E˙=12∑i=1NU-˙i-∑i=1NviTUi=12∑i=1NU-˙i+∑i=1NviT(-∑j=1,j≠iN∇piUij(pi,j)-Ki∇piUiA(pi,j)hhhhhhhhhhhhhhhh-∑j=1,j≠iNαi,j(vi-vj))=∑i=1NviT∇piU-i+∑i=1NviT(-∑j=1,j≠iN∇piUij(pi,j)hhhhhhhhhhhhhhhhhhh-Ki∇piUi,goal(pi,goal)hhhhhhhhhhhhhhhhhhh-∑j=1,j≠iNαi,j(vi-vj))=∑i=1NviT∇pi(∑j=1,j≠iNUi,j+KiUi,goal)+∑i=1NviT(-∑j=1,j≠iN∇piUij-Ki∇piUi,goal+∑i=1NviThh-∑j=1,j≠iNαi,j(vi-vj))=∑i=1NviT(-∑j=1,j≠iNαi,j(vi-vj))=-vT(LC⊗I3)v,
where Lc is UAV swarms system satisfying the Laplacian matrix of the communication conditions. Therefore, the quadratic form is explicitly described as follows:
(17)E˙=-vxTLcvx-vyTLcvy-vzTLcvz≤0.
Consider the following collections: {vi,Pi,j∣E≤C} is a closed set. The following is to verify that it is a compact set, and there is a clear conclusion that Pi,j≤C. Similarly viTvi≤C,∥vi∥≤C1/2, and according to the definition of the potential field we obtain ∥Pi,j∥≤Ui,j-1(C(N-1)). According to the LaSalle invariance principle, the system will converge to the largest invariant set in the area and meet E˙=0. According to E˙=0, when the system enters the steady state, the speed of each individual is equal, and all individuals move to the target position Pgoal, making the overall potential energy minimum.Theorem 2.
Consider the UAV swarms consisting ofN UAVs. The position of individual i is pi. All individuals in the swarms will eventually build up to the spherical region:
(18)∥pi-pc∥≤2λmax(G)(N-1)lCλmin(ATG+GA).Proof.
Consider(19)Pc˙=VcVc˙=1N∑i=1Nv˙i=1N∑i=1N(-∑j=1,j≠iN∇piUij(∥pi-pj∥))-βi∇piUiA(∥pi-pc∥)-∑j=1,j≠iNαi,j(vi-vj)=0,
where Pc(t)=(1)/N∑i=1Npi(t),Vc=(1)/N∑i=1Nαi,jvi.
By making the variable replacementεip=Pi-Pc,εip=vi-vc, we get
(20)ε˙ip=εiv,ε˙ip=(-∑j=1,j≠iN∇piUij(∥pi-pj∥)-βi∇piUiA(∥pi-pc∥)∑j=1,j≠iNαi,j(vi-vj))=-βiεip-N(∑j=1,j≠iNαi,j)εiv-∑j=1,j≠iN∇piUij(∥pi-pj∥).
Then
(21)ζi˙=(ε˙ipε˙iv)=(01-βi-N(∑j=1,j≠iNαi,j))ζi-(0∑j=1,j≠iN∇piUij(∥pi-pj∥)),A=(01-βi-N(∑j=1,j≠iNαi,j)),B=(0∑j=1,j≠iN∇piUij(∥pi-pj∥)).
Since Ki,N>0,
(22)λ(A)=-N(∑j=1,j≠iNαi,j)2±(N(∑j=1,j≠iNαi,j))2-4ki2,-(∑j=1,j≠iNαi,j)Re(λ(A))<0.
The largest and the smallest eigenvalues of symmetric positive definite matrix S are λmax(S) and λmin(S), respectively. The symmetric positive definite matrixSwith appropriate dimensions satisfies the following conclusion [41, 42].
(23)λmin(S)∥ζi∥2≤ζiTSζi≤λmax(S)∥ζi∥2.
Finally, select Lyapunov function(24)Ei=ζiTGζi,GT=G,λ(G)>0.
Time derivative can be obtained:(25)Ei˙=ζiT˙Gζi+ζiTGζi˙=(Aζi-B)TGζi+ζiTG(Aζi-B)=-ζiT(ATG+GA)ζi-2ζiTGB,Uij(∥pi-pj∥)∥pi-pj∥≤lC.
Therefore, according to the above formula we obtain
(26)Ei˙≤-λmin(ATG+GA)∥ζi∥2+2λmax(G)∥ζi∥(N-1)lC.
When
(27)∥ζi∥>2λmax(G)(N-1)lCλmin(ATG+GA),E˙i<0. The system continues to move closer to the population centre. Therefore, eventually the system stabilizes at a known system of
(28)∥ζi∥≤2λmax(G)(N-1)lCλmin(ATG+GA).


## 4. Simulation of System Flocking Formation Behavior

According to the UAV’s physical characteristics, this paper will discretize the time with high frequency. Thus, a UAVi makes its path decision Pi(t+1) at time-step t and will execute an action as the following equation:
(29)Pi(t+1)=Pi(t)+vΔt.The movement of the individual is not only controlled by itself but also affected by the state of other individuals. Therefore, the individual direction of movement at a certain time is not only relative to its direction one moment before, but also relative to the directions of its surrounding individuals’ movements. The influence of all the individuals to the individuali can be described as the following equation:
(30)1∑1/pij∑j=1nθj(t)pij.Then, the speed direction of the UAVi at time (t+1) can be modified as the following equation:
(31)θi(t+1)=αθi(t)+β1∑1/pij∑j=1nθj(t)pij+γarctany0-y1x0-x1.We consider the swarms of 100 UAVs with six degrees of freedom. The weights of the cost function are set toα=0.3,β=0.5,γ=0.2, V=[200,0,0] m/s and m=25kg. Direction is the rand variable from −2 * pi to 2 * pi. The position of the UAVs is the rand variable. The factors of the influence on the flight are wind and airstream. The results for the case of 100 UAVs are shown in Figures 1 and 2.Figure 1
Velocities with respect to time.Figure 2
Trajectories with respect to time.From Figure1, at t=19s, the velocities of the swarms achieve consensus at θ=3.2.Figure2 describes the trajectories with respect to time. The UAV swarms system will eventually be able to form a stable distance between each individual and the same velocity vectors. Collision between individuals is thus avoided.Figures3, 4, and 5 show the Pitch, Roll, and Attack with respect to time. From the simulation results, we can conclude that the UAVs based on the method successfully fly after the adjustment at the initial stage.Figure 3
Angle of the Pitch with respect to time.Figure 4
Angle of the Roll with respect to time.Figure 5
Angle of the Attack with respect to time.Figure6 shows the Sideslip with respect to time. Through the analysis of the Sideslip Angle, we can find that the Angle of the Sideslip is less than 0.5 degrees and tends to zero, to ensure the turning flight control.Figure 6
Angle of the Sideslip with respect to time.


## 5. Conclusion

This paper analyzed current researches and existent problems of UAV swarms. Afterwards, by the theory of stochastic process and supplemented variables, a differential-integral model was established. The existence and uniqueness of the solution of the system were discussed. The flocking control law is given based on artificial potential with system consensus. At last, we analyzed the stability of the proposed flocking control algorithm based on the Lyapunov approach and proved the conclusion that the system in 28 s can converge to the consensus direction of the velocity. And we performed simulation tests to verify the conclusion.



---

*Source: 368369-2013-11-19.xml*