Control of Large Swarms via Random Finite Set Theory

arXiv:1801.07314v3 [cs.SY] 27 Feb 2018

Bryce Doerr1 and Richard Linares2

Abstract— Controlling large swarms of robotic agents has
many challenges including, but not limited to, computational
complexity due to the number of agents, uncertainty in the
functionality of each agent in the swarm, and uncertainty in the
swarm’s configuration. This work generalizes the swarm state
using Random Finite Set (RFS) theory and solves the control
problem using model predictive control which naturally handles
the challenges. This work uses information divergence to define
the distance between swarm RFS and a desired distribution.
A stochastic optimal control problem is formulated using a
modified L22 distance. Simulation results show that swarm
densities converge to a target destination, and the RFS control
formulation can vary in the number of target destinations.

I. INTRODUCTION
Control of large robotic groups or swarms is currently
an area of great interest in controls research. These robotic
groups are usually tiny robots with limited actuators that do
specific tasks in the formation. For example, the group can
use its combined effort to grasp or move in the environment
[1]. Specifically, in NASA applications, swarm control of
space systems including satellites and rovers can be used
for exploration of asteroids and other celestial bodies of
interest [2]. In military applications, swarms of UAVs can
be used for border patrol, search and rescue, surveillance,
communications relaying, and mapping of hostile territory
[3]. In reviewing these applications, swarm control of large
groups of swarming agents are required.
Changing how the swarm is modeled can alleviate computational complexity for different control methods. Previously,
the swarm/potential model using the random finite set (RFS)
formalism was used to describe the temporal evolution of
the probabilistic description of the robot swarm to promote
decentralized coordination [6]. By using a measure-value
recursion of the RFS formalism for the agents of the swarm,
the dynamics of the swarm in time can be determined with
computational efficiency.
Several control techniques have been implemented on
swarms previously. In centralized control, one agent in the
swarm is the designated controller and manages the execution of other agents in the swarm [7]. Thus, it oversees the
other agents’ system processes. This type of control is the
easiest to implement on a swarm of robots, and it solves
the problem of controlling the swarm to a desired state.
Unfortunately, this type of control suffers from two main
1 Bryce Doerr, Graduate Student, Department of Aerospace Engineering
and Mechanics, University of Minnesota, Minneapolis, MN 55455, USA

doerr024@umn.edu
2 Richard Linares, Assistant Professor, Department of Aerospace Engineering and Mechanics, University of Minnesota, Minneapolis, MN 55455,
USA rlinares@umn.edu

problems. As the number of agents of the swarm increases,
the computational workload becomes more difficult. This
is especially difficult when the swarm robots are low cost
robotic agents. Another issue that arises is that centralized
control is not robust against failures of other agents within
the swarm [8]. Since swarm robots on a level of a thousand
agents are low cost units, communication, actuation, and
sensing are performed with less reliability. Thus, centralized
control may not be a viable option for these systems.
Probabilistic swarm guidance has also been used to enable
the swarm to converge to target distributions [9]. Probabilistic swarm guidance controls a swarm density distribution
through distributed control, so that each agent determines
its own trajectory while the swarm converges to the target
distribution. Distributed control is defined as reformulating a
control problem as a set of interdependent subproblems and
solving these subproblems [5]. Probabilistic swarm guidance
solves issues that involve controlling an increasing number
of individual agents which is computationally complex by
controlling the swarm density distribution of the agents.
Using optimal transport, convergence is achieved faster than
a homogeneous Markov chain approach and the cost function
is minimized [9], [10], [11]. Another method is the use of the
inhomogeneous Markov chain [12], [13]. In this approach,
the agents move in a decentralized fashion which is like
the homogeneous Markov chain, but the algorithm allows
communication with other agents to reach and settle at the
target destination. These algorithms are robust to external
disturbances since the algorithms themselves use the swarm
density distribution.
Velocity field generation for swarm control is a nonoptimal decentralized control method for swarms that synthesizes smooth velocity fields as a function of time and position
[14]. With a designated target distribution, the swarm follows
the velocity field using the heat equation to convergence by
using local agent position information to estimate its local
density. The advantage of this method is that the agents
facilitate collision avoidance and move in a smoother manner
than previous Markov chain approaches [10]. Unfortunately,
the use of the heat equation diffuses the agents in a locally
uniform manner to the target density. Therefore, this is a
non-optimal method of controlling the swarm to a target
distribution.
Model predictive control (MPC) has been heavily studied
for nonlinear systems and for applications including spacecraft maneuvering and attitude control [15], [16], [17]. Decentralized MPC has also been studied for thousands of low
cost spacecraft with limited capabilities [18]. MPC computes
the control input by solving a finite horizon problem. This

allows it to optimize the current time step while taking the
future time steps into account.
One last decentralized approach to controlling swarms
is using sequential convex programming [19]. Sequential
convex programming uses multiple iterations to maintain
accuracy of the convex approximations of the constraints
which create more efficient trajectories. Uncertainties in the
trajectories are accounted for when the algorithm is tuned.
Thus, this provides the robustness for the swarm while
pushing the agents to the designated targets while using
this algorithm in combination with MPC in real time. Morgan also used sequential convex programming to do target
assignment (mapping of agents to targets) and trajectory
generation for varying swarm sizes through time [20]. This
method is viable for swarms, but through RFS control, target
assignment is not necessary.
II. R ANDOM F INITE S ET C ONTROL P ROBLEM
F ORMULATION

(1)

The control
sequence is also defined by a RFS in the form

Ut = u1,t , u2,t , · · · , uNtotal(t) ,t and a RFS probability
density given by
p(U = {u1 , u2 , · · · , un })
= p(|U |= n)p({u1 , u2 , · · · , un } | |U | = n),

where b(x), ps (ζ), and f (x|ζ) are the agents’ birth, survival
probability, and target motion models respectively and ζ is
the previous state [22]. The bar on v̄t (x) denotes that the
PHD is conditioned on measurements Zt−1 . For the time
update, the equation is given by
vt (x) = (1 − pd (x))v̄t (x)
X
pd (x)g(zt |x)v̄t (x)
R
,
+
c(z)
+
pd (ζ)g(zt |ζ)v̄t (ζ)dζ
z ∈Z

(4)

t

The control problem for swarming agents is set up by
using RFS theory [21], [22]. This theory addresses the
decentralized formulation for each agent in the formation.
Each agent has the challenge of estimating its local formation
configuration and designing a decentralized control policy
to achieve some local configuration. It is assumed that
each agent within the swarm is identical and using unique
identifiers on each agent is unnecessary. Using RFS theory,
the number of agents and their states is determined from
measurements. The number of agents in the field is denoted
by Ntotal (t) and
 may be varying as agents enter and leave the
field. Xt = x1,t , x2,t , · · · , xNtotal(t) ,t denotes a realization
of the RFS distribution for agents. At any time, t, the RFS
probability density function can be written as
p(X = {x1 , x2 , · · · , xn })
= p(|X|= n)p({x1 , x2 , · · · , xn } | |X| = n).

the cardinality of the agent set that is Poisson distributed
[22]. It is noted that the assumptions made by the PHD filter
are strong assumption for swarming robotics. However, this
is a good starting point for initial proof-of-concept study.
The PHD recursion for a general intensity function, vt (x),
is given by
Z
v̄t (x) = b(x) + ps (ζ)f (x|ζ)v(ζ)dζ,
(3)

(2)

since the number of the agents in the field to be controlled
is also varying. The Bayesian multi-target state estimation
problem generally has no closed-form solution, but when
the form p(|X|= n) and p({x1 , x2 , ..., xn } | |X| = n) are
assumed, approximate solutions can be determined by the
Probability Hypothesis Density (PHD) filter. Using the PHD
filter assumption, the RFS is represented by its first moment (intensity function ν(x, t)) [21]. ν(x, t) represents the
probability of finding an agent in a region of state space.
The PHD is defined over the single agent’s state space and
has the property that the expected number of agents over a
region is simply the integral of the PHD over that region.
Propagation of the PHD can be determined if the agents are
assumed to be independent and identically distributed with

where pd (x), g(zt |x), and c(z) are the probability of detection, likelihood function, and clutter model of the sensor
respectively [22]. By using this recursion, the swarm probabilistic description can be updated. Unfortunately, Eqs. (3)
and (4) do not contain a closed-form solution and the numerical integration suffers from higher computational time as the
states increase. Fortunately, a closed-form solution exists if
it is assumed that the survival and detection probabilities are
state independent, and the intensities of the birth and spawn
RFSs are Gaussian mixtures of the form
ν(x, t) = f (x) =

ν̄(x, t) = g(x) =

Nf
X

(5)

i=1


(i)
wf N x; mif , Pfi ,

Ng
X


wg(i) N x; mig , Pgi ,

(6)

i=1


where w are the weights and N x; mi , P i is the probability
density function of a multivariate Gaussian distribution with
a mean and covariance respectively. Note that closed form
solutions using Gaussian mixtures exist for cases without
PNf (i)
the state independent assumption. Additionally, i=1
wf =
PNg (i)
= N̄total (t) where N̄total (t) is the
Ntotal (t) and
i=1 wg
desired number of agents. The intensity function ν(x, t) is
in terms of the swarm state while ν̄(x, t) is in terms of the
desired state. The swarm intensity function can be propagated
through updates on the mean and covariance as given by
mif,k+1 = Ak mif,k + Bk uif,k ,

(7)

i
i
Pf,k+1
= Ak Pf,k
ATk + Qk ,

(8)

where Qk is process noise. Then given the Gaussian mixture
intensities assumption, a control variable is calculated for
each component uif,k . Additionally, each Gaussian mixture
component can represent many agents since the intensity
function integrates to the total number of agents. Note
that although linear dynamics is used, the RFS-based cost
function is nonquadratic.

Each individual swarm agent will run a local PHD observer to estimate the state of the swarm by modeling
PSfrag replacements
the swarm as a distribution. Thus, using RFS PSfrag
theory, replacements
it is
assumed that the individual swarm agents form an intensity
function that is a Gaussian mixture density in which the
y
means and covariances of the Gaussian mixture are propx
agated and controlled. An optimal control problem is set
(a) Cauchy-Schwarz Divergence
up that tracks a desired swarm formation by minimizing its
control effort in the following objective function
Z T
J(u) =
E {u⊺ Ru} + D(ν(x, t), ν̄(x, t))dt, (9) PSfrag replacements
2.5

400

2

200

1.5

0

1

-200

0.5
60

-400
60

50

40

50

40

40

40

30

20

0

30

20

20

20

y

10

0

0

10
0

x

(b) L22 Distance

1500

1000

500

0

-500

0

where ν̄(x, t) is the desired formation and u is the control
effort for the Gaussian mixture swarm. Both ν(x, t) and
ν̄(x, t) are defined over the complete state space which
include position and velocity parameters. D(ν(x, t), ν̄(x, t))
is the distance between these Gaussian mixtures which has
several closed-form solutions.
III. D ISTRIBUTIONAL D ISTANCE BASED -C OST
A. Cauchy-Schwarz Divergence
The Cauchy-Schwarz divergence is based on the CauchySchwarz inequality for inner products, and it is defined for
two random vectors with probability densities f and g given
by


hf, gi
,
(10)
DCS (f, g) = ln
kf kkgk
where h·, ·i is the L22 inner product over the densities. The
argument of the logarithm is non-negative since probability
densities are non-negative, and it does not exceed one
by the Cauchy-Schwarz inequality. The Cauchy-Schwarz
divergence can be interpreted as an approximation to the
Kullback-Leibler divergence but has a closed-form expression for Gaussian mixtures [24]. This is useful for calculating
the distance between two-point processes represented by intensity functions. The Cauchy-Schwarz divergence between
two Poisson point processes with Gaussian mixture intensities is shown in Eq. (10) when Eqs. (5) and (6) are substituted
inside. Note that in the control formulation used, only ν(x, t)
is assumed to depend on the control u. Therefore, the term
that depends only on ν̄(x, t) is omitted from the objective
function since ν̄(x, t) does not depend on u.
Figure 1(a) shows the surface plot using the CauchySchwarz divergence for 4 densities in the swarm at an initial
time instance. Each initial state has hills while the destination
states have valleys. Thus, each density will repel each other
while moving towards the final state through time. If the
initial state is too large compared to the destination states,
it will take longer for the 4 densities to converge to the
destination values. Also, the repelling effect due to the hills
are relatively small. Thus, the Cauchy-Schwarz divergence
may not be the fastest converging solution for the objective
function minimization.
B. L22 Distance
Alternatively, the distance between two Poisson point processes with Gaussian mixture intensities can be determined

60

40
20

y

0

0

10

x

20

30

40

50

(c) L22 + Quadratic Term

Fig. 1: Surface Plots of the RFS divergences.
by using the L22 distance between the probability densities.
The L22 is given by
Z
2
DL22 (f, g) = (f − g) dx.
(11)
Figure 1(b) shows surface plot using the L22 distance for
a 4-density swarm at an initial time instance. Each initial
state has more defined hills compared to the Cauchy-Schwarz
divergence. Thus, the initial states have a stronger repel
effect from each other. Also, the destination states have
large valleys that create a large attraction effect for each
initial state to move to. Thus, the optimization solution will
be faster in the L22 distance case. Unfortunately, the L22
distance suffers from a similar issue to the Cauchy-Schwarz
divergence. If the initial conditions increase farther away
from the destination states, the optimization may take much
longer due to a flat surface away from the destination states.
C. L22 Distance with Quadratic Term
The issue of convergence remains for the L22 distance when
the initial states are far larger than the destination states. To
achieve faster convergence, an additional term is added to the
L22 distance to shape the gradient descent through a quadratic
term as given by
DL22 mod (f, g) = DL22 (f, g)


Ng Nf
X
X
(12)
(i)
wg(j) wf N (mjg ; mif , Pgi + Pfj ).
− ln 
j=1 i=1

Note that this term is referred as quadratic, although it
may be more appropriate to call it quadratic-like. Figure
1(c) shows the surface plot using the modified L22 distance
for a 4-density swarm in Eq. (12). Compared to the L22
distance, the initial and destination states provide the hills
and valleys necessary to obtain convergence. However, as
the initial states move outwards, the surface map decreases
in a quadratic fashion instead of staying flat. Thus, the
optimization can occur at any point to reach convergence

in a timely manner. To solve the minimization problem for
the objective function, traditional LQR based solutions are
not applicable since the terms inside the objective functions
are nonquadratic [23]. The minimization of the objective
function in discrete time is
min

uk ,k =1,...,T

+

J(u) =

T
X

E{u⊺k Ruk }

k=1

Nf Nf
X
X

j
i
wf,k wf,k N (mjf,k ; mif,k , Pf,k
+ Pf,k
)

Ng Ng
X
X

j
i
wg,k wg,k N (mjg,k ; mig,k , Pg,k
+ Pg,k
)

(j)

(i)

(j)

(i)

j =1 i =1

+

j =1 i =1

−2

Ng Nf
X
X

(j)

V. R ESULTS

(i)

j
i
wg,k wf,k N (mjg,k ; mif,k , Pg,k
+ Pf,k
)

j =1 i =1



− ln 

which measures the distance of the agents to the means of
the densities. As the swarm evolved, this distance determined
which units belong to a given component. Agents were
controlled according to their placement in each density
through an open loop method using the MPC control input
obtained for each density. Although an open loop method
was used, feedback control can be used if the state estimates
are determined from the PHD filter.
The great advantage of using MPC is that the algorithm
can handle nonlinearities in the objective function. Unfortunately, numerically solving an optimization at each timestep is not particularly fast, so it can be difficult to run
optimization problems in real-time to control swarm robots
towards a target state.

Ng Nf
X
X
j=1 i=1



(j) (i)
j
i
wg,k wf,k N (mjg,k ; mif,k , Pg,k
+ Pf,k
) ,

(13)

Subject to : mif,k+1 = Ak mif,k + Bk uif,k ,

(14)

N

where uk = [(u1f,k )T , · · · , (uf,kf )⊺ ]⊺ is the collection of all
control variables.
IV. R ECEDING H ORIZON DYNAMIC P ROGRAMMING

In the simulation experiment, a 4-density swarm on a 2D plane was initialized in a square grid where the densities
1, 2, 3, and 4 were defined counterclockwise starting on
the first quadrant. With the 4-density swarm, four different
test cases were implemented to bring the densities to the
target trajectories and to test the control theory involved
from the control formulation. The first test case shows the
limitations of using only the L22 distance with varying initial
conditions in a square grid and four target destinations
located at (±1,±1). The last three cases used the L22 distance
with a quadratic term and varying target destinations. In
Case 2, four target destinations are located at (±1,±1).
For Case 3, three target destinations are located at (±1,1)
and (-1,-1). Lastly in Case 4, five target destinations are
located at (±1,±1) and (0,0). The receding horizon dynamic
programming was used to determine the control effort and
trajectories for the simulation, but other methods including
iterative LQR can be used too [23].

m2,y

m3,y

m4,y

m1,y

m2,x

m3,x

m4,x

m1,x

An optimal solution based on minimizing the objective
uses MPC or receding horizon control. This uses dynamic
programming techniques to obtain an optimal control solution u from the objective function and has been implemented
on swarm control previously [20].
Conceptually, at a time t, the knowledge of the system
model is used to derive a sequence u(t|t), u(t + 1|t), u(t +
4
4
2|t), · · · , u(t + Tp |t) where Tp is the finite horizon from
2
2
the current state x̄(k) [26]. With the input sequence, the
0
0
0
0.1
0.2
0.3
0.4
0
0.1
0.2
0.3
0.4
state is moved forward in time by the control horizon,
Tc ;
PSfrag replacements
Time (s)
Time (s)
0
4
usually one time-step. Then the same strategy is repeated
-2
2
-4
0
for time t + 1. Tp can be chosen to be either small or large.
0
0.1
0.2
0.3
0.4
0
0.1
0.2
0.3
0.4
As Tp increases, the degrees of freedom in the optimization
Time (s)
Time (s)
0
0
increases which can slow down the algorithm considerably
-2
-2
-4
-4
even though more of the future reference trajectory would
0
0.1
0.2
0.3
0.4
0
0.1
0.2
0.3
0.4
Time (s)
Time (s)
be useful to bring the output closer to the reference. With
4
0
2
-2
a smaller Tp , the computation time will be faster, but the
0
-4
0
0.1
0.2
0.3
0.4
0
0.1
0.2
0.3
0.4
optimization may be more suboptimal. Thus, the swarm may
Time (s)
Time (s)
not converge to the desired configuration.
For the RFS control formulation, a u that controls the
swarm densities through their statistics (mean and covari- Fig. 2: The position responses of density statistics (means)
ance) was found by minimizing the objective as given by are shown for Case 1 and Case 2.2 Black-dashed and redEq. (13) and (14). This was done using MATLAB’s fminunc dotted lines are the responses for L2 distance with (±3,±3))
The bluesolver and the Quasi-Newton algorithm [27]. Then, the and (±1.5,±1.5)) initial conditions respectively.
2
solid
lines
are
the
responses
for
the
L
+Quadratic
term with
2
agents in the swarm were initialized to the closest density
initial
conditions
of
(±3,±3)).
The
green-dotted
lines
are the
using the Mahalanobis distance given by
target
destinations
for
both
cases.
q
DM (x, mif,k ) = (x − mif,k )⊺ P −1 (x − mif,k ), (15)

A. Case 1: L22

5

time = 0.00 s

5

time = 0.05 s

0.08

y

y

y

y

y
y

y

y

For Case 1, four swarm densities are controlled to move
0.07
0
0
towards the target destinations at initial conditions farther
0.06
away (square grid at (±3,±3)) and closer to (square grid
-5
-5
0.05
at (±1.5,±1.5)) the target destinations as shown by mean
-5
0
5
-5
0
5
x
x
PSfrag
replacements
responses given by the black-dashed and red-dotted
lines
in
0.04
time = 0.10 s
time = 0.40 s
5
5
Figure 2 respectively. From the trajectory snapshots given
0.03
by Figure 3(a), initial conditions that are far from the target
0
0
destinations do not have a converging control solution. From
0.02
the surface visualization in Figure 1(b), the general plane
0.01
-5
-5
is flat in areas away from the target destinations and states
-5
0
5
-5
0
5
x
x
of the density. Therefore, optimization using MPC is more
difficult in these flat areas and may not converge to a solution.
Fig. 4: Case 2: Trajectories of the swarm with individual
If the densities are initialized much closer to the target
agents (red dots) and four target destinations (black x’s).
destinations as shown in Figure 3(b), the flatness in the
general plane is minimal, and the optimization step in MPC
converges to a solution. By using the L22 distance, converging C. Case 3: L2 +Quadratic Term, Three Target Destinations
2
control solutions can only be found for initial conditions and
Case
3
illustrates
the effect of three target destinations
target destinations that are close.
on the final trajectories of the four swarm densities. From
Figure 5, it can be visually shown where the swarm densities
time = 0.00 s
time = 0.40 s
5
5
are located relative to the target destinations at each time
0.08
step. The trajectories for density 1 and density 3 reach their
0.07
0
0
target, but densities 2 and 4 reach the third target with
0.06
approximately 0.30 of steady-state error. The results obtained
-5
-5
follow directly from the RFS control theory using the L22 plus
Sfrag replacements
0.05
-5
0
5
-5
0
5
x (a) Initial conditions of (±3,±3) x
quadratic distance term. By using this L22 with a quadratic
0.04
time = 0.00 s
time = 0.40 s
5
5
term in the objective function, the individual densities will
0.03
attract towards the target destinations while repulsing away
0
0
from each other. This can be shown in the surface map shown
0.02
in Figure 1(c) where the hills are areas of repulsion and
0.01
-5
-5
valleys are areas of attraction. Thus, for Case 3, densities 2
-5
0
5
-5
0
5
time = 0.10 s
x (b) Initial conditions of (±1.5,±1.5) x
and 4 are attracted to the same target, but they stay away
from each other which causes both to share the same target
Fig. 3: Case 1: Trajectories of the swarm with individual
at a distance away.
agents (red dots) and four target destinations (black x’s).
5

time = 0.00 s

5

time = 0.05 s
0.12

B. Case 2: L22 +Quadratic Term, Four Target Destinations
y

y

y

y

0
0
0.1
In Case 2, four swarm densities move towards the four
target destinations given by the blue-solid lines (mean re0.08
-5
-5
sponses) in Figure 2. Figure 4 shows the trajectory snapshots
-5
0
5
-5
0
5
x
x
and final states of each of the swarm densitiesPSfrag
duringreplacements
the
0.06
time = 0.10 s
time = 0.40 s
5
5
simulation. The target destinations are plotted as black x’s.
0.04
The red dots are the individual swarm agents that form
0
0
the Gaussian mixture densities. From the figure, all four
0.02
densities converge to the target destinations in approximately
-5
-5
0.17 seconds and approximately 0.03 of steady-state error
-5
0
5
-5
0
5
x
x
between the densities’ position to the target destinations.
In comparison to Case 1, Figure 2 shows that for small Fig. 5: Case 3: Trajectories for the swarm with individual
distances between the initial state and the target destination, agents (red dots) and three target destinations (black x’s).
the L22 distance is sufficient for state convergence, but as the
2
distance increases, the L2 distance diverges away. By adding D. Case 4: L2 +Quadratic term, Five Target Destinations

2

the quadratic term to L22 , the optimization step can directly
determine the minimum for the control solution shown in
Figure 1(c). Therefore, the target destinations attract the
swarm densities at distances that fail for only L22 distance
given by Figure 2.

For Case 4, four swarm densities are attracted to five
target destinations. The trajectory snapshots of the densities
are visually shown relative to the five target destinations
in Figure 6. All densities have a steady state error of
approximately 0.17 which follow the theory as expected.

Since the swarm densities are far from each other, the effects
of repulsion are minimal. Also, the densities are attracted
to the four target destinations that make up a square, but
they are also attracted to the target destination at the origin.
This is due to the minimization of the objective function
that has a L22 and a quadratic term where the individual
densities will attract towards the target destinations. Since
there is an additional target destination of attraction at the
origin, all four swarm densities are affected by the origin
as they are moving towards the 4 square target destinations.
Thus, compared to Case 2 with only 4 target destinations,
the densities in this case will have steady state error due to
attraction of the target destination that has no density.
time = 0.00 s

5

5

time = 0.05 s

0.09

y

y

0.08

0

0
0.07
0.06

-5
-5

Sfrag replacements
5

0

-5
-5

5

x
time = 0.10 s

5

0

5
0.05

x
time = 0.40 s

0.04

y

y

0.03

0

0
0.02
0.01

-5
-5

0

5

-5
-5

x

0

5

x

Fig. 6: Case 4: Trajectories for the swarm with individual
agents (red dots) and five target destinations (black x’s).
VI. C ONCLUSIONS
The objective of this paper was to formulate a method to
control swarming formations by MPC and RFS theory. By
setting up the problem using information divergence to define
the distance between the swarm RFS and a desired distribution, an optimal control problem was found that tracked a
linear system with a nonquadratic objective function using
MPC that approximates the optimal control solution that
minimizes the objective function.
ACKNOWLEDGMENT
The authors wish to acknowledge support by the National
Aeronautics and Space Administration under Contract Number NNX16CP45P issued through the NASA STTR Program
by the Jet Propulsion Laboratory (JPL) and led by Amir
Rahmani at JPL. The authors wish to acknowledge useful
conversations related to satellite technologies with Chuck
Hisamoto, Vaughn Weirens, and Suneel Sheikh of ASTER
Labs, Inc.
R EFERENCES
[1] Kube, C. Ronald, and Hong Zhang. “Collective robotics: From social
insects to robots.” Adaptive behavior 2.2 (1993): 189-218.
[2] Vassev, Emil, Mike Hinchey, and Joey Paquet. “Towards an ASSL
specification model for NASA swarm-based exploration missions.”
Proceedings of the 2008 ACM symposium on Applied computing.
ACM, 2008.
[3] Ryan, Allison, et al. “An overview of emerging results in cooperative
UAV control.” Decision and Control, 2004. CDC. 43rd IEEE Conference on. Vol. 1. IEEE, 2004.

[4] Rubenstein, Michael, Alejandro Cornejo, and Radhika Nagpal.
“Programmable self-assembly in a thousand-robot swarm.” Science
345.6198 (2014): 795-799.
[5] Lunze, Jan. Feedback control of large-scale systems. London: PrenticeHall, (1992): 9-11.
[6] Pace, M., M. Birattari, and M. Dorigo. “The swarm/potential model:
Modeling robotics swarms with measure-valued recursions associated
to random finite sets.” IEEE Transactions on Robotics, page submitted
(2013).
[7] Centralized Control. Ian Sommerville , 2008. Web. 05 May 2017.
[8] Mondada, Francesco, et al. “The cooperation of swarm-bots: Physical
interactions in collective robotics.” IEEE Robotics & Automation Magazine 12.2 (2005): 21-28.
[9] Bandyopadhyay, Saptarshi, Soon-Jo Chung, and Fred Y. Hadaegh.
“Probabilistic swarm guidance using optimal transport.” Control Applications (CCA), 2014 IEEE Conference on. IEEE, 2014.
[10] Acikmese, Behet, and David S. Bayard. “A markov chain approach to
probabilistic swarm guidance.” American Control Conference (ACC),
2012. IEEE, 2012.
[11] Chattopadhyay, Ishanu, and Asok Ray. “Supervised self-organization
of homogeneous swarms using ergodic projections of markov chains.”
IEEE Transactions on Systems, Man, and Cybernetics, Part B (Cybernetics) 39.6 (2009): 1505-1515.
[12] Bandyopadhyay, Saptarshi, Soon-Jo Chung, and Fred Y. Hadaegh. “Inhomogeneous Markov chain approach to probabilistic swarm guidance
algorithm.” 5th Int. Conf. Spacecraft Formation Flying Missions and
Technologies. 2013.
[13] Hadaegh, Fred Y., Soon-Jo Chung, and Harish M. Manohara. “On
development of 100-gram-class spacecraft for swarm applications.”
IEEE Systems Journal 10.2 (2016): 673-684.
[14] Eren, Utku, and Behet Açıkmeşe. “Velocity field generation for density
control of swarms using heat equation and smoothing kernels.” IFACPapersOnLine 50.1 (2017): 9405-9411.
[15] Camacho, Eduardo F., and Carlos Bordons Alba. Model predictive
control. Springer Science & Business Media, 2013.
[16] Di Cairano, S., H. Park, and I. Kolmanovsky. “Model predictive
control approach for guidance of spacecraft rendezvous and proximity
maneuvering.” International Journal of Robust and Nonlinear Control
22.12 (2012): 1398-1427.
[17] Manikonda, V., et al. “A model predictive control-based approach for
spacecraft formation keeping and attitude control.” American Control
Conference, 1999. Proceedings of the 1999. Vol. 6. IEEE, 1999.
[18] Morgan, Daniel, Soon-Jo Chung, and Fred Y. Hadaegh. “Decentralized
model predictive control of swarms of spacecraft using sequential
convex programming.” Advances in the Astronautical Sciences 148
(2013): 1-20.
[19] Morgan, Daniel, Soon-Jo Chung, and Fred Hadaegh. “Spacecraft
swarm guidance using a sequence of decentralized convex optimizations.” (2012): Art-No.
[20] Morgan, Daniel, Soon-Jo Chung, and Fred Y. Hadaegh. “Swarm assignment and trajectory optimization using variable-swarm, distributed
auction assignment and model predictive control.” AIAA guidance,
navigation, and control conference. 2015.
[21] Mahler, Ronald PS. “Multitarget Bayes filtering via first-order multitarget moments.” IEEE Transactions on Aerospace and Electronic
systems 39.4 (2003): 1152-1178.
[22] Vo, B-N., and W-K. Ma. “The Gaussian mixture probability hypothesis
density filter.” IEEE Transactions on signal processing 54.11 (2006):
4091-4104.
[23] Todorov, Emanuel, and Weiwei Li. “A generalized iterative LQG
method for locally-optimal feedback control of constrained nonlinear
stochastic systems.” American Control Conference, 2005. Proceedings
of the 2005. IEEE, 2005.
[24] Hoang, Hung Gia, et al. “The Cauchy Schwarz divergence for Poisson point processes.” IEEE Transactions on Information Theory 61.8
(2015): 4475-4485.
[25] Inaba, Masayuki, and Peter I. Corke. “Linear-quadratic control and
smoothing.” Robotics Research: The 16th International Symposium
ISRR. Switzerland: Springer, 2016. 43-47. Print.
[26] Findeisen, Rolf, and Frank Allgwer. “An introduction to nonlinear
model predictive control.” 21st Benelux Meeting on Systems and
Control. Vol. 11. Eindhoven, The Netherlands: Technische Universiteit
Eindhoven Veldhoven, 2002.
[27] Find minimum of unconstrained multivariable function - MATLAB
fminunc, MathWorks. Web.

