Electronic Journal of Differential Equations, Vol. 2021 (2021), No. 104, pp. 1–24.

ISSN: 1072-6691. URL: http://ejde.math.txstate.edu or http://ejde.math.unt.edu

DYNAMICS OF FLOCKING MODELS WITH TWO SPECIES

QINGJIAN ZHAO, SHAOYUN SHI, WENLEI LI

Abstract. This article studies the flocking behavior of self-organized agents

in two species. First, referring to the work of Olfati-Saber and the classical

Cucker-Smale model, we establish a discrete system describing the flocking dy-
namic of the agents in two species. Second, by using the LaSalle’s invariance

principle, we show that the system with global interaction will achieve uncon-

ditional time-asymptotic flocking, and the system with local interaction has a
time-asymptotic flocking under certain assumptions. Moreover, we investigate

the local asymptotic stability of a class of flocking solutions. Finally, some
numerical simulations and qualitative results are presented.

1. Introduction

The self-organized behavior, illustrating that some agents interact with each
other to reach a position and velocity balance, wildly exist in the natural world and
human society. Scientists have established many models of various phenomena to
describe these behaviors, such as birds flocking [4, 10], fishes swarming [2], complex
network [6] and the price monitor in the market [1, 18]. The corresponding models
have also been applied in the field of systems science and control engineering, such
as the flying control for space flight [26, 27, 34] and the formation control for mobile
robots team [29].

In 1986, Reynolds [28] introduced three heuristic rules showing the basic require-
ments on the collective behavior of self-organized agents, which can be summarized
as follows.

(1) cohesion: attempt to stay close to nearby flockmates;
(2) separation: avoid collisions with nearby flockmates;
(3) alignment: attempt to match velocity with nearby flockmates.

Following the above rules, scientists have established many standard models. One
class of these works is the flocking models, which mainly describe a final state
of alignment under some principles of interaction between agents. In detail, let
xi(t) ∈ Rd and vi(t) ∈ Rd be the position and velocity of the i-th agent in a group
N with n agents, i = 1, . . . , n. Then mathematically, the agents in N achieve
flocking if for any 1 ≤ i, j ≤ n, the following two properties hold:

(1) ‖vj(t)− vi(t)‖ = 0, as t→∞.
(2) there exists a constant C > 0 such that for all t ≥ 0, ‖xj(t)− xi(t)‖ ≤ C.

2010 Mathematics Subject Classification. 34D20, 92D25, 37D10, 65L07.
Key words and phrases. Flocking dynamics; discrete model, LaSalle invariance principle;

invariant manifold.
©2021. This work is licensed under a CC BY 4.0 license.

Submitted August 14, 2021. Published December 30, 2021.

1



2 Q. ZHAO, S. SHI, W. LI EJDE-2021/104

In 2006, Olfati-Saber [24] studied the following flocking model in which the
interactions between agents are defined as the addition of an action function of
distances and a matching function of velocities,

ẋi(t) = vi(t),

ẋi(t) =

n∑
j=1

φ(‖xj(t)− xi(t)‖ − d)~nij + c(vj(t)− vi(t)),
(1.1)

where the interaction function φ is a monotone increasing function with φ(0) = 0,

~nij =
xj(t)−xi(t)
‖xj(t)−xi(t)‖ , and c, d are positive constants called the velocity matching

coefficient and the balance distance, respectively. Under the condition that the
corresponding graph of the agents keeps connected, the author obtained that the
agents achieve flocking.

In 2007, Cucker and Smale [12] considered the model (C-S model):

ẋi(t) = vi(t),

v̇i(t) =

n∑
j=1

aij(x(t))(vj(t)− vi(t)),
(1.2)

with the velocity matching coefficient function

aij =
H

(1 + ‖xj(t)− xi(t)‖2)β
,

where H > 0, β > 0 are parameters. They concluded that the agents reach a final
state as unconditional flocking if β < 1/2, yet for β ≥ 1/2, it is necessary to have
some restrictions on the initial values.

Following the C-S model, scientists have developed various models and gained
many results. In 2008, Ha and Tadmor used the BBGKY hierarchy to turn the
C-S model into a Vlasov type mean-field model and proved that for β = 1

2 , the
system still has unconditional flocking [32]. In 2009, Carrillo extended the C-
S model to a continuous kinetic version model with a Boltzmann-type equation
[3], and proved that the solutions will exponentially concentrate their velocities
to their mean and they will converge towards a translational flocking solution. In
2014, Motsch and Tadmor presented a generalized C-S model in [23] and proved the
flocking of the system under a row-stochastic assumption and some requirements
on the interaction between agents. There are many other interesting results related
to the above models, such as flocking with random influence of white noise [11, 17],
randomly switching topologies [15], stochastic mean-field limit [9], non-symmetric
interaction or leadership [13, 30], collision avoidance of obstacle [5, 35], and time-
delay [7, 14, 21, 25, 31], etc.

One can notice that nearly all the above research works studied only the flocking
behaviors of agents in one species, that is to say, all the agents follow the same inter-
action function and velocity matching coefficient. While, no matter in the natural
world or human society, the interactions and gathering of multi-species agents are
more widespread situations, such as several kinds of birds migrate together, vari-
ous vehicles move together on the road, and peoples live together in a multi-ethnic
society. In these circumstances, it can always be observed that two or more species
of agents with different influence disciplines tend to form a whole cluster and move
in the same way.



EJDE-2021/104 DYNAMICS OF FLOCKING MODELS WITH TWO SPECIES 3

As it was pointed out in [8], the mathematical study of two-species behaviors is
a hot topic in the research area of flocking behaviors. Several continuous models on
two-species flocking have been established and analysed, see [8, 20, 33] for examples.
The main goal of this paper is to establish a discrete model of some self-organized
agents from two species and to investigate the flocking phenomena and related
dynamics of the model. Furthermore, we will present several numerical patterns
of flocking dynamics to illustrate our results and some other phenomena related to
the models.

The outline of this article is as follows: in section 2, we establish the two-species
flocking model. In section 3, the model is considered to have global interaction
functions and is proved to have unconditional flocking. In section 4, we carry out a
local model by adding a restriction on the interaction range. Then we use LaSalle’s
invariance principle and the invariant manifold theory to analyze the stability of
the solutions. Finally, some numerical simulations and qualitative results will be
concluded in section 5.

2. Flocking model for two groups

In this section, we will establish a flocking model of two species based on the
Reynolds three rules and the models (1.1) and (1.2) in section 1. Before the mod-
eling work, let us recall some concepts of graph theory. A graph G is defined as a
pair (V, E) with a set of vertices V = {1, 2, . . . , n} and a set of edges E ⊆ V×V. Let
‖V‖ and ‖E‖ be the numbers of vertices and edges in graph G, respectively, then
it is easy to find that ‖V‖ = n and ‖E‖ ≤ n2. The set of neighbours of vertex i is
defined by

Ni = {j ∈ V : (i, j) ∈ E}.
If for any i, j ∈ V, there has a sequence of vertices {i, v1, . . . , vn, j} such that
v1 ∈ Ni, j ∈ Nvn , and vk+1 ∈ Nvk , k = 1, . . . , n−1, then we say that the graph G is
connected and this sequence is called a walk between i and j. If G is disconnected,
then V can be decompose into a pair of separated sets V1 and V2, where V1∪V2 = V,
V1 ∩ V2 = ∅, and (i, j) /∈ E for any i ∈ V1, j ∈ V2.

Definition 2.1. Let G = (V, E) be a graph. A subset V̂ ⊆ V is called a cluster if

the corresponding graph Ĝ = (V̂, Ê) is a connected graph with Ê = V̂ × V̂ ∩ E .

An n× n matrix A(x) = [aij(x)] can be viewed as a weighted adjacency matrix
of graph G = (V, E), i.e. for any i, j ∈ V, aij(x) 6= 0 if and only if (i, j) ∈ E . The
degree matrix D(x) of G is defined as a diagonal matrix with diagonal elements∑n
i=1 aii(x). Then, matrix L(x) = [lij(x)] induced by

L(x) = D(x)−A(x) =


∑
j 6=1 a1j −a12 . . . −a1n
−a21

∑
j 6=2 a2j . . . −a2n

...
...

. . .
...

−an1 −an2 . . .
∑
j 6=n anj

 , (2.1)

is called the Laplacian matrix of G.
It is easy to find that in the matrix L(x), the elements of each row have a sum

of zero. This fact implies that L(x) has a eigenvalue λ1 = 0 with the related

eigenvector 1n = (1, . . . , 1)
>

. Furthermore, the Laplacian matrix L(x) has the
following properties [16, 22]:



4 Q. ZHAO, S. SHI, W. LI EJDE-2021/104

(1) If A(x) is a non-negative matrix, then the eigenvalues of the corresponding
L(x) satisfy 0 = λ1 ≤ λ2 ≤ · · · ≤ λn.

(2) If A(x) is a non-negative symmetric matrix, then the corresponding L(x)
is a positive semi-definite matrix which satisfies

z>Lz =
1

2

∑
aij 6=0

aij(zj − zi)2, z ∈ Rn. (2.2)

(3) If A(x) is the adjacency matrix of a graph G, then the Fiedler number

λ2 = min
z⊥1n

z>Lz

‖z‖2
> 0,

if and only if the graph G is connected.

In the follow-up post, we consider two groups of agents which are divided ac-
cording to the species. The two groups are indicated by the sets N1, N2 with
N1 ∩ N2 = ∅ and N1 ∪ N2 = N . For the sake of convenience, here we give some
notations which will be used. Let
ni: the number of agents in Ni, i = 1, 2;
n: the total number of all the agents in N , i.e., n = n1 + n2;
xl(t) ∈ Rd: the position of the l-th agent in N1 at time t;
ul(t) ∈ Rd: the velocity of the l-th agent in N1 at time t;
yi(t) ∈ Rd: the position of the i-th agent in N2 at time t;
vi(t) ∈ Rd: the velocity of the i-th agent in N2 at time t.
Therefore, we can use the pairs (xl, ul) ∈ N1 and (yi, vi) ∈ N2 as the l-th agent in
N1 and i-th agent in N2, respectively.

Based on the Newton’s second law, we illustrate the interactions of the agents
in two groups by the kinetic model

ẋl = ul,

u̇l = f1l + g1l , l = 1, . . . , n1,

ẏi = vi,

v̇i = f2i + g2i , i = 1, . . . , n2,

(2.3)

here, for k = 1, 2 and l ∈ {1, . . . , nk}, fkl represents the influence acting on the
l-th agent in Nk from other agents in Nk and gkl represents that from all agents in
N/Nk.

First, according to the cohesion and separation rules, when the agents achieve
flocking, the distance between every two agents will tend to reach a proper value.
This ideal value under the flocking state is called the balance distance. In the model
of two groups, we introduce three different balance distances to show the distinction
between the two groups. In details, for agents in N1 and N2, the balance distances
are respectively denoted by d1, d2, and for agents from different groups, the balance
distance is d0 which should be not less than d1, d2. Hence, we utilize the idea of
Ofati-Saber [24] presented in model (1.1) and give a collective potential function V
as follows

V (x, y) = V1(x) + V2(y) + V0(x, y), (2.4)

with

V1(x) =
1

2

n1∑
k=1

n1∑
l=1

ψ(‖xk − xl‖ − d1),



EJDE-2021/104 DYNAMICS OF FLOCKING MODELS WITH TWO SPECIES 5

V2(y) =
1

2

n2∑
j=1

n2∑
i=1

ψ(‖yj − yi‖ − d2),

V0(x, y) =
1

2

n1∑
k=1

n2∑
i=1

ψ(‖xk − yi‖ − d0),

where x = (x1, . . . , xn1
), y = (y1, . . . , yn2

), and ψ(z) : R → R is a nonnegative
function with its strict minimum at the original point. This property of ψ(z)
reflects the attractive and repulsive behaviors between agents and these behaviors
ensure that the agents can always tend to their balance distances.

Next, to depict the alignment rule, it is natural to apply the velocity matching
idea of the C-S model [12] as it is presented in (1.2). Specifically, we can give the
expression of the influence functions,

f1l =

n1∑
k=1

φ(‖xk − xl‖ − d1)~e(xk, xl) +

n1∑
k=1

a1(xk, xl)(uk − ul),

g1l =

n2∑
j=1

φ(‖yj − xl‖ − d0)~e(yj , xl) +

n2∑
j=1

a0(yj , xl)(vj − ul),

f2i =

n2∑
j=1

φ(‖yj − yi‖ − d2)~e(yj , yi) +

n2∑
j=1

a2(yj , yi)(vj − vi),

g2i =

n1∑
k=1

φ(‖xk − yi‖ − d0)~e(xk, yi) +

n1∑
k=1

a0(xk, yi)(uk − vi),

where the interaction function φ(z) is the derivation of ψ(z), ~e(z1, z2) is a unit
vector along the line of the two agents, expressed as

~e(z1, z2) =
z1 − z2
‖z1 − z2‖

,

and the functions a1, a2, a0 are the velocity matching coefficients.
Now, conclusively, we present the system of agents in two groups which will be

investigated in this article,

ẋl = ul,

u̇l =

n1∑
k=1

φ(‖xk − xl‖ − d1)~e(xk, xl) +

n1∑
k=1

a1(xk, xl)(uk − ul)

+

n2∑
j=1

φ(‖yj − xl‖ − d0)~e(yj , xl) +

n2∑
j=1

a0(yj , xl)(vj − ul),

ẏi = vi,

v̇i =

n2∑
j=1

φ(‖yj − yi‖ − d2)~e(yj , yi) +

n2∑
j=1

a2(yj , yi)(vj − vi)

+

n1∑
k=1

φ(‖xk − yi‖ − d0)~e(xk, yi) +

n1∑
k=1

a0(xk, yi)(uk − vi).

(2.5)

Remark 2.2. System (2.5) naturally corresponds to a graph G = (N , EN ) with a
block adjacency matrix induced by the velocity matching coefficients.



6 Q. ZHAO, S. SHI, W. LI EJDE-2021/104

According to the definition of flocking given by Cucker and Smale [12], we carry
out a new definition of two-groups flocking.

Definition 2.3. For agents in the two groups N1 and N2, they have a time-
asymptotic flocking if and only if for any l = 1, . . . , n1 and i = 1, . . . , n2, the
vector (xl(t), ul(t), yi(t), vi(t)) satisfies

(1) ‖ul(t)− vi(t)‖ = 0, as t→∞,
(2) there exists a constant C such that ‖xl(t)− yi(t)‖ ≤ C for all t ≥ 0.

Same as the case of one-group flocking, the above definition implies that if agents
in the two groups achieve flocking, their velocities are asymptotically equal, and in
the whole process, there is no agent get away from the groups.

For the convenience of studying, we then consider the translation invariance of
system (2.5). We translate the system as:

x̃l = xl −Xc, ỹi = yi −Xc,

ũl = ul − Vc, ṽi = vi − Vc,
where

Xc =
1

n

( n1∑
l=1

xl +

n2∑
i=1

yi

)
, Vc =

1

n

( n1∑
l=1

ul +

n2∑
i=1

vi

)
.

Noting that Ẋc = Vc and V̇c = 0, which means

Xc = Vc(0)t+Xc(0), Vc ≡ Vc(0),

we can obtain that the translation of the system has invariance. Thus, without loss
of generality, we may assume that Vc = Xc = 0 and this assumption can be realized
by the above translation. In other word, we only need to investigate system (2.5)
restricted on the invariant manifold M =M1 ×M2 with

M1 =
{

(x, y) :

n1∑
l=1

xl +

n2∑
i=1

yi = 0
}
, M2 =

{
(u, v) :

n1∑
l=1

ul +

n2∑
i=1

vi = 0
}
,

where u = (u1, . . . , un1
), v = (v1, . . . , vn2

).
With the definition of two-groups flocking, it is easy to discover that if the

system restricted on M achieves flocking, the velocities will asymptotically equal
to 0. Thus the flocking problem of system (2.5) on Rdn turns to be the stability
problem of fixed points on the manifold M.

3. Collective dynamic under global interaction

In this section, we consider the case when system (2.5) has global interaction,
which means any two agents in N interact with each other at any time. We first
present several assumptions on the interaction function and velocity matching co-
efficients.

(A1) The interaction function φ(z) : R → R is a C1-smooth monotonically in-
creasing function with φ(0) = 0 and there exists a constant C0 such that
for all z ∈ R, |φ(z)| ≤ C0.

(A2) The potential function

ψ(z) =

∫ z

0

φ(s) ds (3.1)

satisfies ψ(z)→∞ as z →∞.



EJDE-2021/104 DYNAMICS OF FLOCKING MODELS WITH TWO SPECIES 7

(A3) The velocity matching coefficients a1 ≡ c1, a2 ≡ c2, a0 ≡ c0, where c1, c2 ≥
0 and c0 > 0.

An example fulfilling (A1) and (A2) is the action function presented in [24]. It is
an s-type function defined by

φ(z) =
a+ b

2
σ(z + c) +

a− b
2

, (3.2)

where σ(z) = z/
√

1 + z2, b ≥ a > 0, and c = |a− b|/(2
√
ab). The corresponding ψ

is

ψ(z) =
a+ b

2
(
√

1 + (z + c)2 −
√

1 + c2) +
a− b

2
z.

By choosing the parameters as a = 2, b = 3, the graphs of φ(z) and ψ(z) can be
illustrated as following figures.

Figure 1. Functions φ(z) (left), and ψ (right)

Theorem 3.1. Let (A1)–(A3) hold for system (2.5). Then the agents in groups
N1 and N2 have a time-asymptotic flocking.

Proof. The structural energy of system (2.5) is indicated by a Hamiltonian

H(x, u, y, v) = V (x, y) +K(u, v), (3.3)

where V (x, y) is the potential energy induced by the differences of positions which
is defined by (2.4) and

K(u, v) =
1

2

n1∑
l=1

‖ul‖2 +
1

2

n2∑
i=1

‖vi‖2 (3.4)

is the kinetic energy induced by velocities.
Let (x(t), u(t), y(t), v(t)) be the solution of system (2.5) with the initial value

(x(0), u(0), y(0), v(0)) ∈M. With the assumptions (A1) and (A3), by differentiat-
ing H(x, u, y, v) along the solution with the respect of time, we have

Ḣ =
1

2

∑
k 6=l

φ(‖xk − xl‖ − d1)〈~e(xk, xl), uk − ul〉+

n1∑
l=1

〈 d
dt
ul, ul〉

+
1

2

∑
i 6=j

φ(‖yj − yi‖ − d2)〈~e(yj , yi), vj − vi〉+

n2∑
i=1

〈 d
dt
vi, vi〉

+

n1∑
l=1

n2∑
i=1

φ(‖xl − yi‖ − d0)〈~e(xl, yi), ul − vi〉

=c1

n1∑
k=1

n1∑
l=1

〈ul, uk − ul〉+ c2

n2∑
j=1

n2∑
i=1

〈vi, vj − vi〉+ c0

n2∑
i=1

n1∑
l=1

〈ul, vi − ul〉



8 Q. ZHAO, S. SHI, W. LI EJDE-2021/104

=c1

n1∑
k=1

n1∑
l=1

〈ul, uk − ul〉+ c2

n2∑
j=1

n2∑
i=1

〈vi, vj − vi〉+ c0

n2∑
i=1

n1∑
l=1

〈ul, vi − ul〉

=− 1

2
c1
∑
k,l

‖uk − ul‖2 −
1

2
c2
∑
j,i

‖vj − vi‖2 −
1

2
c0
∑
i,l

‖vi − ul‖2 ≤ 0,

and

Ḣ(x, u, y, v) = 0 ⇐⇒ ul = vi, l = 1, . . . , n1, i = 1, . . . , n2. (3.5)

Therefore,

H(x(t), u(t), y(t), v(t)) ≤ H0 := H(x(0), u(0), y(0), v(0)), t ≥ 0.

Let

Ω0 = {(x, u, y, v) ∈M : H(x, u, y, v) ≤ H0}.

Then Ω0 is a positive invariant set of system (2.5). Furthermore, Ω0 is bounded.
In fact, for any (x, u, y, v) ∈ Ω0, it is clear that K(u, v), V1(y), V2(x) are all less
than H0, so ‖ul‖ and ‖vi‖ are bounded, and by (A2), we deduce that ‖xk−xl‖ and
‖yj − yi‖ are bounded. As the fact that (x, u, y, v) ∈ M, we obtain that ‖xl‖ and
‖yi‖ are bounded. Therefore, Ω0 is a bounded positive invariant manifold.

From the LaSalle’s invariance principle [19], (x(t), u(t), y(t), v(t)) converges to
the largest invariant set in

S = {(x, u, y, v) ∈ Ω0 : Ḣ(x, u, y, v) = 0}.

By (3.5), we obtain that for any l = 1, . . . , n1, i = 1, . . . , n2,

‖ul(t)− vi(t)‖ → 0, as t→∞.

From this and the boundedness of positions, it follows that the agents have a time-
asymptotic flocking. �

Theorem 3.1 shows that the two groups of agents with global interaction fulfilling
the assumptions will unconditionally achieve flocking. However, in the real world,
one agent may not pay attention to all other agents at the same time. Therefore,
we will develop a more realistic model in the next section.

4. Collective dynamic under local interaction

In this section, we investigate the flocking system with local interaction, which
means the interaction between any two agents is considered to have finite cut-
off at a certain distance. As it is proposed in the above section, we used three
different distances d1, d2, and d0 to represent the balance distances of agents in
N1, agents in N2, and agents from different groups. In the same way, to realize
the cut-off of interaction between agents, we introduce three relevant non-influence
distances r1, r2, and r0 which are called the interaction radii. Similarly, r0 is not less
that r1 and r2. In addition, for the sake of convenience, we reorder the positions
and velocities of agents in N as q = (q1, . . . , qn) = (x1, . . . , xn1

, y1, . . . , yn2
) and

p = (p1, . . . , pn) = (u1, . . . , un1
, v1, . . . , vn2

).



EJDE-2021/104 DYNAMICS OF FLOCKING MODELS WITH TWO SPECIES 9

With the above arguments, we carry out the following local interaction functions
of agents in the two groups,

φij =


φ1(‖qj − qi‖) = ρ1(‖qj − qi‖)φ̂(‖qj − qi‖ − d1), i, j = 1, . . . , n1,

φ2(‖qj − qi‖) = ρ2(‖qj − qi‖)φ̂(‖qj − qi‖ − d2), i, j = n1 + 1, . . . , n,

φ0(‖qj − qi‖) = ρ0(‖qj − qi‖)φ̂(‖qj − qi‖ − d0), otherwise,

(4.1)

where φ̂ is a function fulfilling the assumptions (A1) and (A2), and

ρα(s) = ρ
( s
rα

)
, α = 0, 1, 2,

are three bump functions with ρ(z) a C1−smooth, uniformly bounded, and mono-
tonic decreasing function satisfying that ρ(z) = 0, for all z ≥ 1 and ρ(z) = 1, for
all z ≤ 1− ε, here 0 < ε� 1 is a sufficiently small parameter. One suitable choice
is the following function introduced in [24],

ρ(z) =


1, z ∈ (0, h),
1
2

(
1 + cos(π z−h1−h )

)
, z ∈ [h, 1],

0, otherwise,

where h = 1 − ε. The use of these bump functions ensures that the interactions
between the agents smoothly vanish as their distances go beyond the interaction
radii.

Then, on the basis of the above discussion, we establish the following system to
describe the dynamic of agents in two groups under local interaction,

q̇i = pi,

ṗi =

n∑
j=1

φij~e(qj , qi) +

n∑
j=1

aij(pj − pi)
(4.2)

where

aij = a(‖qj − qi‖) =


c1ρ

1(‖qj − qi‖), i, j = 1, . . . , n1,

c2ρ
2(‖qj − qi‖), i, j = n1 + 1, . . . , n,

c0ρ
0(‖qj − qi‖), otherwise,

(4.3)

with c1, c2, and c0 are positive constants.
Similar to the remark in section 2, the above system (4.2) naturally corresponds

to a graphGN (q) = (N , EN (q)) with (i, j) ∈ EN (q) if and only if aij = a(‖qj−qi‖) 6=
0, i, j = 1, . . . , n. Then, the structural energy of system (4.2) is

Ĥ(q, p) = V̂ (q) +K(p),

with

V̂ (q) =
1

2

n1∑
i,j=1

ψ1(‖qj − qi‖) +
1

2

n∑
i,j=n1+1

ψ2(‖qj − qi‖)

+

n1∑
j=1

n∑
i=n1+1

ψ0(‖qj − qi‖),
(4.4)

where

ψα(z) =

∫ z

dα

φα(s) ds, α = 0, 1, 2. (4.5)



10 Q. ZHAO, S. SHI, W. LI EJDE-2021/104

We remark that ψα(z) has the following properties:

ψα(z) ≡ hα := ψα(rα), z ≥ rα,
ψα(z) < hα, dα ≤ z < rα, α = 0, 1, 2.

(4.6)

Now we can state the following conditional flocking result for the local model.

Theorem 4.1. Let φij and aij in system (4.2) be given by (4.1) and (4.3), respec-
tively. If for any t ≥ 0, N is a cluster with connected graph GN = (N , EN ), then
the agents have a time-asymptotic flocking.

Proof. Let (q(t), p(t)) be the solution of system (4.2) with (q(0), p(0)) ∈ M. By

differentiating Ĥ(q, p) along the solution with the respect of time, we have

˙̂
H =

1

2

n1∑
i,j=1

φ1(‖qj − qi‖)〈~e(qj , qi), pj − pi〉

+
1

2

n∑
i,j=n1+1

φ2(‖qj − qi‖)〈~e(qj , qi), pj − pi〉

+

n1∑
i=1

n∑
j=n1+1

φ0(‖qj − qi‖)〈~e(qj , qi), pj − pi〉+

n∑
i=1

〈 d
dt
pi, pi〉

=− 1

2

n1∑
i,j=1

c1ρ
1(‖qj − qi‖)〈pj − pi, pj − pi〉

− 1

2

n∑
i,j=n1+1

c2ρ
2(‖qj − qi‖)〈pj − pi, pj − pi〉

− 1

2

n1∑
i=1

n2∑
j=n1+1

c0ρ
0(‖qj − qi‖)〈pj − pi, pj − pi〉

=− 1

2

n∑
i,j=1

aij‖pj − pi‖2 ≤ 0,

(4.7)

therefore,

Ĥ(q(t), p(t)) ≤ Ĥ0 := Ĥ(q(0), p(0)). (4.8)

We claim that if
˙̂
H(q, p) = 0, then pj − pi = 0, i, j = 1, . . . , n. For any agents

(qj , pj), (qi, pi) ∈ N , by the connectivity of GN (q) = (N , EN (q)), there exists a
walk (k1, . . . , km) which goes through agents (qk1 , pk1), . . . , (qkmpkm) in N with
(qk1 , pk1) = (qj , pj), (qkm , pkm) = (qi, pi), and akiki−1

> 0, i = 1, . . . ,m. By (4.7),
we can deduce that pk1 − pk2 = pk2 − pk3 = · · · = pkm−1 − pkm = 0, i.e., pj − pi = 0.

Let Ω̂ = Ω̂0 ∩ Σ with

Ω̂0 = {(q, p) ∈M : Ĥ(q, p) ≤ Ĥ0}, Σ = {(q, p) ∈M : {q1, . . . , qn} is a cluster}.

Since for any t > 0, GN is connected, combining (4.8), we know that Ω̂ is a

positive invariant manifold of system (4.2). Furthermore, for any (q, p) ∈ Ω̂, we

have K(p) < Ĥ0 and ‖p‖ is bounded, and {q1, . . . , qn} is a cluster, which leads to

‖q‖ is bounded. Thus, Ω̂ is a bounded positive invariant manifold of (4.2).



EJDE-2021/104 DYNAMICS OF FLOCKING MODELS WITH TWO SPECIES 11

Then, from LaSalle’s invariance principles, the solution (q(t), p(t)) of (4.2) start-

ing in Ω̂ converges to the largest invariant set in

S = {(q, p) ∈ Ω̂| ˙̂H(q, p) = 0}.
From (4.7), we can deduce that for i, j = 1, . . . , n,

‖pj(t)− pi(t)‖ → 0, as t→∞.
The boundedness of ‖q‖ for all time t can be deduced by the connectivity of GN (q),
consequently, the agents have a time-asymptotic flocking. �

Remark 4.2. Because of the dissipation of Ĥ, the assumption that N is a cluster
at any time could be replaced by proposing a constraint on the initial value of the

energy function, for example, Ĥ0 ≤ min{hα, ψα(0), α = 0, 1, 2}.

Besides the above conditional flocking consideration, we investigate the stability
of a class of flocking solutions which can be, in fact, entirely viewed as an invariant
manifold of system (4.2). We first present the following propositions.

Proposition 4.3. For q ∈ M1, if GN (q) is connected, then q ∈ B := {q ∈
M1| ‖q‖ ≤ n2r0}.

Proof. If GN (q) is connected, then max{‖qj − qi‖, i, j = 1, . . . , n} ≤ nr0, thus

‖qi‖ =
1

n
‖nqi‖ =

1

n
‖(n− 1)qi −

∑
k 6=i

qk‖ ≤
1

n

∑
k 6=i

‖qk − qi‖ ≤ (n− 1)r0.

Consequently, ‖q‖ ≤ n(n− 1)r0 < n2r0 and this leads to q ∈ B. �

Proposition 4.4. For any q ∈ M1, if GN (q) is disconnected, then there exist a

q∗ ∈M1 such that GN (q∗) is connected and V̂ (q∗) < V̂ (q).

Proof. Since GN (q) is disconnected, the agents in N can be decomposed into two
sets V1, V2, where V1 is a cluster, as well, V1 and V2 are a pair of separated sets.
Let Vij = Vi ∩ Nj , i, j = 1, 2. Then V1 = V11 ∪ V12, V2 = V21 ∪ V22, and we note
that

(qi, pi) ∈ Vi, (qij , pij) = ((qij1 , p
ij
1 ), . . . , (qijnij , p

ij
nij )) ∈ Vij , i, j = 1, 2,

where nij is the number of agents in Vij . Then, by (4.4), we derive

V̂ (q) = V 1(q1) + V 2(q2) + n11n21h1 + n12n22h2 + n12n21h0 + n11n22h0,

where

V 1(q1) =
1

2

n11∑
i,j=1

ψ1(‖q11j − q11i ‖) +
1

2

n12∑
i,j=1

ψ2(‖q12j − q12i ‖) +

n11∑
j=1

n12∑
i=1

ψ0(‖q11j − q12i ‖),

V 2(q2) =
1

2

n21∑
i,j=1

ψ1(‖q21j − q21i ‖) +
1

2

n22∑
i,j=1

ψ2(‖q22j − q22i ‖) +

n21∑
j=1

n22∑
i=1

ψ0(‖q21j − q22i ‖).

As V1 and V2 are separated sets, we have

δ1 := d(V11,V22)− r0 > 0, δ2 := d(V11,V21)− r1 > 0,

δ3 := d(V12,V21)− r0 > 0, δ4 := d(V12,V22)− r2 > 0 .

For A1, A2 ⊂ N , we define

d(A1, A2) = min{‖qj − qi‖, (qi, pi) ∈ A1, (qj , pj) ∈ A2}.



12 Q. ZHAO, S. SHI, W. LI EJDE-2021/104

Without loss of generality, we suppose that δ1 = min{δ1, δ2, δ3, δ4}. Then there
exist (q11i0 , p

11
i0

) ∈ V11, (q22j0 , p
22
j0

) ∈ V22, such that

‖q11i0 − q
22
j0 ‖ − r0 = δ1 > 0. (4.9)

Let q∗ = (q1∗, q2∗) with

q1∗ = q1 − ~e(q11i0 , q
22
j0 )δ, q2∗ = q2 + ~e(q11i0 , q

22
j0 )

n11 + n12
n21 + n22

δ,

where δ is a positive constant satisfying

(1 +
n11 + n12
n21 + n22

)δ = δ1 + ε

for some sufficiently small ε > 0. Then we have

‖q11∗j − q11∗i ‖ = ‖q11j − q11i ‖, i, j = 1, . . . , n11,

‖q12∗j − q12∗i ‖ = ‖q12j − q12i ‖, i, j = 1, . . . , n12,

‖q21∗j − q21∗i ‖ = ‖q21j − q21i ‖, i, j = 1, . . . , n21,

‖q22∗j − q22∗i ‖ = ‖q22j − q22i ‖, i, j = 1, . . . , n22;

therefore,
V 1(q1∗) = V 1(q1), V 2(q2∗) = V 2(q2). (4.10)

Meanwhile, we infer that

‖q11∗i0 − q
22∗
j0 ‖ =

∥∥(‖q11i0 − q22j0 ‖ − δ − n11 + n12
n21 + n22

δ
)
~e(q11i0 , q

22
j0 )
∥∥

= ‖q11i0 − q
22
j0 ‖ − (δ1 + ε)

= r0 − ε.

(4.11)

Then as δ2 ≥ δ1,

‖q11∗i − q21∗j ‖ = ‖q11i − q21k −
(
δ +

n11 + n12
n21 + n22

δ
)
~e(q11i0 , q

22
j0 )‖

≥ ‖q11i − q21k ‖ − (δ1 + ε)

≥ δ2 + r1 − (δ1 + ε)

≥ r1 − ε, i = 1, . . . , n11, k = 1, . . . , n21.

Similarly,

‖q22∗j − q12∗l ‖ ≥ r2 − ε, j = 1, . . . , n22, l = 1, . . . , n12,

‖q21∗k − q12∗l ‖ ≥ r0 − ε, k = 1, . . . , n21, l = 1, . . . , n12.

By (4.9), for any i 6= i0 and j 6= j0,

‖q11∗i − q22∗j ‖ ≥ ‖q11i − q22j ‖ − (δ1 + ε) ≥ ‖q11i0 − q
22
j0 ‖ − (δ1 + ε) = r0 − ε. (4.12)

Combining (4.11)-(4.12), with the properties shown in (4.6), we have
n11∑
i=1

n22∑
j=1

ψ0(‖q22∗j − q11∗i ‖) < n11n22h0,

n11∑
i=1

n21∑
j=1

ψ1(‖q21∗j − q11∗i ‖) ≤ n11n21h1,

n12∑
i=1

n22∑
j=1

ψ2(‖q22∗j − q12∗i ‖) ≤ n12n22h2,
n12∑
i=1

n21∑
j=1

ψ0(‖q21∗j − q12∗i ‖) ≤ n12n21h0.

Adding the above expression and (4.10), we have V̂ (q∗) < V̂ (q).



EJDE-2021/104 DYNAMICS OF FLOCKING MODELS WITH TWO SPECIES 13

If GN (q∗) is still disconnected, we can reuse the above process to get an new q∗

such that GN (q∗) is connected and V̂ (q∗) < V̂ (q). �

Let M = Mq ×Mp with

Mq =
{
q ∈M1 : V̂ (q) = min

y∈M1

{V̂ (y)}
}
, Mp = {p ∈M2|pi = 0}.

We then state the following properties of Mq.

Lemma 4.5. Mq is a nonempty compact subset of M1.

Proof. For any q ∈ M1 with GN (q) disconnected, by Proposition 4.3, there exist

a q∗ ∈ M1 such that GN (q∗) is connected and V̂ (q∗) < V̂ (q). By Proposition 4.1,

we know that q∗ ∈ B. Thus, the minimum of V̂ (q) in M1 is just the minimum in
B.

Then from (4.4), V̂ (q) is a continuous function which must have minimum in B.
Thus, Mq is not empty. Meanwhile, Mq ⊆ B, therefore Mq is compact. �

By Lemma 4.5, M is nonempty. Moreover, for any (q, p) ∈ M , it is obviously
that (q, p) is a flocking solution of system (4.2), thus, M is an invariant manifold
of (4.2) and dim(M) = dim(Mq).

To analyze the local stability of M , we ought to investigate the Jacobian matrix
of the system (4.2) calculated on (q, p) ∈M . Noting that

J =

(
0 I
G −L

)
, (4.13)

where G = −∇2V̂ (q), L = L⊗I with L is the Laplacian matrix defined by (2.1) and
I is unit matrices with corresponding order, we have the following observations.

Lemma 4.6. ker(L) ⊂ ker(G).

Proof. Noting that ker(L) ⊂ ker(G) is just equivalent to that if Lξ = 0 for some
ξ ∈ Rd, then Gξ = 0. Meanwhile, one can induce by (2.1) that if Lξ = 0, then
ξ = ξ01d with ξ0 ∈ R. Hence, we only need to prove that G1d = 0.

We observe that G := [Gij ]n×n, where Gij = −∇j∇iV̂ (q) is a d × d matrix.

Therefore, G1d = 0 if and only if

n∑
j=1

Gij = O, i = 1, . . . , n, (4.14)

where O is a d× d zero matrix.
To prove (4.14), we rewrite V̂ (q) as

V̂ (q) =
1

2

n1∑
k,l=1

U1(‖qk − ql‖2) +
1

2

n∑
i,j=n1+1

U2(‖qj − qi‖2)

+

n1∑
k=1

n∑
i=n1+1

U0(‖qk − qi‖2),

where Uα(s) = ψα(s1/2), s ≥ 0, α = 0, 1, 2.



14 Q. ZHAO, S. SHI, W. LI EJDE-2021/104

For i ∈ {n1 + 1, . . . , n}, k ∈ {1, . . . , n1} and j ∈ {n1 + 1, . . . , n}, we have

∇iV̂ (q) = ∇i
n∑

j=n1+1

U2(‖qj − qi‖2) +∇i
n1∑
k=1

U0(‖qk − qi‖2)

= 2

n∑
j=n1+1

∂

∂qi
U2(‖qj − qi‖2)(qj − qi) + 2

n1∑
k=1

∂

∂qi
U0(‖qk − qi‖2)(qk − qi),

∇2
i V̂ (q) = 2

( n∑
j=n1+1

∂

∂qi
U2(‖qj − qi‖2)

)
Id + 2

( n∑
k=1

∂

∂qi
U0(‖qk − qi‖2)

)
Id

+ 4

n∑
j=n1+1

∂2

∂qi∂qj
U2(‖qj − qi‖2)‖qj − qi‖2

+ 4
n∑
k=1

∂2

∂qi∂qk
U0(‖qk − qi‖2)‖qk − qi‖2,

∇j∇iV̂ (q) = −2
( ∂

∂qi
U2(‖qj − qi‖2)

)
Id − 4

∂2

∂qi∂qj
U2(‖qj − qi‖2)‖qj − qi‖2,

∇k∇iV̂ (q) = −2
( ∂

∂qi
U0(‖qk − qi‖2)

)
Id − 4

∂2

∂qi∂qk
U0(‖qk − qi‖2)‖qk − qi‖2,

which leads to
n∑
j=1

Gij = Gii +

n1∑
k=1

Gik +

n∑
j=n1+1

Gij = O.

In the similar way, (4.14) also holds for i ∈ {1, . . . , n1}, k ∈ {1, . . . , n1} and j ∈
{n1 + 1, . . . , n}. This completes the proof. �

Proposition 4.7. dim(ker(G)) ≥ dim(Mq).

Proof. For any smooth curve Γ(s) ∈Mq with Γ(0) = q ∈Mq, we have

∇V̂ (Γ(s)) = 0, ∇2V̂ (Γ(s))
d

ds
Γ(s) = 0,

this means that d
dsΓ(s)

∣∣
s=0

, the tangent vector of Mq at q, belongs to ker(G).
Therefore, the dimension of ker(G) is not less than that of the tangent space of Mq.
Thus dim(ker(G)) ≥ dim(Mq). �

Now we can state and prove a result about the local stability of two-groups
flocking solutions under local interaction.

Theorem 4.8. Let φij and aij in system (4.2) be given by (4.1) and (4.3), re-
spectively. Assuming that dim(ker(G)) = dim(Mq), then the manifold M is locally
asymptotically stable.

Proof. It is not difficult to see that (u, v)> ∈ ker(J) if and only if

v = 0,

Gu− Lv = 0,
(4.15)

or equivalently, u ∈ ker(G) and v = 0, so dim(ker(J)) = dim(ker(G)). Combining
this with the assumption that dim(ker(G)) = dim(Mq), we obtain

dim(ker(J)) = dim(Mq) = dim(M).



EJDE-2021/104 DYNAMICS OF FLOCKING MODELS WITH TWO SPECIES 15

Now we show that if λ 6= 0 is an eigenvalue of J , then Re(λ) < 0. Let (ξ, η)> be
the eigenvector corresponding to λ, then

η = λξ,

Gξ − Lη = λη,

which leads to
λ2ξ + λLξ − Gξ = 0.

premultiplying by ξ>, we solve the eigenvalue λ as

λ =
−ξ>Lξ ±

√
(ξ>Lξ)2 + 4‖ξ‖2ξ>Gξ

2‖ξ‖2
.

From Lemma 4.6 and the property of L in (2.2), we conclude that ξ>Lξ > 0.
Moreover, since q ∈ Mq is a minimum point, ξ>Gξ ≤ 0. As mentioned above, we
get Re(λ) < 0.

Combining this with Lemma 4.5, we state that M is a normally hyperbolic
invariant manifolds onM, which induce that M is locally asymptotically stable. �

Remark 4.9. 1. Theorem 4.8 leads to the flocking solutions in manifold M being
locally asymptotically stable.

2. The assumption dim(ker(G)) = dim(Mq) can happen in many situations. For

example, when n1 = n2 = 1, φα(z) = ραz, and dα < rα ≤
√

2dα, α = 0, 1, 2, we
have dim(ker(G)) = dim(Mq) = d.

5. Numerical simulations

We have proved in the above sections that the systems of two groups with global
and local interaction could both achieve flocking under some proper assumptions. In
this section, some numerical results will be presented to show the above conclusions
visually. In addition, we would like to remark that the relationship between the
parameters and the final configuration of the flocking solutions is also a notable
problem, and we will display a few numerical results to show the related phenomena
in subsection 5.3.

In all the following figures, the positions and the velocities of agents in N1 are
showed by red stars and red arrows, while inN2 are showed by blue squares and blue
arrows. For the two groups of agents with global interaction, they follow system
(2.5) in section 3, where φ is chosen as (3.2). For agents with local interaction, they

follow system (4.2) in section 4 with φ̂ chosen as (3.2). Both the parameters in φ

and φ̂ are chosen as a = 1, b = 7. For the velocity matching coefficients, we find
that their differences do not have much effect on the final configurations, thus for
convenience, we choose c0 = c1 = c2 = 1.

5.1. One group versus two groups. We present four sets of figures showing the
time-asymptotical flocking of one group and two groups of agents under global and
local interactions, respectively. The number of agents is n = 30 with n1 = n2 = 15.
It is clear that if the balance distances are chosen as d0 = d1 = d2, the two-groups
system turns into a one-group system. The initial positions of all the four sets of
figures are chosen at random from [0, 160]2 and the initial velocities are chosen at
random from [0, 10]2.

Figures 2–5 show that for random initial data fulfilling the connective require-
ment of the local model, both two systems can get a flocking consensus. Meanwhile,



16 Q. ZHAO, S. SHI, W. LI EJDE-2021/104

Figure 2. One-group system under global interactions with d0 =
d1 = d2 = 40.

when considering the patterns, there shows to have some differences among these
figures: for the one-group system, the agents mix together at the final states, while
for the two-groups system, the agents present a layering phenomenon as they tend
to get close with agents in their own group. Furthermore, we can also find that
under local interactions, the agents have a clear lattice structure, while under global
interactions, they do not have a regular pattern.

5.2. Global interaction versus local interaction. In this subsection, we present
three sets of figures to show the effects on the final flocking states from connectivity
of the initial data.

In Figures 6–8, the initial velocities are chosen at random from the set [−5, 15]2

and the initial positions are given in the corresponding captions. The balance dis-
tances are chosen as d1 = d2 = 40, d0 = 80, and in the local model, the interaction
radii are chosen as r1 = r2 = 56, r0 = 100. The flocking states of local and global
systems are separately presented in parts (b) and (c) of each figure.



EJDE-2021/104 DYNAMICS OF FLOCKING MODELS WITH TWO SPECIES 17

Figure 3. Two-group system under global interactions with d1 =
d2 = 40 and d0 = 60.

From the numerical results, we see that under global interactions, no matter the
initial data is connected or not, the agents always tend to form a whole cluster and
achieve flocking. While under local interactions, the connectivity of initial data
could affect the flocking of the agents and cause different final states: Figure 6
shows that the initial data are totally connected and the agents achieve flocking
with d(N1,N2) = d0 = 80; Figure 7 shows that N1 and N2 are a pair of separated
sets in the initial state and the agents do not achieve flocking all along the time;
Figure 8 shows that the initial data of N2 are not connected, and in the final state,
the agents in N2 separated into two parts, one part form a cluster with N1 while
the other part get away from the cluster. We remark that for some disconnected
initial data, the agents do not achieve flocking under local interactions.

5.3. Interlaced versus separated. As it has been referred in subsection 5.1, the
change of the balance distances and interaction radii could ultimately affect the
patterns of the two groups. Hence in this subsection, we will use three sets of
figures to perform the effect on the patterns caused by d0 and r0. Here, the agents



18 Q. ZHAO, S. SHI, W. LI EJDE-2021/104

Figure 4. One-group system under local interactions with d0 =
d1 = d2 = 40 and r0 = r1 = r2 = 56.

are considered under the local interactions, and three sets of initial data are chosen
to show this phenomenon.

To describe the interlacing state of the two groups N1,N2, we introduce the
following criterion: if d(N1,N2) = d0 and for both of the two groups, there exist
three agents A,B,C from one groups and one agent O from the other group, such
that O inside 4ABC, then N1,N2 are interlaced. Otherwise, N1,N2 are strongly
separated.

For Figures 9–11, the initial velocities of agents are chosen at random from
the set [0, 20]2 and the initial positions are given in their captions. The balance
distances and interaction radii inside the two groups are chosen as d1 = d2 = 40
and r1 = r2 = 56. Then we gradually increased d0 and r0 to find the variation
of the configuration. In each set of figures, we present part (b) as a one-group
contrast, part (c) as an interlaced pattern and part (d) as a strongly separated
pattern. Their relative d0 and r0 are placed in the caption.



EJDE-2021/104 DYNAMICS OF FLOCKING MODELS WITH TWO SPECIES 19

Figure 5. Two-group system under local interactions with d1 =
d2 = 40, d0 = 60 and r1 = r2 = 56 and r0 = 84.

(a) Initial status (b) Local model (c) Global model

Figure 6. The initial positions of agents in N1 and N2 are chosen
at random from [0, 200]2 and [250, 450]2.

For initial data in Figures 9–11, we can always find a d0 such that the two groups
are interlaced, and another d0 large enough such that d(N1,N2) = d0 and the two



20 Q. ZHAO, S. SHI, W. LI EJDE-2021/104

(a) Initial status (b) Local model (c) Global model

Figure 7. The initial positions of agents in N1 and N2 are chosen
at random from [0, 200]2 and [300, 450]2.

(a) Initial status (b) Local model (c) Global model

Figure 8. The initial positions of agents in N1 and N2 are chosen
as random from [0, 250]2 and [300, 550]2.

groups are strongly separated. Thus we would like to remark that the larger the
difference between d0 and d1, d2, the stronger the tendency of agents in different
groups to be separated.

Acknowledgments. This work was supported by the NSFC grant (No. 11771177,
12171197), by the China Automobile Industry Innovation and Development Joint
Fund (No. U1664257), by the Program for Changbaishan Scholars of Jilin Province
and Science and Technology Development Project of Jilin Province
(No. YDZJ202101ZYTS141, 20190201132JC).

References

[1] H. Bae, S. Cho, S. Lee, S. Yun; A particle model for the herding phenomena induced by
dynamic market signals, Journal of Statistical Physics, 177 (2019), 365–398.

[2] E. Ben-Jacob; Bacterial self-organization: co-enhancement of complexification and adaptabil-

ity in a dynamic environment, Philosophical Transactions of the Royal Society, 361 (2003),
1283–1312.

[3] J. A. Carrillo, M. Fornasier, J. Rosado, G. Toscani; Asymptotic flocking dynamics for the

kinetic cucker smale model, SIAM Journal on Mathematical Analysis, 42 (2009), 218–236.
[4] J. A. Carrillo, M. Fornasier, G. Toscani, F. Vecil; Particle, kinetic, and hydrodynamic models

of swarming, Mathematical Modeling of Collective Behavior in Socio-Economic and Life
Sciences, Series: Modelling and Simulation in Science and Technology, Birkhäuser, 2010,

pp. 297–336.



EJDE-2021/104 DYNAMICS OF FLOCKING MODELS WITH TWO SPECIES 21

(a) Initial status (b) d0 = 40 and r0 = 56

(c) Interlaced: d0 = 54 and r0 = 76 (d) Strongly separated: d0 = 70 and r0 = 98

Figure 9. The number of agents are n1 = n2 = 30. The agents
have a strict defined initial positions shown in figure(a) where the
initial positions of the two groups have one column of intersect.

[5] D. Chang, S. Shadden, J. Marsden, R. Olfati-Saber; Collision avoidance for multiple agent

systems, Proceedings of the IEEE Conference on Decision and Control, 1 (2003), 539–543.
[6] S. Chen, Y. Qiu, J. Liu, S. Nie; Fast flocking algorithm of multi-agent network via community

division, Control and Decision, 33 (2018), 1523–1526.

[7] Y. Choi, J. Haskovec; Cucker-smale model with normalized communication weights and time

delay, Kinetic & Related Models, 2 (2017), 1033–1039.
[8] Y. Choi, B. Kwon; Two-species flocking particles immersed in a fluid, Communications in

Information and Systems, 13 (2013), 123–149.
[9] Y. Choi, S. Salem; Cucker-smale flocking particles with multiplicative noises: Stochastic

mean-field limit and phase transition, Kinetic & Related Models, 12 (2019), 573–592.

[10] F. Cucker, C. Huepe; Flocking with informed agents, Maths in Action, 1 (2008), 1–25.
[11] F. Cucker, E. Mordecki; Flocking in noisy environments, Journal De Mathematiques Pures

Et Appliquees, 89 (2008), 278–296.



22 Q. ZHAO, S. SHI, W. LI EJDE-2021/104

(a) Initial status d0 = 40 and r0 = 56

(c) Interlaced: d0 = 48 and r0 = 67 (d) Strongly separated: d0 = 80 and r0 = 113

Figure 10. The number of agents in the two groups are n1 =
n2 = 30. The initial positions of the agents in N1 are chosen at
random in the set [0, 150]2 and of the agents in N2, in the set
[100, 250]2.

[12] F. Cucker, S. Smale; On the mathematics of emergence, Japanese Journal of Mathematics,

2 (2007), 197–227.
[13] F. Dalmao, E. Mordecki; Cucker-smale flocking under hierarchical leadership and random

interactions, SIAM Journal on Applied Mathematics, 71 (2009), 1307–1316.

[14] J. Dong, S. Ha, D. Kim; Emergent behaviors of the kuramoto model with a time delay on a

general digraph, SIAM Journal on Applied Dynamical Systems, 19 (2020), 304–328.
[15] J. Dong, S. Ha, D. Kim; Emergence of mono-cluster flocking in the thermomechanical cucker-

smale model under switching topologies, Analysis and Applications, 19 (2021), 305–342.
[16] C. Godsil, G. F. Royle; Algebraic graph theory, Graduate Texts in Mathematics, vol. 207,

Springer-Verlag, New York, 2001.

[17] S. Ha, K. Lee, D. Levy; Emergence of time-asymptotic flocking in a stochastic cucker-smale
system, Communications in Mathematical Sciences, 7 (2009), 453–469.

[18] M. O. Jackson; Social and economic networks, Princeton University Press, 2010.



EJDE-2021/104 DYNAMICS OF FLOCKING MODELS WITH TWO SPECIES 23

(a) Initial status (b) d0 = 40 and r0 = 56

(c) Interlaced: d0 = 48 and r0 = 67 (d) Strongly separated: d0 = 60 and r0 = 84

Figure 11. The number of the agents are n1 = 5 and n2 = 55.
The initial positions of agents N1 are chosen at random from
[50, 150]2 and of agents in N2 are chosen at random from [0, 250]2.

[19] J. P. LaSalle; Some extensions of liapunov’s second method, Ire Transactions on Circuit

Theory, 7 (1961), 520–527.
[20] Y. Lei and H. Zhao; Negative sobolev spaces and the two-species vlasov-maxwell-landau sys-

tem in the whole space, Journal of Functional Analysis, 267 (2014), 3710–3757.

[21] I. I. Matveeva; Exponential stability of solutions to nonlinear time-varying delay systems of

neutral type equations with periodic coefficients, Electronic Journal of Differential Equations,
2020 (2020), No. 20, 1-12.

[22] B. Mohar; The laplacian spectrum of graphs, Graph Theory, Combinatorics, and Applica-
tions, 2 (1991), 871–898.

[23] S. Motsch, E. Tadmor; Heterophilious dynamics enhances consensus, SIAM Review, 56

(2014), 577–621.
[24] R. Olfati-Saber; Flocking for multi-agent dynamic systems: Algorithms and theory, IEEE

Transactions on Automatic Control, 51 (2006), 401–420.



24 Q. ZHAO, S. SHI, W. LI EJDE-2021/104

[25] R. Olfati-Saber, R. M. Murray; Consensus problems in networks of agents with switching

topology and time-delays, IEEE Transactions on Automatic Control, 49 (2004), 1520–1533.

[26] D. Paley, N. E. Leonard, R. Sepulchre, D. Grunbaum, J. K. Parrish; Oscillator models and
collective motion: spatial patterns in the dynamics of engineered and biological networks, Die

Makromolekulare Chemie, 192 (2007), 2293–2305.

[27] L. Perea, P. Elosegui, G. Gomez; Extension of the cucker-smale control law to space flight
formations, Journal of Guidance Control and Dynamics, 32 (2009), 527–537.

[28] C. W. Reynolds; Flocks, herds and schools: A distributed behavioral model, ACM SIGGRAPH

Computer Graphics, 21 (1987), 25–34.
[29] K. Saulnier, D. Saldaa, A. Prorok, G. Pappas, V. Kumar; Resilient flocking for mobile robot

teams, IEEE Robotics and Automation Letters 2 (2017), 1039–1046.

[30] J. Shen; Cucker-smale flocking under hierarchical leadership, SIAM Journal on Applied
Mathematics, 68 (2007), 694–719.

[31] S. Su, G. Zhang; Global stability of traveling waves for delay reaction-diffusion systems
without quasi-monotonicity, Electronic Journal of Differential Equations, 2020 (2020), No.

46, 1-18.

[32] E. Tadmor, S. Ha; From particle to kinetic and hydrodynamic descriptions of flocking, Kinetic
& Related Models 1 (2008), 415–435.

[33] J. Tello, M. Winkler; Stabilization in a two-species chemotaxis system with a logistic source,

Nonlinearity 25 (2012), 1413–1425.
[34] C. Viragh, G. Vasarhelyi, N. Tarcai, T. Szrenyi, G. Somorjai, T. Nepusz, T. Vicsek; Flocking

algorithm for autonomous flying robots, Bioinspiration & Biomimetics 9 (2013), 025012.

[35] Q. Xiao, H. Liu, Z. Xu, Z. Ouyang; On collision avoiding fixed-time flocking with measurable
diameter to a cucker smale type self-propelled particle model, Complexity 11 (2020), 1–12.

Qingjian Zhao

School of Mathematics, Jilin University, Changchun 130012, China
Email address: zhaoqj17@mails.jlu.edu.cn

Shaoyun Shi

School of Mathematics & State key laboratory of automotive, simulation and control,
Jilin University, Changchun 130012, China

Email address: shisy@jlu.edu.cn

Wenlei Li (corresponding author)

School of Mathematics, Jilin University, Changchun 130012, China

Email address: lwlei@jlu.edu.cn


	1. Introduction
	2. Flocking model for two groups
	3. Collective dynamic under global interaction
	4. Collective dynamic under local interaction
	5. Numerical simulations
	5.1. One group versus two groups
	5.2. Global interaction versus local interaction
	5.3. Interlaced versus separated
	Acknowledgments

	References