ABSTRACT
Title of dissertation: Optimality, Synthesis and a
Continuum Model for Colletive Motion
Udit Halder, Dotor of Philosophy, 2019
Dissertation direted by: Professor P. S. Krishnaprasad
Dept. of Eletrial & Computer Engineering
It is of importane to study biologial olletives and apply the wisdom so
arued to modern day engineering problems. In this dissertation we attempt to
gain insight into olletive behavior where the main ontribution is twofold. First,
a `bottom-up' approah is employed to study individual level ontrol law synthesis
and emergene thereby of olletive behavior. Three dierent problems, involv-
ing single and multiple agents, are studied by both analytial and experimental
means. These problems arise from either a pratial viewpoint or from attempts
at desribing biologially plausible feedbak mehanisms. One result obtained in
this ontext for a double agent senario is that under a partiular onstant bearing
pursuit strategy, the problem exhibits ertain features ommon with the Kepler
two body problem. Laboratory demonstrations of the solutions to these problems
are presented. It is to be noted that these types of individual level ontrol prob-
lems an help understand and onstrut building bloks for group level behaviors.
The seond approah is `top-down' in nature. It treats a olletive as a whole
and asks if its movement minimizes some kind of energy funtional. A key goal
of this work is to develop wave equations and their solutions for a natural lass
of optimal ontrol problems with whih one an analyze information transfer in
oks. Controllability arguments in innite dimensional spaes give strong sup-
port to onstrut solutions for suh optimal ontrol problems. Sine the optimal
ontrol problems are innite dimensional in the state spae and one annot simply
expet Pontryagin's Maximum Priniple (PMP) to apply in suh a setting, the
work has required are and attention to funtional analyti onsiderations. In this
work, it is shown that under a ertain assumption on nite o-dimensionality of a
reahable set, PMP remains valid. This assumption is then shown to hold true for
the ase of a spei ensemble of agents, eah with state spae as the Heisenberg
group H(3). Moreover, analysis of optimal ontrols demonstrates the existene
of traveling wave solutions in that setting. Synhronization results are obtained
in a high oupling limit where deviation from neighbors is too ostly for every
agent. The ombination of approahes based on PMP and alulus of variations
have been fruitful in developing a solid new understanding of wave phenomena in
olletives. We provide partial results along these lines for the ase of a ontinuum
of planar agents (SE(2) ase).
Finally, a dierent top-down and data-driven approah to analyze olletive be-
havior is also put forward in this thesis. It is known that the total kineti energy
of a ok an be divided into several modes attributed to rigid-body translations,
rotations, volume hanges, et. Flight reordings of multiple events of European
starling oks yield time-signals of these dierent energy modes. This approah
then seeks an explanation of kineti energy mode distributions (viewed as ok-
sale deisions) by appealing to tehniques from evolutionary game theory and
optimal ontrol theory. We propose the notion of ognitive ost that alulates
a suitably dened ation funtional and measures the ost to an event, resulting
from temporal variations of energy mode distributions.
Optimality, Synthesis and a
Continuum Model for Colletive
Motion
by
Udit Halder
Dissertation submitted to the Faulty of the Graduate Shool of the
University of Maryland, College Park in partial fulllment
of the requirements for the degree of
Dotor of Philosophy
2019
Advisory Committee:
Professor P. S. Krishnaprasad, Chair/Advisor
Professor Steve Marus
Professor Nuno Martins
Professor André Tits
Professor Nikhil Chopra, Dean's Representative
© Copyright by
Udit Halder
2019
Dediated to my parents
ii
Aknowledgments
I would like to take this opportunity to thank my advisor Prof. P. S. Krish-
naprasad for all the help and support he has generously provided over the past
ve years. He has given me the freedom to think deeply about hallenging and
fundamental problems. It is his onstant enouragement and enthusiasm for si-
ene over ountless hours of enlightening disussions that positively propelled the
work presented in this thesis. His persistent emphasis on attention to detail in
every aspet of researh has helped me grow as a researher. It has been an ab-
solute pleasure having him as my thesis advisor. I would also like to thank Prof.
Steve Marus, Prof. Nuno Martins, Prof. André Tits and Prof. Nikhil Chopra for
being in my dissertation ommittee and reviewing this doument. The insightful
omments from the members of the dissertation ommittee have helped improve
the quality of this thesis. Collaboration with Dr. Uroš Kalabi¢ at the Mitsubishi
Eletri Researh Laboratory has been very helpful. I also ollaborated with Dr.
Eri Justh on a entral problem of my thesis and I thank him for all his feedbak
and wise opinions.
Over the past ve years, I have been fortunate to interat with my olleagues in
the Intelligent Servosystems Laboratory. I thank Ms. Vidya Raju, Dr. Biswadip
Dey, Dr. Yunlong Huang, Prof. Kevin Galloway, Mr. Brent Shlotfeldt, and Mr.
Kenneth Miltenberger for all the stimulating disussions we have had. I enjoyed
a lot while ollaborating with Dr. Biswadip Dey, Mr. Brent Shlotfeldt, and Ms.
Vidya Raju at several stages of my PhD. Speial thanks to Vidya who has been
iii
on a same boat as me from the start of my graduate study. I would also like to
thank my friends outside the lab without whose ompany and support my life as
a graduate student would have been inomplete. Biswadip, Vidya, Dipankar, Sid-
dharth, Soham, Debdipta, Abhishek, Rajdeep, Ankit, Proloy, Soumyadip, Agniv,
Jiaul, Kushal and others  thank you all for making this journey a memorable
one.
Finally, I would like to thank my parents who have been a soure of onstant
inspiration and enouragement. I may not have aomplished this without their
support.
This researh was supported by the Air Fore Oe of Sienti Researh
under AFOSR Grant FA9550-10-1-0250, an AFOSR FY2012 DURIP Grant No.
FA2386-12-1-3002, the ARL/ARO MURI Program Grant No. W911NF-13-1-
0390, through the University of California Davis (as prime), the ARL/ARO Grant
No. W911NF-17-1-0156, through the Virginia Polytehni Institute and State
University (as prime), and by Northrop Grumman Corporation.
iv
Table of Contents
List of Tables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . vii
List of Figures. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .viii
List of Abbreviations and Notations.. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . ix
1 Introdution . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1
1.1 Mathematial Bakground . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5
1.1.1 Self Steering Partile Model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5
1.2 Experimental Setup. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7
I Synthesis of Colletive Motion: Bottom-up Approahes 9
2 Feedbak Laws for Colletive Motion .. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
2.1 Introdution . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
2.2 Steering for Beaon Pursuit under Limited Sensing . . . . . . . . . . . . . . . . . . . 11
2.2.1 Traking a Moving Beaon . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13
2.2.2 Dynamis Restrited to the Invariant Manifold . . . . . . . . . . . . . . . 15
π
2.2.3 Speial Case: α = . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17
2
2.2.4 The Limited Field of View (FOV) Problem . . . . . . . . . . . . . . . . . . . 23
2.2.5 Implementation . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 26
2.2.6 Assoiated Lagrangian . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30
2.3 Biomimeti Algorithms for Coordinated Motion . . . . . . . . . . . . . . . . . . . . . . 33
2.3.1 Mutual Motion Camouage (MMC) . . . . . . . . . . . . . . . . . . . . . . . . . . . 34
2.3.2 Topologial Veloity Alignment (TVA) . . . . . . . . . . . . . . . . . . . . . . . . 36
2.3.3 Implementation Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 38
3 Optimal Steering of Agents on a Plane . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 45
3.1 Introdution . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 45
3.2 Optimal Steering of a Uniyle . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 47
3.2.1 Optimal Control Solution . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49
3.2.2 Charaterizing the Types of Motion . . . . . . . . . . . . . . . . . . . . . . . . . . . 51
3.2.3 On Time-optimality . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59
3.3 Optimal Control of a Colletive of Agents . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 60
3.4 Conluding Remarks. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 64
v
II Analysis of Colletive Motion: Top-down Approahes 65
4 Continuum Floking and Control . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66
4.1 Motivation.. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66
4.2 A Control System on a Loop Group.. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 69
4.3 Controllability.. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 71
4.4 Optimal Control Problems . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 76
4.4.1 Calulus of Variations Approah . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 76
4.4.2 Maximum Priniple Approah . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 79
4.5 Speial Case : G = H(3) . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 84
4.5.1 Controllability . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 85
4.5.2 Equations of Optimal Control . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 87
4.5.3 Behavior of Optimal Control . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 92
4.5.4 Strong Coupling Limit, χ → ∞ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 98
4.5.5 Simulation Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 101
4.6 A Continuum of Agents on the Plane . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 105
4.6.1 Equations of Optimal Control . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 106
4.6.2 Strong Coupling Limit, χ → ∞ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 110
4.6.3 Simulation Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 111
4.7 Disussion and Sope of Future Researh . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 118
5 Cognitive Cost of Floking: A Geometri and Hamiltonian Per-
spetive . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .120
5.1 Introdution . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 120
5.2 Floking Data . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 122
5.3 Data Smoothing . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 123
5.3.1 A linear generative model . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 124
5.3.2 Data smoothing in the Eulidean setting . . . . . . . . . . . . . . . . . . . . . . 125
5.4 Energy Modes . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 127
5.5 Generative model on the 1-simplex and the data-smoothing problem131
5.6 Data Fitting Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 136
5.7 Disussion . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 138
6 Conlusions and Diretions for Future Researh . . . . . . . . . . . . . . . . . .148
A An Optimal Control Problem in an Innite Dimensional Setting152
A.1 Introdution . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 152
A.2 Maximum Priniple .. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 156
Bibliography.. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .167
vi
List of Tables
5.1 Details of aptured oking events . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 123
5.2 Hamiltonian Signature . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 138
vii
List of Figures
1.1 Mobile robot based experimental platform (Pioneer 3 DX) with
two-wheel dierential and aster. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7
2.1 Illustration of salar shape variables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13
2.2 Phase portrait . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18
2.3 Nulllines of (2.5) with u = 1, ν1 = 0.5, ν2 = 1. . . . . . . . . . . . . . . . . . . . . . . . . 20
2.4 Robot (Pioneer 3 DX) with Kinet mounted, and the orange one
used as the beaon. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 27
2.5 Implementation results for the limited FOV problem .. . . . . . . . . . . . . . . . 29
2.6 Implementation results for MMC ontrol law. . . . . . . . . . . . . . . . . . . . . . . . . . 40
2.7 Implementation results for TVA ontrol law .. . . . . . . . . . . . . . . . . . . . . . . . . . 41
3.1 Intersetion of Hamiltonian and Casimir surfaes . . . . . . . . . . . . . . . . . . . . . 51
3.2 Sample trajetories for dierent values of c, η and s1, s2 . . . . . . . . . . . . . . 53
4.1 Numerial solution of (4.61) for experiment 1 . . . . . . . . . . . . . . . . . . . . . . . . . 99
4.2 Evolution of x3 for experiment 1 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 100
4.3 Numerial solution of (4.61) for experiment 2 . . . . . . . . . . . . . . . . . . . . . . . . . 101
4.4 Evolution of x3 for experiment 2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 102
4.5 Numerial solutions of (4.61) for experiment 3 . . . . . . . . . . . . . . . . . . . . . . . . 103
4.6 Numerial solution of (4.99) and state evolution for experiment 1. . . 112
4.7 State evolution for experiment 2. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 114
4.8 State evolution for experiment 3. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 115
4.9 Numerial solution of (4.99) and state evolution for experiment 4. . . 116
4.10 Numerial solution of (4.99) and state evolution for experiment 5. . . 117
5.1 Hamiltonian signatures . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 136
5.2 Optimal energy alloation for event 1 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 140
5.3 Optimal energy alloation for event 2 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 141
5.4 Optimal energy alloation for event 3 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 142
5.5 Optimal energy alloation for event 4 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 143
5.6 Optimal energy alloation for event 5 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 144
5.7 Optimal energy alloation for event 6 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 145
5.8 Optimal energy alloation for event 7 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 146
5.9 Optimal energy alloation for event 8 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 147
viii
List of Abbreviations and Notations
ROS Robot Operating System
CB Constant Bearing pursuit
FOV Field-of-View
MMC Mutual Motion Camouage
TVA Topologial Veloity Alignment
PMP Pontryagin's Maximum Priniple
ri Position of agent i
{xi,yi} (planar) Natural Frenet Frame for agent i
tr Trae (of a matrix)
1n n× n identity matrix
am z Amplitude z
sn z Ellipti sine z, Jaobian ellipti funtion
n z Ellipti osine z, Jaobian ellipti funtion
dn z Delta amplitude z, Jaobian ellipti funtion
TmM Tangent spae at m to manifold M
G A matrix Lie group
g Lie algebra of G
S1 Cirle group
G ∞ 1Loop group C (S ;G)
L ∞ 1Loop algebra C (S ; g)
H(3) Heisenberg group
SE(2) Speial Eulidean group
h(3) Lie algebra of H(3)
se(2) Lie algebra of SE(2)
χ Coupling onstant
ix
Chapter 1
Introdution
The last few deades have witnessed an inrease in researh eorts towards un-
overing mehanisms behind olletive motion [Nagy et al., 2010; Ballerini et al.,
2008a; Cavagna et al., 2010; Inada and Kawahi, 2002℄ and pursuit behavior [Ol-
berg et al., 2000; Mizutani et al., 2003; Ghose et al., 2006; Chiu et al., 2010℄
in nature. Ranging from sh shools to bird oks, olletive behavior is seen
abundantly in nature. The onept of safety in numbers is used in aomplishing
variety of goals, from foraging food to evading predators. Reent improvements
in data olletion and proessing tehnology has enabled researhers to study
these natural oks in more detail than ever before [Ballerini et al., 2008a,b℄.
The driving question then beomes to answer how loal interations between in-
dividual agents in the olletive give rise to group level ohesion and synhrony.
Although several attempts have been made to understand these behaviors [Cuker
and Smale, 2007; Mora and Bialek, 2011; Bialek et al., 2012; Young et al., 2013;
Attanasi et al., 2014℄, the individual level mehanisms responsible for emergene
1
of olletive behavior remain mostly elusive to researhers. It is therefore a sig-
niant goal of this thesis to pursue suh questions.
This thesis is distintively divided into two parts, where we take two dier-
ent approahes to understand olletive behavior. The rst approah is alled a
`bottom-up' approah, i.e. instead of studying the ok as a whole, we onen-
trate on dynamis of individual agents and analyze simple interation laws among
small number of agents. Studying these interations are important sine they an
be used as a building blok for group level motion. In a 1995 paper [Visek et al.,
1995℄, a novel disrete time self-driven partile model was rst introdued to ad-
dress self-ordered motion in a system of partiles. The onept of self-steering
partiles was developed in the following deades [Justh and Krishnaprasad, 2004,
2006; Reddy et al., 2006; Mishiati and Krishnaprasad, 2010, 2012; Galloway
et al., 2013℄. We undertake the self-steering partile model under gyrosopi on-
trol [Justh and Krishnaprasad, 2003, 2004℄ as the basi model for individuals in
the ok. This model desribes a trajetory of an individual as a urve, desribed
by the natural Frenet frame equations [Bishop, 1975℄ in the Eulidean spae; and
the driving ontrols are given by speed and urvature of the urve. We show in
Chapter 2, 3 that even in the single agent or double agent ase, interesting mo-
tion patterns an be synthesized from arefully seleting these ontrol inputs. The
ontrol inputs an be generated from an underlying optimal ontrol problem or by
applying some biologially plausible feedbak strategies. Parallel to the quest of
2
mathemetial modeling, some groups in the robotis ommunity have performed
suessful implementation of various ontrol strategies [Thurrowgood et al., 2014;
Vásárhelyi et al., 2014℄, and thereby demonstrated the power of a bio-inspired ap-
proah towards synthesizing olletive motion. Our work is similar in spirit, and
provides indoor demonstrations of problems raised in Chapter 2. Some of these
problems were oneptualized from a pratial perspetive and arry engineering
value.
The other approah to study olletive behavior is what an be alled as `top-
down' view. Instead of speifying agent level ontrol laws, the idea is to infer
those laws from solving a bigger problem that investigates the ok as a whole.
Existing literature employs several methods suh as optimal ontrol [Justh and
Krishnaprasad, 2015b,a℄, statistial physis [Mora and Bialek, 2011; Bialek et al.,
2012℄ et. It is the framework of optimal ontrol [Justh and Krishnaprasad, 2015b℄
that we undertake and extend in this thesis. It has been observed from empir-
ial data [Ballerini et al., 2008a℄ that interation among starlings in the ok is
loal, i.e. eah bird interat with six/seven neighbors during ight. Taking in-
spiration from this idea, the entral onept of [Justh and Krishnaprasad, 2015b℄
is to set up an optimal ontrol problem whih penalizes ontrols of individual
agents oupled with mismath in ontrol with its `neighbors'. The neighbors are
determined by a previously dened interation graph. We then let the number
of agents in the ok to go to innity in order to propose a ontinuum model for
3
oking. Various ontinuum models have been studied for olletives [Kudrolli
et al., 2008; Topaz et al., 2006; Zhang et al., 2010℄. These models study a set of
partial dierential equations that desribe spatio-temporal evolution of the ok
density. Our approah is dierent in the sense that the system dynamis an be
seen as an ordinary dierential equation in an appropriate innite dimensional
Lie group setting. The oupling between birds are introdued through the mis-
math term in the ost funtional. A natural question of ontrollability of suh
a system is addressed by using a generalized Chow-Rashevsky theorem for in-
nite dimensional systems. This enables us to formulate the underlying optimal
ontrol problem in an innite dimensional setting in whih the usual Pontryagin's
maximum priniple fails in general without further assumptions. In Chapter 4,
we invoke a maximum priniple atered for this spei setting. A spei ex-
ample of ontinuum of nonholonomi integrators is also studied in detail. This
an be viewed as a ontinuum version of single agent Heisenberg ase [Justh and
Krishnaprasad, 2016℄. It has been found that optimal ontrol solutions possess
a traveling wave harater, whih might enable information transfer in the ok.
In addition to the Heisenberg ase, we provide optimal ontrol equations in the
ase of a ontinuum ok of planar agents. Synhronization results and numerial
simulations are presented for both the ases.
In Chapter 5, we present another `top-down' approah to the oking problem.
This approah is data-driven in nature. Kinemati energy modes of European
4
starling oks are represented on a simplex whih is then subjeted to desription
as trajetory of some evolutionary game dynamis. Solution of this data-tting
problem on the simplex results in ontrol inputs that are interpreted as modulation
of tness assoiated with the energy modes. We note that in ontrast to Chapter
4, where the ontrol inputs were individual agent-level (or `low-level') ontrols,
the ontrols obtained by this data-driven approah are ok-level (or `high-level')
ontrols. The ok is oneptualized to apply these ontrols to optimally alloate
its kineti energy among dierent modes.
1.1 Mathematial Bakground
1.1.1 Self Steering Partile Model
We desribe the partile model that is the underlying generative model in all our
subsequent analysis throughout this thesis. The trajetory of a single agent an
be desribed by a funtion r : [0, T ] → R3, for some T > 0. We assume r(t) to
be a regular urve, i.e. ṙ(t) =6 0, ∀t ∈ [0, T ]. Let s be the ar length parameter,
∫ t
i.e. s(t) = ‖ṙ(σ)‖ dσ. Under the regularity assumption, s(t) is monotonially
0
inreasing and invertible funtion of time. We an then reparametrize the urve
r(t) by the ar length parameter s and the evolution equations an be expressed
in terms of well known Fernet-Serret frames. However, this way of representation
requires thrie dierentiability of the urve and need the urvature of the urve to
be stritly positive. To overome these diulties, we take an alternate approah
5
for framing the urve, known as the Natural Frenet frame [Bishop, 1975℄. This
approah requires only twie dierentiability and is well dened even when the
seond derivative vanishes.
In 3D, to any point on the urve r(t), we attah an orthonormal moving frame
{x(t),y(t), z(t)}. The unit vetor x(t) is tangent to the urve and points toward
the heading of an individual. The unit vetors {y(t), z(t)} are hosen in the plane
normal to x(t). The evolution of these vetors are given by the frame equations,
ṙ(t) = ν(t)x(t)
ẋ(t) = ν(t)(u(t)y(t) + v(t)z(t))
(1.1)
ẏ(t) = −ν(t)u(t)x(t)
ṙ(t) = −ν(t)v(t)x(t),
where ν(t) is the speed (‖ṙ(t)‖) and (u(t), v(t)) are alled natural urvatures of
the trajetory [Justh and Krishnaprasad, 2005℄. In a planar setting, we have the
frame {x(t),y(t)} and the evolution equations are written as,
ṙ(t) = ν(t)x(t)
ẋ(t) = ν(t)u(t)y(t) (1.2)
ẏ(t) = −ν(t)u(t)x(t).
We an therefore treat ν and u variables as ontrol inputs to steer the individual
on the plane, ν as the veloity input and u as the urvature ontrol input.
6
Figure 1.1: Mobile robot based experimental platform (Pioneer 3 DX) with two-
wheel dierential and aster.
1.2 Experimental Setup
We provide a omprehensive desription of the laboratory set up in the Intelli-
gent Servosystems Lab, University of Maryland. All the laboratory experiments
presented in this thesis are done under this setup. Our experimental test-bed is
omprised of Pioneer 3 DX wheeled robots from Adept MobileRobots [Pioneer℄.
These ompat, dierential-drive mobile robots are equipped with reversible DC
motors, high-resolution motion enoders and 19m wheels, and the onboard om-
putation is done via a 32-bit Renesas SH2-7144 RISC miroproessor, inluding
the P3-SH miroontroller with ARCOS. The sensors on the robot inlude eight
forward-faing ultrasoni (sonar) sensors. ARIA [ROS-ARIA℄ provides an inter-
fae for ontrolling and reeiving data from the robot, and ommuniation with
the robot for sending ontrol ommands (forward veloity and turning rate) is
done via 802.11-b/g/n networking. The width of the robot is 380 mm and it has
a swing radius of 260 mm.
7
Algorithm implementation (i.e, feedbak law omputation) has been done in
C++ using ROS [ROS℄, along with ROS-ARIA [ROS-ARIA℄, as the interfaing
robotis middleware. The experiments have been arried out in a laboratory en-
vironment equipped with a sub-millimeter aurate Vion motion apture system
[Vion℄. We use a Dell workstation to run ROS, and this omputer is onneted
to the Vion server via a dediated Ethernet onnetion.
The Vion system aptures the motion of the robots and sends out the position
and heading data to the omputer running ROS. The ontrol law program listens
to this data, and transmits the individual veloities and turning rates. Both of
these operations are arried out at a frequeny of 25 Hz. The ontrol law program
omputes the ontrols aording to the strategy that is spei to the problem
onsidered. The omputed veloity and urvature ontrol variables ν(t), u(t) an
be translated to the turning rate ω(t) (in degrees/se) as:
( )
180
ω(t) = ν(t)u(t). (1.3)
π
8
Part I
Synthesis of Colletive Motion:
Bottom-up Approahes
9
Chapter 2
Feedbak Laws for Colletive Motion
2.1 Introdution
Colletive motion plays a ruial role in modern day robotis and engineering. It
is beoming ommonplae for a group of unmanned, remote ontrolled vehiles
to be deployed to aomplish goals ranging from searh and resue to surveil-
lane. For the swarm of robots to funtion in a harmonious manner, it is very
important to ontrol them arefully. Natural olletives are indeed an inspira-
tion in this endeavor. On the other hand, a thorough study of those olletives
remain inomplete without understanding agent level interation laws. In this
hapter therefore, we will build models for olletive motion from `bottom-up',
i.e. from individual level ontrol strategy to ok level synhrony through in-
teration among agents. This hapter presents a range of theoretial results as
well as laboratory demonstrations of ontrol laws that we propose or has been
proposed before. This hapter has two main ontributory setions, most of whih
10
are taken verbatim from their respetive publiations. We start with a problem
that studies a partiular dyadi interation under the setting of a pursuit strategy
alled onstant bearing pursuit [Halder et al., 2016℄. We gain interesting insight
from this problem that onnets to the Kepler two body problem. The obtained
result of this problem is then used to solve a problem arising from a pratial robot
maneuvering senario. The last problem is purely experimental [Halder and Dey,
2015℄ whih demonstrate another dyadi interation strategy potentially useful
for surveillane, and a oking strategy involving many agents.
2.2 Steering for Beaon Pursuit under Limited Sens-
ing
In this setion, we will try to understand simple dyadi pursuit strategies (i.e.
strategies based on pairwise interations), and exploit them as building bloks
for synthesis of omplex motion patterns for olletives. In [Galloway et al.,
2009, 2013℄, using symmetry priniples and nonlinear dynamis, a spei strat-
egy, known as onstant bearing yli pursuit, is shown to produe a rih variety
of behaviors for appropriate hoies of parameters (bearing angles). In [Justh
and Krishnaprasad, 2006℄ a biologially plausible feedbak ontrol law is inves-
tigated that exeutes motion amouage, a type of stealthy pursuit assoiated
with visually-guided ight in insets (e.g. hoveries and dragonies). Stealth
arises from nulling opti ow in the visual eld of the target of pursuit, thereby
11
inreasing the hane of suess in prey apture or territorial battle against a on-
spei. This type of dyadi interation is also exploitable in oordinated motion,
for instane see [Mishiati and Krishnaprasad, 2010℄.
The present work is similar to [Justh and Krishnaprasad, 2004℄ in motivation.
We onsider a problem of two agents moving in a plane with onstant (not ne-
essarily idential) speeds and, one of them is free i.e. it assumes any open loop
steering (urvature) ontrol, while the other pursues it. The free agent may be
onstrued as a beaon and the pursuer's task is to reah a safe viinity of the
beaon and irulate around it. In the interesting ase when the beaon is sta-
tionary, but the pursuer has a sensor with limited eld of view (FOV) to detet
the beaon, the irling law proposed in [Justh and Krishnaprasad, 2004℄ may be
foiled. One goal of this work is to devise a prinipled approah (ontrol algorithm)
for this problem that opes with sensor limitation. We do this via a two-step pro-
ess. We rst analyze a slightly dierent problem of traking a (slowly) moving
beaon assuming that: (a) the beaon trak is of onstant urvature (i.e. on a
straight line or on a irle); and (b) the sensor on-board the pursuer has no FOV
limitation. For the hoie of a onstant bearing pursuit feedbak ontrol law, one
obtains a rih dynamis. The phase portrait in turn suggests the seond step  a
feedbak law modiation that is appliable to the setting of stationary beaon,
and limited FOV. In this ase, one needs an additional ingredient  an estimator
using odometry to trak the beaon when it has fallen out of the FOV. The idea
here is to use the odometry-based estimate in the feedbak law as if it is exat
12
PSfrag replaements
y1
x1 x2
κ1
y2
κ2
r1 r12 = r1 − r2
r2
Figure 2.1: Illustration of salar shape variables (ρ, κ1, κ2) used to parametrize
the shape spae.
(a type of ertainty equivalene) and use diret observation of the beaon when
it is re-aptured in the FOV. Suh an estimator is used in an implementation of
beaon traking in a laboratory test-bed with a range amera (Kinet [Kinet℄) as
the sensor mounted on a mobile ground robot. A high preision (marker-based)
motion apture system (Vion [Vion℄) is used to determine ground truth and
analyze the performane. In addition, the unonstrained traking of the moving
beaon problem is revisited and it is shown that the resulting dynamis an be
identied with motion of a harged partile in an eletromagneti eld. Moreover,
at a partiular value of the beaon urvature, the ombined dynamis is exatly
same as the Kepler problem of two bodies.
2.2.1 Traking a Moving Beaon
Let us onsider two agents moving on a plane, eah abiding the self steering par-
tile equations of motion. We assume that both their speeds ν1, ν2 are onstants.
It is possible to represent the dynamis of the system of two agents by the help
13
of salar shape variables ρ, κ1, κ2 (Fig. 2.1) as
ρ̇ = −ν1 cosκ1 − ν2 cosκ2
− 1κ̇1 = ν1u1 + (ν1 sin κ1 + ν2 sin κ2) (2.1)
ρ
− 1κ̇2 = ν2u2 + (ν1 sin κ1 + ν2 sin κ2).
ρ
This is rather a straightforward alulation. We view agent 1 as a slowly moving
beaon to whih agent 2 pays attention. Let us make the following assumption.
A-1: The speed of agent 1 is less than the speed of agent 2, i.e. ν1 < ν2.
We pik the feedbak ontrol law for agent 2 as follows:
( ) ( )
− r21 1 r21u2 = µ̃ R(α)y2 · | | −
⊥
r ν | · ṙr | |r | 21 , (2.2)21 2 21 21
⊥
for some µ̃ > 0. Here we denote a = R(π/2)a, for any vetor a in the plane of
motion, R(·) is the planar rotation matrix. Note that this ontrol law is a standard
onstant bearing (CB) pursuit law [Galloway et al., 2013℄ with parameter α. The
feedbak ontrol law an be expressed in terms of the salar shape variables as
1
u2 = µ̃ sin(κ2 − α) + (ν1 sin κ1 + ν2 sin κ2)
ν2ρ
The losed loop dynamis of (2.1) then takes the form
ρ̇ = −ν1 cosκ1 − ν2 cos κ2
1
κ̇1 = −ν1u1 + (ν1 sin κ1 + ν2 sin κ2) (2.3)
ρ
κ̇2 = −µ̃ν2 sin(κ2 − α).
A fundamental result [Galloway et al., 2013℄ for the CB strategy tells us that
under the ation of the ontrol law (2.2), the manifold
Mα + 1 1CB = {(ρ, κ1, κ2) ∈ R × S × S : κ2 = α}
14
is an attrative invariant manifold for all initial onditions exept κ2(0) = α +
π. The invariane follows diretly from the losed loop dynamis (2.3). The
attrativeness an be proved by dening Λ(t) = − cos(κ2(t)− α). Thus Λ(0) =6 1
implies Λ(t) → −1 as t → ∞ or equivalently κ2 onverges to α.
2.2.2 Dynamis Restrited to the Invariant Manifold
At this stage, we are ready to make another assumption:
A-2: We onsider the urvature of the beaon to be onstant, i.e. u1 = u, for
some u ∈ R onstant.
Now we fous our analysis on the dynamis on the invariant manifold (alled
M + 1Shape = R × S ) whih may be expressed as
ρ̇ = −ν1 cos κ1 − ν2 cosα
(2.4)
1
κ̇1 = −ν1u+ (ν1 sin κ1 + ν2 sinα).
ρ
It is of interest to haraterize the solutions of the restrited dynamis (2.4) on
the invariant manifold. Note that given ν1, ν2, u and α, (2.4) might have at most
(ρ∗, κ∗
ν sinκ∗+sinα
two equilibrium points 1), with cosκ
∗ = − cosα ρ∗ 11 and = , whereν νu
ν = ν1we denote < 1. Existene of suh equilibrium points is guaranteed if ν ≥
ν2
| cosα| and (ν sin κ∗1 + sinα)u > 0. Linearizing (2.4) around suh an equilibrium
point gives the Jaobian matrix
 
0 ν sin κ∗
 1
ν  2 ,
 
−νu − cosα
ρ∗ ρ∗
15
with assoiated Eigenvalues
( )
√
− cosα± cos2 α− 4ν2ρ∗u sinκ∗
λ = ν 12 .
2ρ∗
Depending on α, following ases will arise
I. α ∈ (−π,−π/2) ∪ (π/2, π]
In this region, cosα < 0, whih makes the equilibrium points unstable. So
all trajetories tend to blow up in (ρ, κ1) plane.
II. α ∈ (−π/2, π/2)
Here, cosα > 0 ∗, then (loally) stable equilibrium exists if u sin κ1 > 0,
otherwise (ρ∗, κ∗1) is unstable whih leads to eventual ollision. We note
that [Davis, 1962℄ (pages 119-125) studies the same problem with α = 0.
For the α = 0 ase, the existene onditions of equilibrium read ν ≥ 1 and
u sinκ∗1 > 0, whih in turn means we will have a stable equilibrium only
when ν ≥ 1. The urrent problem an be viewed as a generalization of that
onsidered in [Davis, 1962℄.
III. α ∈ {π/2,−π/2}
In this ase, however, the dynamis (2.4) produes a rih behavior whih we
analyze next. We only provide the analysis for α = π/2 ase, the other ase
being analogous.
16
π
2.2.3 Speial Case: α =
2
We rewrite (2.4) in this partiular ase,
ρ̇ = −ν1 cosκ1
(2.5)
1
κ̇1 = −ν1u+ (ν2 + ν1 sin κ1).
ρ
All the trajetories of (2.5) are losed and we will prove this in the same way as
exploited in [Mishiati and Krishnaprasad, 2012℄. We rst introdue the following
denitions and a theorem due to Birkho.
Denition 2.1 (Involution). A dieomorphism F : M → M from a manifold
M onto itself is said to be an involution if F 6= idM , the identity dieomorphism
2
and F = idM , i.e. F (F (m)) = m, ∀m ∈ M .
Denition 2.2 (F-reversibility). A vetor eld X dened over a manifold M
is said to be F-reversible if there exists an involution F suh that F∗(X) = −X,
i.e. F maps orbits of X to orbits of X , reversing the time parametrization. Here
(F∗(X))(m) = (DF )F−1(m)X(F
−1(m)), ∀m ∈ M is the push-forward of F . We
all F the reverser of X .
Theorem 2.2.1 (G.D. Birkho, [Birkho, 1915℄). Let X be a F-reversible vetor
eld on M and ΣF the xed-point set of the reverser F . If an orbit of X through
a point of ΣF intersets ΣF at another point, then it is periodi.
See [Mishiati and Krishnaprasad, 2012℄ for a detailed proof of this theorem.
Based on these denitions and the theorem, we propose the following.
17
90   4 90   8 90   4
120 60 120 60 120 60
  3   6   3
150   2 30 150   4 30 150   2 30
  1   2   1
180 0 180 0 180 0
210 330 210 330 210 330
240 300 240 300 240 300
270 270 270
(a) u = −2 (b) u = 0 () u = +2
Figure 2.2: Phase portrait (polar plot) of the system dynamis restrited to the
invariant manifold (2.5) with dierent values of u, keeping ν1 = 0.5, ν2 = 1 xed
for all three ases.
Theorem 2.2.2. (i) The quantity,
1
E(ρ, κ1) = ρ(ν2 + ν1 sin κ1)− ν1uρ2 = E(ρ(0), κ1(0)) (2.6)
2
is onserved along any trajetory of (2.5).
(ii) Every solution of (2.5) is periodi.
Proof. (i) Denote, χ = ρ(ν2 + ν1 sin κ1), then
dχ
= ρ̇(ν2 + ν1 sin κ1) + (ρν1 cosκ1)κ̇1
dt
dρ
= ν1uρ
dt
=⇒ 1χ(t) = ν1uρ2(t) + c,
2
where c = χ(0)− ν1uρ2(0)/2 = onstant, whih, in turn implies
1
E(ρ, κ1) = ρ(ν2 + ν1 sin κ1)− ν1uρ2 = E(ρ(0), κ1(0)).
2
18
(ii) Step-1: Vetor eld dened by (2.5) is F-reversible with reverser F (ρ, κ1) =
(ρ, π − κ1).
2
Clearly, F is an involution sine F (ρ, κ1) = (ρ, κ1). Next,
(F∗(X))(ρ, κ1) = (DF )F−1(ρ,κ )X(F
−1(ρ, κ1))1
 
1 0
 
=  X(ρ, π − κ1)
 
0 −1
= −X(ρ, κ1).
Hene, X is F-reversible.
π
Step-2: Fixed point set of F is given by ΣF = {(ρ, κ1) : ρ > 0, κ1 = ± }. So2
π
every orbit of (2.5) rossing κ1 = ± line twie is periodi. Now, depending on2
the value of u, dierent ases will arise.
(a) u ≤ 0 : In this ase, we note that the assumption ν2 > ν1 is suient
to guarantee monotoniity of κ̇1, in partiular κ̇1 > 0 for all time. Hene, any
trajetory originating from any point on the κ1 = ±π/2 line (exluding the origin)
will travel ounter lokwise until it hits the line again when κ1 gets inremented
by an amount of π radian (see Fig. 2.2a, 2.2b). Note that the onserved energy, E
prohibits any trajetory that starts with positive energy to go to the origin (with
zero energy).
(b) u > 0 : Beause κ̇1 an assume any sign under this ase, we need a more
serious argument for this ase. To determine the nullline κ̇1 = 0, we ompute
19
PSfrag replaements
90
  4
120 60
A   3
150   2 30
  1
III
IV
I
II
180 0
O
B
NK
210 330
240 300
270
Figure 2.3: Nulllines of (2.5) with u = 1, ν1 = 0.5, ν2 = 1.
ρ = ν2+ν1 sinκ1 and ν2 > ν1 ensures existene of a valid ρ for eah value of κν u 1. The1
ρ̇ nullline is simply the κ1 = ±π/2 line. This is illustrated in Fig. 2.3. It is
( )
ν2+ν1 π
immediate that in this ase we have two equilibrium points, A = ,+ and
ν1u 2
( )
= ν2−ν1 ,−πB . The trajetories starting from either of those points are learly
ν1u 2
periodi.
Depending on the nulllines, the whole spae an be divided into four regions as
shown in Fig. 2.3 and those regions are haraterized as
Region I : ρ̇ < 0, κ̇1 < 0, Region III : ρ̇ > 0, κ̇1 > 0,
Region II : ρ̇ > 0, κ̇1 < 0, Region IV : ρ̇ < 0, κ̇1 > 0.
Now imagine trajetories starting on the line segment OA, exluding both
points. Sine κ̇1 < 0, they will move into region III, whih an produe two
outomes:
(b1) It hits the OB line (exluding both points).
20
(b2) It exits region III through NK and goes into region II. Now the onstany
of E gives the following observation, if a trajetory of (2.5) rosses NK at
(1) (2)
a ertain angle κ1 , then it must ross it again at a dierent angle κ1 and
these angles are symmetri about κ1 = π/2. This means trajetories must
not enter region III from region II. This leads all the trajetories in region
II to hit the κ1 = +π/2 line beyond point A.
Next, onsider trajetories starting at the boundary of region I and II with
κ1 = π/2. κ̇1 < 0 gives rise to lokwise motion into region I. Again, we need to
analyze two senarios:
(b3) The trajetories reah boundary between region I and II with κ1 = −π/2.
(b4) They enter region IV through NK. Similar argument as in ase (b2) an
be employed to prove they must reah boundary between region III and IV
with κ1 = π/2.
Now, trajetories starting on OB line (exluding both points) must move in
region IV and hene must hit OA line (exluding both points).
Finally, trajetories starting on the boundary between region I and II with
κ1 = −π/2 must go into region II and must eventually hit the boundary between
region I and II with κ1 = π/2 following a lokwise path (and without entering
region III).
This ompletes the proof. 
Remark 2.1. In a speial ase u = 0, (2.6) reads exatly as the polar equation of
1
an ellipse = c̃(1 + ν cos θ), where the origin is plaed at one fous of the ellipse,
ρ
21
the angle θ is measured from the origin with respet to the major axis of the ellipse
ν = ν1and is the eentriity of the ellipse (reall A-1). It is of interest to note
ν2
that the same ellipti equation omes from the analysis of the Kepler two-body
problem [Goldstein et al., 2001℄.
A result of the dynamis (2.5) regarding time period is immediate in the light
of Proposition 2.2.2.
Corollary 2.2.1. Every orbit of (2.4) has a period
∫ ρmax dρ
T = 2 √ , (2.7)
ρ 2 1min ν1 − ( ν1uρ+ E0 − ν 22 ρ 2)
where ρmin and ρmax are solutions of the pair of equations
1
ρ(ν ± ν )− ν uρ22 1 1 = E(ρ(0), κ1(0)) =: E0 (2.8)
2
2πν2E0
In partiular, for the ase u = 0, the time period beomes, T = 3 .
(ν2 2 22−ν1 )
Remark 2.2. Note that for an admissible value of E0, the pair of equations
(2.8) has only two solutions. For the speial ase, u = 0, we know the losed loop
trajetories are desribed by the ellipses ρ = E0/(ν2 + ν1 sin κ1) with semi-major
( )
a = 1axis, (ρmin + ρ
E0 1 1 E0ν2
2 max
) = +
2 ν +ν ν − =ν ν2−ν2 . Then, from Corollary2 1 2 1 2 1
( )
2 4π2 3
2.2.1, we nd T = a .
ν2E0
Remark 2.3. The ondition ν2 > ν1 is neessary for the existene of periodi
orbits for the ase u = 0 while it is merely a suient ondition for other values
of u.
22
2.2.4 The Limited Field of View (FOV) Problem
In this setion, we will desribe the problem of traking a stationary beaon by
a ontrolled agent equipped with a sensor (e.g. a depth amera like Kinet) with
limited eld of view. As opposed to the problem disussed in the previous setion,
this problem is inspired by an implementation point of view. The bakbone model
(2.1) of the system stays the same. Agent 2 is supposed to sense the position of
the beaon relative to its own position and use the sensed quantities to determine
the ontrol ation. Using the shape variables, it has aess to the pair (ρ, κ2) (refer
to Fig. 2.1), with the limitation that |κ2| ≤ κmax < π/2, whih we all the eld of
view onstraint. Although various feedbak ontrol laws have been proposed (for
e.g. [Justh and Krishnaprasad, 2004℄) to enirle a stationary beaon, permanent
loss of the target (beaon) from the eld of view annot be avoided by those laws.
More preisely, the limited FOV problem boils down to enirling the beaon while
being able to sense it (at least) periodially.
Putting beaon speed, ν1 = 0 in (2.1) and ignoring κ1 dynamis, the equivalent
shape spae equations an be redued to
ρ̇ = −ν2 cosκ2
(2.9)
− 1κ̇2 = ν2u2 + ν2 sin κ2.
ρ
Remark 2.4. From (2.9), it is guaranteed that under the eld of view onstraint,
the attempt of enirling the beaon would eventually make the beaon perma-
nently invisible (as long as a irular orbit around the beaon is onsidered).
Moreover, from ρ dynamis, meeting the onstraint |κ2| ≤ κmax < π/2 for all time
23
will lead to denite ollision.
One the beaon goes out of the FOV, the only hoie for the vehile would be
to eiently estimate the position of the beaon and apply the ontrol based on
those estimates. Aepting the fat mentioned in remark 2.4, one an only try to
design the ontrol u2 in suh a way that the ontrol law provides some promise to
bring the beaon bak in the eld of view after losing it. The following proposition
is meant to serve that purpose.
Proposition 2.2.1. The feedbak ontrol law given by
uFOV
µ
2 = u0 − , u0 ≤ 0, µ > ν2, (2.10)ρν2
guarantees the periodi return of the beaon to the eld of view under ideal esti-
mates.
Proof. With the feedbak ontrol (2.10), the losed loop system beomes,
ρ̇ = −ν2 cosκ2
(2.11)
κ̇2 = −
1
ν2u0 + (µ+ ν2 sin κ2).
ρ
Notiing that (2.11) is equivalent to (2.5), the laim follows diretly from Theorem
2.2.2. From the polar phase portrait (Fig. 2.2), we see that the ondition u0 ≤ 0
◦
is required for the angle variable κ2 to go through a full 360 rotation whih is
essential in order to bring the beaon bak in the FOV. 
As we will disover next, the ondition on the parameter u0 an be relaxed
to inlude positive values as well. From Theorem 2.2.2, the ondition µ > ν2 is
24
only suient for any value of u0. It beomes neessary for the partiular ase of
u0 = 0.
Although the feedbak ontrol law (2.10) produes all periodi orbits and it
inherently takes are of ollision avoidane problem (see Fig. 2.2), it laks in the
freedom of driving the vehile to a partiular desired orbit (f. the irular orbits
with desired radius as in [Justh and Krishnaprasad, 2004℄). Moreover, the periodi
orbits are not orbitally asymptotially stable whih makes them suseptible to
disturbanes. To overome these shortomings, we propose the following.
Proposition 2.2.2. Let Ed denote the admissible value of the desired energy, i.e.
there exists a periodi orbit with E(ρ, κ2) = Ed. Here E(ρ, κ2) is as in (2.6) with
the pair (ν1, ν2) interpreted as (ν2, µ) in present ontext. Then the ontrol law
u = uFOV + uAD
µ
2 2 2 = u0 − + kd(E(ρ, κ2)−Ed) cosκ2, (2.12)ρν2
with kd > 0 makes the orbit with energy Ed asymptotially stable with region of
attration given by MShape \ {(ρ, κ2) : cosκ (µ+ν2 sinκ2)2 = 0, ρ = , u0 > 0}, whereν2u0
M = R+ 1Shape × S .
Proof. Note that the trajetories of (2.9) with ontrol (2.12) will no longer be peri-
AD
odi beause of the inlusion of the extra u2 term. Sine E(ρ, κ2) is a ontinuous
2
funtion of both ρ and κ2, it sues to prove that the quantity (E(ρ, κ2) − Ed)
is monotonially dereasing. We obtain,
d
(E(ρ, κ2)−E 2d) = −2kdρν22(E(ρ, κ )− E )22 d cos2 κ2.dt
25
Clearly, with kd > 0, d(E(ρ, κ2)−Ed)2/dt ≤ 0, for all (ρ, κ2) suh that ρ > 0. The
2
largest invariant subset of {(ρ, κ2) : d(E(ρ, κ2)−Ed) /dt = 0} is indeed {(ρ, κ2) :
E(ρ, κ2) = Ed}, whih, in turn establishes the statement of the proposition. 
In the light of Theorem 2.2.2, we see that the restrition u0 ≤ 0 in Proposition
2.2.1 an be relaxed. In partiular, for u0 > 0, one only has to hoose Ed suh
that k2 ompletes full 360
◦
rotation (for e.g. one might pik Ed = E(ρ,−π/2),
µ−ν2
with ρ > , see Fig. 2.2()).
ν2u0
2.2.5 Implementation
In this laboratory implementation, we hose to use the newest Kinet model,
whih was reated for Mirosoft's Xbox One. The Kinet primarily funtions as
a motion-sensing input devie, enabling players to interat with video games in
exiting ways. To aomplish this, the devie is equipped with several sensors
inluding an RGB sensor, 3D Depth Sensor, as well as Multi-array Mirophones.
The Kinet's RGB sensor has a 70.6 degree horizontal eld of view, and a 60
degree vertial eld of view (see Fig. 2.4). The Kinet operates at a rate of 30
Hz, and has an eetive range between 0.5 meters, and 4.5 meters where auray
is reliable. Despite it's original use ase as a video game ontroller, the Kinet
has been studied reently as a sensor for many robotis appliations, inluding
autonomous vehiles [Oliver et al., 2012℄ and healthare [Nghiem et al., 2012℄.
In this experiment, the Kinet RGB amera ats as a primary sensor for de-
termining the distane and relative heading of the beaon (i.e. ρ, κ2), whih in
26
Figure 2.4: Robot (Pioneer 3 DX) with Kinet mounted, and the orange one
used as the beaon.
our experiment was an orange one. OpenCV [OpenCV℄ is used to perform a
simple blob detetion algorithm that alulates the entroid of the one in pixel
oordinates, and then uses the Kinet's oordinate mapping feature to transform
the result into physial, or amera spae. To take advantage of these API features,
we mount a laptop running the Windows operating system onto the robot, and
utilize a ustom TCP/IP server to stream the oordinates bak to the robot on-
trol station. The ontrol station is a Dell omputer running ROS [ROS℄, and the
algorithm implementation is done using the MATLAB ROS toolbox [MATLAB℄.
Finally, the Vion motion apture system [Vion℄ is used to trak the motion of
the robot and beaon in the lab oordinate spae to obtain ground truth results
of the implementation.
2.2.5.1 Estimation
In order to suessfully implement the proposed ontrol law, the robot has to be
able to eiently determine the beaon position relative to its own position during
27
the periods of time when the beaon is not in FOV. Equivalently, the estimation
problem then is to integrate the losed loop shape dynamis ((2.9) with ontrol
(2.12)) given some initial ondition, whih in this ase would be the last known
(ρ, κ2) value before the robot loses sight of the beaon. Sine there is a onserved
quantity assoiated with (2.11) (Proposition 2.2.2), a mid-point based update
rule performs better than the naive Euler rule [Austin et al., 1993℄. Denoting the
estimate of (ρ, κ2) by (ρ̂, κ̂2) and the disrete time step by ∆t, the update rule
may be impliitly expressed as follows:
(
ρ̂n+1 − ρ̂n κ̂n + κ̂n+1)
= −ν 2 22 cos
∆t 2
(2.13)
(
κ̂n+1 − κ̂n 2ν κ̂n + κ̂n+1)2 2 = − 2ν2un2 + sin 2 2 ,∆t ρ̂n+1 + ρ̂n 2
n n n 0 0
where u2 = u2(ρ̂ , κ̂2 ) as in (2.12) and (ρ̂ , κ̂2) is the last suessful measurement
of the beaon position. We then solve the nonlinear equations (2.13) numerially
(using MATLAB's fsolve) to produe the neessary estimates whenever the beaon
is not in the eld of view of the sensor. This proedure an be summarized in
Algorithm 1.
2.2.5.2 Experimental Results
To demonstrate our solution to the limited eld of view problem, we onstruted
an experiment for whih the robot sees the orange one and attempts to enir-
le it using the desribed ontrol mehanisms. The result is a trajetory that
periodially brings the one bak in its FOV so that the robot an fulll its net
28
Algorithm 1: Steering Law Computation for Limited FOV Problem
Data: ρ : measured distane to beaon, κ2 : measured angle to beaon;
ρ̂ : estimate of ρ, κ̂2 : estimate of κ2
Parameters: u0, µ, kd, Ed
begin
while not stopped do
if beaon visible then
Compute u2 = u2(ρ, κ2) using equation (2.12)
else
0 0 last last
Initialize: (ρ̂ , κ̂2) ←− (ρ , κ2 )
Calulate (ρ̂, κ̂2) from (2.13)
Determine u2 = u2 (ρ̂, κ̂2) using (2.12)
Expe rimental
 90   2000 Ideal
Beacon
2000 120 60Robot
  1500
) 1500
150   1000 30
m
m 1000   500
n
i
( 500 180 0
PSfrag replaements 0
−500 210 330
−1000
 240 300
−2000 −1000 0 1000 2000
x 270(in mm)  
(a) Robot Trajetory (b) Phase Portrait
Figure 2.5: Implementation results as reorded by Vion (with ν2 = 200 mm/se,
u −10 = 0 mm , µ = 1000 mm/se, kd = 5 × 10−9 −3 6mm se and Ed = 1.44 × 10
2
mm /se).
29
y
enirling of the beaon. The ground truth data was obtained using Vion. We
ran the experiment using a robot speed ν2 = 200 mm/se, and the parameters
were u0 = 0
−1
mm , µ = 5ν −9 −32 = 1000 mm/se, kd = 5 × 10 mm se and Ed
was taken to be energy orresponding to the orbit whih maintains the minimum
6 2
distane of 1200 mm from the beaon, in other words Ed = 1.44× 10 mm /se.
Implementation results are shown in Fig. 2.5. The ground truth polar plot
an be seen in omparison to the desired ellipse (sine u0 was taken to be 0) in
Fig. 2.5b. The mid-point rule estimation method results in a robust ontroller
that ahieves the desired trajetory although it is slightly loser to the beaon
than the theory predits. The error between these two orbits is observed (∼ 200
mm) to be within the size of the robot (∼ 400 mm).
2.2.6 Assoiated Lagrangian
Here we will re-visit the problem of traking a moving beaon as onsidered in
Setion 2.2.1. It is of spei interest to ask whether the system dynamis admits
some Lagrangian formulation. Without loss of generality, at this stage we take
ν2 = 1 and denote ν1 = ν (note that A-1 translates to ν < 1). Writing r = r1−r2,
we see that on the invariant manifold, the feedbak ontrol law (2.2) takes the
form
r · ṙ⊥
u2 = − | .r|2
30
Now it is a straightforward exerise to see that the baseline vetor satises the
following seond order ODE:
r̈ = ν2uy1 − u2y2
( )( )
2 r · ṙ⊥ r= ν uy1 + | . (2.14)r|2 |r|
Note that here we used the fat that on the invariant manifold the bearing angle
of the seond agent to the beaon, κ2 = π/2 and hene the y2 vetor is aligned
( )
r
with the baseline vetor r, or in other words y2 = | | .r
Proposition 2.2.3. On every level set of E, (E(ρ, κ1) = E0, as in (2.6)) the
two dimensional system (2.14) is atually the Euler-Lagrange equation of the La-
grangian funtion (of the type kineti energy−potential energy)
L(r, ṙ) = K(r, ṙ)− V (r, ṙ)
( )
1 | |2 − −E0 − 1= ṙ | | νu|r| − A(r) · ṙ , (2.15)2 r 2
1 ⊥
where A is dened as A(r) := − νur .
2
Proof. Note that the quantity, E0 = ρ(1 + ν sin κ) − 1νuρ2 is onserved and2
r · ṙ⊥ = −ρ(1 + ν sin κ). From here, we an rewrite (2.14) as
( )( )
r̈ = νu(νy1 −
E0 1 r
y2) + −| | + νur 2 2 |r|
( )
E0 1
= (ṙ×B)−∇ −| | − νu|r|r 2
= (ṙ×B) + E(r), (2.16)
where we introdue B := −νuẑ, a stati magneti eld in the diretion per-
pendiular to the plane of motion of the agents, ẑ being the unit vetor in
31
that diretion. Also, we all E(r) := −∇Φ(r), the eletri eld and Φ(r) :=
( )
− E0 + 1| | νu|r| , the orresponding eletrostati potential.r 2
Now, the equation (2.16) resembles with that of the equation of motion of a
harged partile in an eletromagneti eld and one an nd obvious similarity
of the right hand side with that of Lorentz fore law. Then a standard result
[Goldstein et al., 2001℄ in the theory of eletromagnetism gives the Lagrangian
formulation of (2.16). Sine ∇ · B = 0, B an be written as url of a magneti
vetor potential A(t, r), i.e. B = ∇×A. Also, the eletri eld E an be writ-
ten as E = −∇Φ(r) − ∂A . Mathing this with (2.16), we an see that A is a
∂t
vetor valued funtion of r only. It is a straightforward exerise to show that
A = −1(r × B) = −1νur⊥ satises B = ∇ × A. Then the Lagrangian whih
2 2
generates (2.16) is given by
L(r, ṙ) = K(r, ṙ)− V (r, ṙ)
1
= |ṙ|2 − (Φ(r)−A(r) · ṙ)
2
1
= |ṙ|2 E0 1 1+ | | + νu|r| − νur
⊥ · ṙ
2 r 2 2

Remark 2.5. In essene, Proposition 2.2.3 reveals that with open loop onstant
urvature ontrol of one agent and with feedbak ontrol (2.2), namely onstant
bearing pursuit law with parameter α = π/2, of the other, (on every level set
of an invariant manifold) the oupled system behaves exatly the same way as a
harged partile in a stati eletromagneti eld. Moreover, the magneti eld, B
is dependent on the onstant urvature (u) of the rst agent. Thus u = 0 implies
32
the absene of the magneti eld and the agents are subjet to the eletrostati
eld only. The eld E has familiar inverse square form and hene the results
for the speial ase u = 0 agree with Kepler's laws whih are obtained from a
two-body problem subjet to Newtonian gravitational fore.
2.3 Biomimeti Algorithms for Coordinated Mo-
tion
In this setion, we will report implementation of two feedbak ontrol strategies
on our laboratory test-bed. The rst of these two strategies is alled mutual mo-
tion amouage (MMC) [Mishiati and Krishnaprasad, 2012℄. Existing literature
on dragonies [Corbet, 1999℄ provides qualitative analysis of territorial battles,
wherein the trajetories display spiraling motion onsistent with the theoretial
preditions [Mishiati and Krishnaprasad, 2011℄. This partiular bio-inspired on-
trol algorithm inherits an appealing overage property through the mehanism of
spae lling urves, and our implementations are able to reprodue overage pat-
terns similar to the predited ones.
Although there has been a long history of ontrol algorithms for oking, al-
most every model of olletive motion predits diusive transport of information.
But, ontrary to the existing models, reent ndings [Attanasi et al., 2014℄ from
starling oks suggest that diretional information within a ok propagates with
an almost onstant speed, and this linear growth of information an be explained
33
by models with wave-like aspets. In [Dey, 2015; Halder and Dey, 2015℄ a on-
trol strategy alled topologial veloity alignment (TVA) was introdued, whih
onforms to this riterion and an explain how information about loal neighbors
an inuene the agents in a ok to align their headings in a single ommon
diretion. Hene it seems reasonable to use TVA strategy for olletive motion
synthesis. Furthermore, our implementation results in real robots have shown
that redution in neighborhood size and external perturbation (similar to preda-
tor attak) an split a ok into smaller subgroups.
2.3.1 Mutual Motion Camouage (MMC)
Here we onsider the mutual interation between two agents eah applying the
same pursuit law, while pereiving the other one as a target. As the dynamis of
MMC in a planar setting has been studied earlier [Mishiati and Krishnaprasad,
2012℄, we just reiterate some key results in order to have a omprehensive frame-
work. Allowing dierent speeds for the agents, we begin with the following sym-
metry:
u1ν1 = u2ν2 = u. (2.17)
Then the dynamis of the relative motion vetors, namely r = r1 − r2, g = ṙ =
ν x ⊥ ⊥1 1 − ν2x2 and h = g = ṙ , an be expressed as
ṙ = g
ġ = uh (2.18)
ḣ = −ug.
34
Now we introdue three salar shape variables dened as ρ = |r|, γ = (r · g)/|r|
and λ = (r ·h)/|r|. Then, aording to [Justh and Krishnaprasad, 2006℄, we have
( ) ( )
− ru = µ ⊥ r|r| · ṙ = −µ | | · h = −µλ, (2.19)r
where, µ > 0 denotes the feedbak gain. As shown earlier, the dynamis of relative
motion (2.18) an be redued to yield a seond order dynamis given by
ρ̇ = γ
(2.20)
( )
γ̇ = (1/ρ− µ) δ2 − γ2 ,
where, δ = |g| = |h| is onserved along any trajetory of (2.18). As detailed in the
original work [Mishiati and Krishnaprasad, 2012℄, individual trajetories an be
reonstruted from the solutions of (2.20). Moreover, the solutions of the redued
dynamis (2.20) onstitute level sets for another onserved quantity, dened as
E(ρ, γ) = ρ2(δ2 − γ2)e−2µρ = E(ρ0, γ0). (2.21)
However, the absene of damping in the redued dynamis (2.20) has poten-
tial to deteriorate the performane of the original MMC law (2.19). A modied
feedbak law, with an added dissipative term to neutralize any deviation from the
predited trajetories, an be expressed as
( )
utot = u+ udis = −µλ+ kdλγ E(ρ, γ)−Ed , (2.22)
where Ed is set as the initial value of the onserved quantity E(ρ, γ), i.e. Ed =
E(ρ0, γ0). Previous work [Mishiati and Krishnaprasad, 2010℄ has shown that this
modied ontrol law (2.22) with kd > 0 makes the periodi orbit (with energy Ed)
35
orbitally asymptotially stable, and the orresponding domain of attration is
haraterized by {(ρ, γ) : ρ > 0,−δ < γ < δ, (ρ, γ) 6= (1/µ, 0)}.
2.3.2 Topologial Veloity Alignment (TVA)
Here we formalize the strategy of topologial veloity alignment (TVA) [Halder
and Dey, 2015; Dey, 2015℄, and assume that eah member in a group of n-agents
uses this strategy to move together while keeping its heading parallel to the neigh-
borhood enter of mass veloity. Letting Ni denote the neighborhood of the i-th
agent, the enter of mass (COM) veloity of this neighborhood is given by
1 ∑
v = ν x , (2.23)
COM |N | j ji j∈Ni
where |Ni| represents the number of neighbors inuening the i-th agent. Next,
by assuming that v does not vanish to zero, we dene the diretion of the
COM
enter of mass motion as
v
x = COMN | | . (2.24)i v
COM
It should be noted that xN is not well-dened over a thin set in the state spae.i
As the hane of getting into the thin set is very small, we an overlook this
situation for all pratial purposes. Now we introdue a ontrast funtion
1
Θi = (xN − xi) · (xN − xi) = 1− xi · xN , (2.25)
2 i i i
as a quantitative measure for the misalignment between the heading of an agent
and the diretion of motion of its neighborhood enter of mass. Clearly, this
ontrast funtion (Θi) assumes its minimum value (= 0) whenever the i-th agent's
36
veloity is aligned with its neighborhood enter of mass veloity, and inreases
monotonially with inrease in the misalignment between them. Thus, Θi an be
interpreted as a measure of departure from our goal of ahieving alignment.
Next, by assuming a non-zero veloity for the neighborhood enter of mass
(v =6 0), the TVA law is given by [Dey, 2015℄
COM
( · )xN yi iui = µ , (2.26)
νi
where µ > 0 denotes a positive gain, and yi arries its usual meaning. Alterna-
tively, lateral aeleration for this hoie of ontrol laws (2.26) an be expressed
as
( [ ] )
lat
ai = µνi xN − xN · xi xi , (2.27)i i
and this provides a physial intuition behind (2.26) as the lateral aeleration
is proportional to the projetion of the normalized veloity of its neighborhood
enter of mass onto the transverse of its own diretion of motion.
Remark 2.6. Earlier works [Justh and Krishnaprasad, 2003, 2004℄ have onsid-
ered a very similar form of this ontrol law with three omponents for attra-
tion (while the agents are far away), repulsion (to avoid ollision) and veloity
alignment. However, the TVA ontrol law onsiders only veloity alignment, but
extends the sope from a planar setting to a three dimensional environment. More-
over it relaxes the assumption on uniform speed of the olletive by allowing the
agents non-uniform and time-varying speed proles. This relaxation plays an im-
portant role in the ontext of applying this ontrol law to a group of heterogeneous
agents.
37
It was shown in [Halder and Dey, 2015; Dey, 2015℄, for a two agent system it
is possible to show that the TVA strategy (2.26) aligns the veloity of the agents,
only if their veloity vetors were not exatly opposite to eah other initially. [Dey,
2015℄ also provides further analytial results for a general n agent system with a
yli interation senario. As analysis of an n-agent system with neighborhood
dened as the set of K-nearest neighbors poses hard hallenges, we propose an
algorithmi way (Algorithm 2) to implement TVA in a real system. We bring in
an additional neighbor into onsideration whenever v beomes zero. Clearly,
COM
this provides a way to avoid ill-posedness assoiated with v being zero beause
COM
non-zero speeds of individual agents ensure that onsidering an extra neighbor will
make an otherwise zero v non-zero.
COM
2.3.3 Implementation Results
We present the implementation results of the two ontrol laws in our roboti test-
bed. In this setion, we are presenting results for whih the speeds of all the
individual agents are same, i.e. νi = νj , ∀i, j. Though it should be noted that
both ontrol laws an be implemented with dierent speeds.
2.3.3.1 Implementation of MMC
Here we will show some implementation results for MMC, and demonstrate the
eetiveness of using a dissipative ontrol term (2.22). Our analysis also inludes
a omparison between the observed trajetories and trajetories obtained from
38
Algorithm 2: Topologial Veloity Alignment
Data: Initial Time - tinitial; Final Time - tfinal; Sampling Interval - ∆;
n
Number of Agents - n; Initial Position and Orientation - {gi}i=1;
Neighborhood Size - K
begin
Initialize: tcurrent ←− tinitial ;
for i = 1 to n do
Initialize: State - Xi ←− gi ;
while tcurrent ≤ tfinal do
for i = 1 to n do
Dene: Ni - the set of K-nearest neighbors ;
Compute: Neighborhood enter of mass veloity - v ;
COM
if v = 0 then
COM
Dene: Ni - the set of K + 1-nearest neighbors ;
Compute: Neighborhood enter of mass veloity - v ;
COM
Compute: Steering ontrol - ui;
n
Implement: Steering Control - {ui}i=1 ;
n
Update: State - {Xi}i=1 ;
Update: Time - tcurrent ←− tcurrent +∆ ;
theoretial preditions, obtained via integrating the redued system dynamis
(2.20). Considering the prersene of a onserved quantity (2.6) in the system,
we used the method desribed in [Austin et al., 1993℄ for integration instead of
general ODE solver, whih otherwise would not be able to keep the quantity
E(ρ, γ) k konstant to its initial value. Then, from the updated values of ρ and γ ,
k k k
we reonstrut the trajetories along with their frame vetors, i.e. ri , xi , yi ,
with i = 1, 2 kand k denoting the time indies. At eah time instane t , the error
39
600 600 300
400 400
200
200 )200

e
0 0 s 100
/
) )
PSfrag replaements -200 PSfrag replaements PSfrag replaements-200 m
m m
0
m -400 -400 mm
n
n n
i
i -600 i -600 -100
(
( (
-800 -800
-1000 -200-1000
-1200 -1200 -300
-500 0 500 1000 -500 0 500 1000 500 1000 1500 2000
x (in mm) x (in mm) ρ (in mm)
(a) Trajetory (Experiment) (b) Trajetory (Simulation) () Phase plot
Figure 2.6: Performane of the system signiantly improved upon addition of
−1
the dissipative ontrol term (with µ = 0.001 mm , ν1 = ν2 = 200 mm/se and
kd = 1× 10−15 −6 3mm se )
(ek k ki ) is omputed as: ei = |ri,expt − rki,ideal|.
The plots of a sample run using the modied MMC feedbak law (2.22) are
shown in Fig 2.6. This modied ontrol law has been applied with kd = 1 ×
10−15 −6 3 −1mm se , and the parameters µ, ν1 and ν2 are seleted as 0.001 mm ,
200 mm/se and 200 mm/se, respetively. The resulting performane is quite
satisfatory as shown in Fig 2.6a (refer [YouTube℄ for implementation video). We
have also observed that the error is bounded (∼ 250 mm) within the size of the
robots (∼ 400 mm).
2.3.3.2 Implementation of TVA
We implemented the TVA ontrol law (2.26) in a 2 dimensional setting (i.e. vi(t) is
ignored). As the implementation is in disrete time, we followed Algorithm 2 in our
implementation in order to avoid the singular ase of |v | = 0. To demonstrate
COM
the performane of TVA ontrol law, we designed three dierent experiments (refer
40
y
y
γ
1000
2000 1000
500
0
1000 0
−500 −1000
0
−1000
−2000
−1000 −1500
−3000
−2000
−2000 −4000
−2500
−3000 −2000 −1000 0 1000 2000 −2000 −1000 0 1000 2000 −5000 −4000 −3000 −2000 −1000 0 1000 2000
x (in mm) x (in mm) x (in mm)
(a) Trajetories: 8 agents, ok- (b) Trajetories: 8 agents, split- () Trajetories: 6 agents, pertur-
ing ting bation
4.5 1.4 4 5
4
1.2
3.5 4
1
3
2.5 0.8
3
2
2 0.6 2
1.5
0.4
1 1
PSfrag replaements PSfrag replaements 0.2 PSfrag replaements
0.5
0 0 0 0
0 10 20 30 40 50 0 5 10 15 20 0 20 40 60 80
Time (in sec) Time (in sec) Time (in sec)
(d) Contrast Funtion: 8 agents, (e) Contrast Funtion: 8 agents, (f) Contrast Funtion: 8 agents,
oking splitting perturbation
Figure 2.7: Robot trajetories and ontrast funtions of TVA for (i) Experiment
1 [Fig (a),(d)℄ with 8 agents, demonstrates oking behavior (K = 3); (ii) Ex-
periment 2 [Fig (b),(e)℄ with 8 agents, desribes the splitting behavior due to low
neighborhood size (K = 1), and (iii) Experiment 3 [Fig (),(f)℄ with 6 agents,
shows that perturbation an ause a swarm to split, the trajetory of the per-
turbing agent is not shown. (µ = 1 Hz and νi = 60 mm/se is kept xed for all
experiments.)
[YouTube℄ for implementation videos). In these experiments, the sonar sensors
on the robots were ativated to sense any obstale in the diretion of motion of
the robots and if any robot an sense suh an obstale, it will simply apply a
sat
maximum turning rate (ω ) to avoid ollision. The sonars are programmed to
41
y (in mm)
Θ(t)
y (in mm)
Θ(t)
y (in mm)
Θ(t)
Number of nearest neighbors
detet an obstale only in lose proximity (∼ 300 mm) of the robots. In all our
sat
experiments ω is taken to be 50 rad/se, forward speeds of all of the robots
are kept onstant at 60 mm/se and the value of the parameter µ is hosen to
be 1 Hz. A system with eight agents is onsidered and we apply same TVA law
to all of them. The neighborhood size is taken to be three (i.e. K = 3). The
robots are initially plaed in arbitrary positions and diretions. The footprints of
∑
the robots and the orresponding ontrast funtion, Θ(t) = i Θi(t) is plotted
against time in Fig 2.7a, 2.7d. The initial and nal diretions of the robots
are shown using arrows and the nal positions of the robots are denoted using
dots. It an be seen from Fig 2.7d that the ontrast funtion deays to zero very
quikly whih indiates perfet veloity alignment within the swarm. Next we
dereased the neighborhood size and made it one (K = 1), so that every robot
only `ommuniates' with its losest neighbor. We hose the initial positions in
suh a way that they may form sub-lusters instead of moving as a single swarm.
This behavior is alled `splitting ' in a swarm. From Fig 2.7b, we an learly see
that the swarm of eight robots gradually split from eah other and form three
dierent lusters. It is to be noted that even if all the agents are not going in
the same diretion, the ontrast funtion still onverges to zero (Fig 2.7e). This
happens beause eah of the robots are aligned with their nearest neighbors and
hene eah of the individual ontrast funtions (Θi(t)) are zero. This experiment
may explain the splitting phenomenon observable in nature. Lastly, we ombined
the above two experiments, and onduted an experiment using six robots in a
42
swarm and another robot as a predator. A separate omputer was used for manual
ontrol of the `predator' robot.
At the beginning, neighborhood size is kept at K = 3, suh that the `om-
muniation' graph among the robots stays onneted and they move as an entire
swarm in a ommon diretion. When the swarm omes lose to the predator, the
neighborhood size is dereased to one. As we are not using any onboard visual
sensing and the sonar sensing is done only in very lose region (∼ 300 mm), the
hange in neighborhood size is made manually. From Fig 2.7f, we an see that
the hange in neighborhood size takes plae at around 20 seonds and we an
also see a tiny jump in the ontrast funtion at that time. The predator then
slowly approahes to one of the agents in the swarm, whih abiding to its ollision
avoidane rule, turns to avoid the predator. In Fig 2.7, the trajetories of the
agents are drawn in dashed lines before the ourrene of this event and in solid
lines afterwards. The trajetory of the predator robot in not shown in the g-
ure. After reating the initial perturbation, the predator is slowly moved through
the swarm ausing some subsequent disturbanes. These perturbations reate a
notieable impat in the swarm. As the attaked agent turns, its neighbor also
tries to align itself with that agent and so does its neighbor. This goes on until
the ommuniation graph beomes disonneted and a split in the swarm is then
observed [YouTube℄ just like in Experiment 2. As we an see in Fig 2.7, the
swarm is divided in two groups after the attak of the predator. The jumps in
plot of the ontrast funtion in Fig 2.7f symbolize the perturbations aused by
43
the external agent. The ontrast funtion eventually onverges to zero after the
members are aligned with their neighbors within eah subgroup.
44
Chapter 3
Optimal Steering of Agents on a Plane
3.1 Introdution
The kinemati uniyle model is often used in path-planning for ground vehiles,
sine the onguration of a ground vehile an often be represented by a point
in a plane that is onstrained to move in the diretion of the urrent heading
[Bellaïhe et al., 1998; LaValle, 2006℄. The state of this system an be represented
as an element of the speial Eulidean group SE(2), where the ontrol inputs are
a urvature input whih ontrols the rate of hange of the heading angle, and a
veloity input whih ontrols the rate of hange of the uniyle position in the
diretion of the heading angle.
Given the urrent onguration of the uniyle and a desired future ong-
uration, an admissible path for moving the uniyle from an initial to a nal
onguration an be determined via the minimization of some ost funtional.
There an be many variations to this problem depending on the hosen ost fun-
45
tional. A muh elebrated problem is the problem of Euler's elastia [Euler, 1744℄
where the minimum urvature path joining two given ongurations on the plane
is onsidered. Partiularly elegant and other well-known variations inlude the
minimum-time solutions of [Dubins, 1957℄ and [Reeds and Shepp, 1990℄. Optimal
paths of the minimum-time problem onsist only of straight-line and irular-ar
segments whih, when pathed together, reate disontinuities in the path ur-
vature and ause potential diulty in implementation sine abrupt hanges in
urvature are hard to trak. Proposed modiations that alleviate this problem
enfore that the urvature stay ontinuous, e.g., [Fraihard and Sheuer, 2004℄,
yet it is also possible to penalize the total urvature along the path in the expe-
tation that the optimal urvature will be ontinuous. [Halder and Kalabi, 2017℄
takes the latter approah, onsidering the minimization of the urvature along a
path onneting initial and nal uniyle ongurations with free nal time.
In this hapter, we will present a problem that penalizes both the urvature
and speed ontrols in maneuvering a uniyle from initial to desired nal ong-
uration. This helps both the urvature and speed ontrols to be smooth along an
optimal trajetory. These optimal trajetories losely resemble to those obtained
in [Halder and Kalabi, 2017℄. Our solution is obtained using geometri optimal
ontrol, where the neessary onditions for optimality are obtained via the Pon-
tryagin Maximum Priniple (PMP), and Lie-Poisson redution [Krishnaprasad,
1993; Ohsawa, 2013℄. Using geometri optimal ontrol on SE(2) to nd solutions
to path-planning problems has also been onsidered by [Sussmann and Tang, 1991;
46
Krishnaprasad, 1993; Agrahev and Sahkov, 2004; Dey and Krishnaprasad, 2014;
Justh and Krishnaprasad, 2015b℄.
Studying a uniyle on the plane is important sine often a olletive of N
agents exhibits single agent behavior under synhronization. In the later part of
this hapter, we will present a framework for analyzing suh a olletive. This
framework is based on [Justh and Krishnaprasad, 2015b℄, where every agent was
assumed to have onstant speed. Here we will, however, onsider a ost funtion
that penalizes both speed and urvature ontrols. In addition to the individual
ontrol osts, there is one ost that is attributed to the `mismath in ontrol' of
an agent with its `neighbors'. The neighbors of any agent is ditated by a xed
graph of interation. The strength of suh interation is aptured by a oupling
parameter. It is shown that in extreme ases (no oupling and high oupling)
the optimal olletive ontrols are diretly assoiated with optimal ontrols for a
single agent problem. This framework is what we use in the later hapters of this
thesis where we onsider optimal ontrol problems for a ontinuum of agents.
3.2 Optimal Steering of a Uniyle
In this setion, we onsider minimizing the urvature and speed ontrol osts of
a path in SE(2) onneting an initial uniyle onguration g0 with its desired
nal onguration gT at time T . We formulate this optimization as a geometri
optimal ontrol problem and derive the neessary onditions using PMP and Lie-
Poisson redution. From the neessary onditions, we show that there are two
47
onstants of motion: the Hamiltonian and the Casimir. We show that there are
three possible families of solutions depending on the values of Casimir and the
Hamiltonian. In the rst ase, the motion onsists of segments of a U-turn; in
the seond ase, the motion onsists of segments of parallel parking trajetories;
in the third ase, the motion onsists of straight lines or asymptoti approahes
thereto.
Consider the uniyle kinemati equations,
ẋ(t) = v(t) cos θ(t), (3.1a)
ẏ(t) = v(t) sin θ(t), (3.1b)
θ̇(t) = u(t), (3.1)
where (x(t), y(t)) ∈ 2R is the position of the uniyle on the Cartesian plane,
θ(t) ∈ S1 is the heading of the uniyle, v(t) is the uniyle speed ontrol, and
u(t) is the steering ontrol, equal to the rate of hange of the heading θ(t).
The onguration of the uniyle an be represented as an element of the
matrix Lie group SE(2). Let g(t) ∈ SE(2) where,
 
cos θ(t) − sin θ(t) x(t)
 
 
 
g(t) = sin θ(t) cos θ(t) y(t) . (3.2)
 
 
 
0 0 1
Then the equations (3.1) an be written in left-invariant form,
ġ(t) = g(t)ξ(u(t), v(t)), (3.3)
where,
ξ(u(t), v(t)) = u(t)X1 + v(t)X2, (3.4)
48
and,
   
0 −1 0 0 0 1
   
   
   
X 1 = 1 0 0 , X2 = 0 0 0 . (3.5)
   
   
   
0 0 0 0 0 0
The matries X1 and X2 are elements of the Lie algebra se(2). Together with
X3 = [X1, X2] = X1X2 −X2X1, X1 and X2 form a basis for se(2).
Without loss of generality, we an assume that g0 = I3, sine g(t) an always
−1
be redened aording to g(t) := g0 g(t). Given a nal time T > 0 and a nal
state gT , ost funtion to be minimized is,
∫
1 T ( )
u(t)2 + v(t)2 dt. (3.6)
2 0
3.2.1 Optimal Control Solution
In order to solve the problem we form the pre-Hamiltonian,
H = 〈 1p, gξ(u, v)〉 − (u2 + v2), (3.7)
2
where p(t) ∈ SE(2)∗ is the adjoint variable. To simplify the Hamiltonian, we
perform Lie-Poisson redution, introduing the variable µ(t) ∈ se(2)∗ satisfying
the translation to identity,
〈µ, ξ(u, v)〉 = 〈p, gξ(u, v)〉.
As an element of the dual spae, µ(t) an be represented as µ(t) = µ X♭1 1+µ
♭
2X2+
µ X♭ {X♭ ♭3 3, where 1, X2, X♭3} are the basis vetors dual to {X1, X2, X3}.
49
The pre-Hamiltonian therefore beomes,
1
H(µ, u, v) = 〈µ, ξ(u, v)〉 − (u2 + v2),
2
= 〈µ1X♭1 + µ2X♭2 + µ3X♭3, uX1 + vX2〉
− 1(u2 + v2),
2
1
= µ1u+ µ2v − (u2 + v2).
2
∗ ∗
Aording to the PMP, the optimal ontrol (u , v ) satises,
H(µ∗, u∗, v∗) = max H(µ∗, u, v). (3.8)
(u,v)∈R2
Therefore the optimal ontrols are given by,
u∗1 = µ1, (3.9)
u∗2 = µ2, (3.10)
The redued Hamiltonian is therefore,
1 ( )
H = H(µ, u∗1, u
∗
2) = µ
2 2
2 1
+ µ2 , (3.11)
and the dynamis of the µi variables are given by [Krishnaprasad, 1993℄,
   
µ̇1 0 −µ   3 µ2
   
    ∂h
µ̇  =  µ 0 0  , (3.12)
 2  3  ∂(µ1, µ2, µ3)
   
   
µ̇3 −µ2 0 0
µ̇1 = −µ2µ3, (3.13a)
µ̇2 = µ1µ3, (3.13b)
µ̇3 = −µ1µ2. (3.13)
50
(a) c = 0.5 (b) c = 1.5 () c = 1
Figure 3.1: M (blak) plotted as an intersetion of C (red) and H (blue) for h = 1
and three values of the Casimir c
Let,
c = µ22 + µ
2
3. (3.14)
This variable is alled the Casimir and it is a onstant of motion, implying that
√
the variables µ2 and µ3 evolve on irle of radius c. Along with the Casimir,
the Hamiltonian (3.11) is also a onserved quantity of motion. For onveniene of
subsequent alulations, we will work with the following saled Hamiltonian,
h = 2H = µ21 + µ
2
2. (3.15)
3.2.2 Charaterizing the Types of Motion
Aording to (3.14) and (3.15), the dynamis (3.13) evolve on the manifold M
where,
M = C ∩ H,
C = {(µ1, µ2, µ3) ∈ 3R : c = µ22 + µ23}, (3.16)
H = {(µ1, µ2, µ3) ∈ 3R : h = µ2 21 + µ2}.
51
The manifold M is one-dimensional and equal to the intersetion of the ylinders
C and H. The shape of M is determined by the Casimir c, while the shape of H
is determined by h. The motion of µ evolves on M. Due to ontinuity, it must
evolve on a onneted omponent of M, so it is important to onsider the types
of possible intersetions, of whih there are three orresponding to three dierent
ases: in Case 1, C is stritly smaller than H; in Case 2, C is stritly larger than
H; in Case 3, C is equal in size to H. To perform a ase-by-ase ategorization of
M √, we note that, aording to (3.14), the variable µ2 is restrited to ± c and,
√ √ √
aording to (3.15), µ2 is restrited to ± h. In Case 1, c < h, so |µ2| ≤ c < h.
Therefore the motion evolves on a onneted omponent of M where µ1 does not
2
hange sign sine µ1 = h − µ22 ≥ h − c > 0. Similarly, in Case 2, c > h, so
√
|µ2| ≤
√
h < c. The motion evolves on a onneted omponent of M where µ3
2 2
does not hange sign sine µ3 = c − µ2 ≥ c − h > 0. In Case 2, c = h, so the
two onstraints agree at the extremes. Instead of having M as a one onneted
omponent, it atually has four disonneted omponents. These omponents
√
meet eah other asymptotially at the extremes µ1 = µ3 = 0, µ2 = ± h. See
Fig. 3.1 for a visualization of the three ases. In the following, we study the three
types of motion in further detail.
Case 1: 0 < c < h
The equations (3.13) admit the following expliit solutions by means of Jaobian
52
6 6 6
5 4 4
4 2 2
3 0 0
2 -2 -2
1 -4 -4
0 -6 -6
-0.5 0 0.5 1 1.5 2 -2 0 2 4 -1 0 1 2
(a) c = 0.5 (solid blak) and 0.75 (b) c = 1.5 (solid blak) and 2 () c = 1
(dot-dashed red) (dot-dashed red)
Figure 3.2: Extremal trajetories for dierent ases (solid blak)  inreasing c
produes the dot-dashed red urve; inreasing η produes the solid blue urve;
hanging s1 or s2 produes the dashed green urve. h = 1 for all three ases.
ellipti funtions [Davis, 1962; Byrd and Friedman, 1971℄,
√ (√ )
µ1 = s1 h dn h(t+ η), k , (3.17a)
√ (√ )
µ2 = s1 c sn h(t+ η), k , (3.17b)
√ (√ )
µ3 = c n h(t+ η), k , (3.17)
2 c
where the modulus k of the ellipti funtions is given by k = . The parameters
h
s1 ∈ {1,−1} and η ∈ R do not depend on c and h, but on the endpoint onstraints.
4K√(k)
Note that these solutions are periodi with the period given by T1 := , whereh
K(k) denotes the omplete ellipti integral of the rst kind.
Proposition 3.2.1. Let 0 < c < h. Assume η = 0, then any extremal trajetory
53
is given by,
s ( (√ ))1
x(t) = 1− dn ht, k (3.18a)
k
1 (√ ( (√ ) ))
y(t) = ht− E am ht, k , k (3.18b)
k
(√ ) ( (√ )) ( (√ ))
θ(t) = s1 · am ht, k = s1 cos−1 n ht, k = s1 sin−1 sn ht, k ,
(3.18)
where E(·, ·) denotes the inomplete ellipti integral of the seond kind.
Proof. We have from (3.9) and (3.17a),
√ (√ )
u1 = µ1 = s1 h dn ht, k .
Integrating the equation θ̇ = u1 gives
(√ ) ( (√ )) ( (√ ))
θ(t) = s1 · am ht, k = s cos−1 ht, k = s sin−11 n 1 sn ht, k .
Sine the optimal speed ontrol is given by (3.10) and (3.17b),
√ (√ )
u2 = µ2 = s1 c sn ht, k ,
we may now integrate x and y dynamis to obtain position variables as funtions
of time. We have,
∫ t
x(t) = u2 cos(θ(t))dt
0
∫ t √ (√ ) (√ )
= s1 c sn ht, k n ht, k dt
0 √
− c 1
[ (√ )]t
= s1√ · dn ht, k
h k2 t=0
s ( (√ ))1
= 1− dn ht, k (3.19)
k
54
Similarly,
∫ t
y(t) = u2 sin(θ(t))dt
0
∫ t√ (√ )
= c 2sn ht, k dt
√0
c 1 [√ ( (√ ) )]t
= √ · ht− E am ht, k , k
h k2 t=0
1 (√ ( (√ ) ))
= ht−E am ht, k , k (3.20)
k

Studying the extremal trajetories in η = 0 ase is important sine any other
trajetory an be expressed by means of these trajetories after a suitable trans-
lation and rotation. This is demonstrated in the following proposition.
Proposition 3.2.2. Assume 0 ≤ c < h. Let us denote an extremal trajetory
belonging to η = 0 ase by (x0(t), y0(t), θ0(t)) aording to (3.18). Then any
extremal trajetory generated by (3.17) an be expressed as,
     
x(t) x (t+ η) x (η)
   0   0 
  = R(−ψ) −   , (3.21)
     
y(t) y0(t+ η) y0(η)
θ(t) = θ0(t + η)− θ0(η), (3.22)
where R(·) is the planar rotation matrix and ψ is dened as,
(√ )
ψ := s · am h η = θ0(η). (3.23)
55
Proof. Indeed, equation (3.22) is readily ahieved by integrating (3.17a). Further,
notie that
√ (√ ) ( (√ ) )
ẋ(t) = s1 c sn h(t+ η), k cos s1 · am h(t+ η), k − ψ
√ ( (√ ) (√ )
= c s1 · cos(ψ) · sn h(t+ η), k n h(t+ η), k
(√ ))
+ sin(ψ) · 2sn h(t + η), k ,
and,
√ (√ ) ( (√ ) )
ẏ(t) = s1 c sn h(t + η), k sin s1 · am h(t+ η), k − ψ
√ ( (√ ) (√ )
= c −s1 · sin(ψ) · sn h(t+ η), k n h(t + η), k
(√ ))
+cos(ψ) · 2sn h(t+ η), k .
Compatly, we have
     
ẋ(t) cos(ψ) sin(ψ) ẋ0(t + η)     
  =   ·   . (3.24)
     
ẏ(t) − sin(ψ) cos(ψ) ẏ0(t + η)
Integration of (3.24) yields (3.21). 
Case 2: c > h
This ase admits solutions analogous to those of Case 1 and we will proeed in a
similar way. To start, we write down solutions to (3.13) as,
√
(√ )
µ1 = h n c(t+ η), k , (3.25a)
√
(√ )
µ2 = s2 h sn c(t + η), k , (3.25b)
√ (√ )
µ3 = s2 c dn c(t+ η), k , (3.25)
56
2 h
with the modulus k = . s2 ∈ {1,−1} and η ∈ R are similar parameters as inc
ase 1. These solutions are periodi as well, with period given by T := 4K√(k)2 .c
Proposition 3.2.3. Let h > 0 and c > h. Assume η = 0, then any extremal
trajetory is given by,
( (√ ))
x(t) = s2 · k 1− n c t, k (3.26a)
(√ ( (√ ) ))
y(t) = s2 c t− E am c t, k , k (3.26b)
θ(t) = cos−
( (√ )) ( (√ ))
1
dn c t, k = sin−1 k sn c t, k , (3.26)
where E(·, ·) denotes the inomplete ellipti integral of the seond kind.
Proof. The proof is similar to that of Proposition 3.2.1 and uses elementary inte-
grals of Jaobi ellipti funtions. 
Analogous to Proposition 3.2.2, we have the following result.
Proposition 3.2.4. Assume h > 0 and c > h. Let us denote an extremal traje-
tory belonging to η = 0 ase by (x0(t), y0(t), θ0(t)) aording to (3.26). Then any
extremal trajetory generated by (3.17) an be expressed as,
     
x(t) x0(t+ η) x0(η)     
  = R(−ψ) −   , (3.27)
     
y(t) y0(t+ η) y0(η)
θ(t) = θ0(t + η)− θ0(η), (3.28)
where R(·) is the planar rotation matrix and ψ is dened as,
− ( (√ )) ( (√ ))ψ := cos 1 c η, k = sin−1dn k sn c η, k = θ0(η). (3.29)
57
Proof. The proof is essentially the same to that of Proposition 3.2.2. Note the
dierent denition of ψ in (3.29) 
Case 3: c = h
2
This ase is transitional between ase 1 and ase 2. Putting the modulus k =
c = 1 in either of the equations (3.17) or (3.25), we get the following solutions
h
√ (√ )
µ1 = s1 c seh c(t + η) , (3.30a)
√ (√ )
µ2 = s2 c tanh c(t + η) , (3.30b)
√ (√ )
µ3 = s3 c seh c(t + η) , (3.30)
with s1, s2, s3 ∈ {1,−1} and η ∈ R. We readily obtain the following result.
Proposition 3.2.5. Let c = h > 0. Assume η = 0, then any extremal trajetory
is given by,
( (√ ))
x(t) = s2 1− seh c t , (3.31a)
(√ (√ ))
y(t) = s1s2 c t− tanh c t , (3.31b)
( (√ ))
θ(t) = s1 tan
−1 sinh c t . (3.31)
Proof. Indeed, integrating θ̇ = u1 = µ1 yields (3.31). Sine u1 = µ2, we get
√ (√ ) (√ )
ẋ = s2 c seh c t tanh c t ,
√ (√ )
ẏ = s1s2 c tanh
2 c t ,
whih in turn gives (3.31a)(3.31b) after integration. 
58
In the same spirit as before, we write any general extremal trajetories in terms
of these trajetories.
Proposition 3.2.6. Assume c = h > 0. Let us denote an extremal trajetory
belonging to η = 0 ase by (x0(t), y0(t), θ0(t)) aording to (3.31). Then any
extremal trajetory generated by (3.30) an be expressed as,
     
x(t) x (t+ η) x (η)
   0   0 
  = R(−ψ) −   , (3.32)
     
y(t) y0(t+ η) y0(η)
θ(t) = θ0(t + η)− θ0(η), (3.33)
where R(·) is the planar rotation matrix and ψ is dened as,
ψ := s tan−
( (√ ))
1
1 sinh c η = θ0(η). (3.34)
This ase onsists of two types of solutions: a straight line solution, orrespond-
ing to the subase where the nal ondition lies on the x-axis, and an asymptoti
solution, whih asymptotially approahes a straight line with slope cotψ. See
Fig. 3.2 for a graphi.
3.2.3 On Time-optimality
From an engineering perspetive, it seems appealing to onsider a similar problem
where nal time T is free. Therefore we want to reah from initial onguration
g 2 20 to the nal onguration gT in a minimal time so that the ontrol ost (u +v )
is minimized along the optimal trajetory. It also makes sense to add a penalty
59
to the time it takes for the uniyle to omplete the maneuver. The ost funtion
an be expressed as
∫
1 T
min (a+ u(t)2 + v(t)2)dt, (3.35)
T,u(·),v(·) 2 0
for some time-penalty parameter a > 0. Note that without the penalty on time, a
solution ould orrespond to a prohibitively large nal time whih is not desired
from a pratial viewpoint. This time-optimal version is just a speial ase of
what we have onsidered in setion 3.2. To see this, we ompute the Hamiltonian
1 2 2
(3.11) as, H = (µ1 + µ2 − a). Sine time-optimality requires the Hamiltonian to2
2 2
be identially zero, we have speial ase of h = µ1 + µ2 = a (.f. (3.15)). This
√
also gives the bounds of the optimal ontrols u1 = µ1, u2 = µ2 to be within ± a.
We an, therefore, use the parameter a to set a desired bound on the ontrols
that is permitted by physial onstraints. A losely related problem was studied
in [Halder and Kalabi, 2017℄, where the speed ontrol v was assumed to be of
onstant magnitude, and the minimum time problem assoiated with minimum
urvature path was onsidered.
3.3 Optimal Control of a Colletive of Agents
Now we onsider a olletive of N agents moving on the plane. Motion of eah
agent an be modeled by the uniyle dynamis (3.1). As seen before, this dynam-
is an be equivalently expressed as a ontrolled dynamis in SE(2), ġk = gkξk(uk),
60
where
ξk(uk) = uk1X1 + uk2X2, k = 1, 2, ..., N. (3.36)
We suppose these agents interat among themselves direted by a xed adjaeny
matrix A = [a ] ∈ RN×Nij . aij = 1 if agent i and agent j interat and aij = 0
otherwise. The nature of the interation will be made preise shortly. LetD be the
degree matrix, i.e. the diagonal matrix where the i-th diagonal entry represents
the number of agents the i-th agent interat with. Then the graph Laplaian is
dened by B := D − A. In this setup, we seek to minimize
∫ T
L = L(ξ1(u1(t)), ..., ξN(uN(t)))dt
0
( )
∫
1 T
N N N
∑ ∑∑
= |ξ 2 2k| + χ akj |ξk − ξj| dt, (3.37)
2 0 k=1 k=1 j=1
for some onstant χ ≥ 0 and the xed endpoint onditions gk(0) = gk0, gk(T ) =
√
gkT , k = 1, ..., N T. Note that we used the trae norm |ξ| = tr(ξ ξ). The param-
eter χ is alled a oupling onstant sine it ats as a weight to the seond term
in the ost funtional (3.37). Without the oupling term, this problem simplies
to solving N opies of the single agent problem as onsidered in Se. 3.2. The
oupling term penalizes agent k through the `mismath in ontrol' with the agents
that it is interating with (i.e. nonzero entries of k-th row of the matrix A). This
type of ost funtional is aimed to apture the `allelomimeti behavior' or the
tendeny to opy neighbors in a natural olletive. [Justh and Krishnaprasad,
2015b℄ studies a very similar problem where the speeds of eah agent is assumed
to be onstant. Here the speed ontrols (uk2, k = 1, 2, ..., N) are to be determined
61
by solving the optimal ontrol problem (3.37).
Sine the underlying optimal ontrol problem is essentially the same, we will
use the results from [Justh and Krishnaprasad, 2015b℄ to derive rst order opti-
mality onditions using the Pontryagin's Maximum Priniple and the Lie-Poisson
∗
redution tehnique. Denote µk ∈ se (2) by,
3
∑
µ = µ X♭k ki i , (3.38)
i=1
[ ]T
♭ ♭ ♭
where {X1, X2, X3} is dual basis to {X1, X2, X3}. If we dene µ̃k = µk1 µk2 µ ,k3
rst order optimality (PMP) yields
   
u µ̃
 1  1
   
   
.
  = Ψ . . . , (3.39)
. .
   
   
   
uN µ̃N
where
Ψ = ((IN + 2χB)⊗ −1I2) = (IN + 2χB)−1 ⊗ I2. (3.40)
∗ N
Here ⊗ denotes the Kroneker produt. The redued Hamiltonian in (se (2))
takes the form
 
µ̃
 1
[ ]  
1  
h = · · · Ψ .µ̃T µ̃T  . . (3.41)2 1 N . 
 
 
µ̃N
The Lie-Poisson redued dynamis is then expressed as follows. Dene µ =
62
[ ]T
µT · · · µT . Then,1 N
µ̇ = Λ(µ)∇h, (3.42)
where Λ(µ) = −diag(Ω1, · · · ,ΩN), with
 
0 µk3 −µ k2
 
 
Ωk = −µ  , k = 1, ..., N. (3.43)
 k3 0 0 
 
 
µk2 0 0
[ ]T [ ]T
Also we have ∇h = (∇h) · · · (∇h) , (∇h)k = ∂h ∂h1 N 0 , k =∂µk1 ∂µk2
1, ..., N . Note that along with the Hamiltonian (3.41), there are N Casimirs
whih are onstants of motion. The Casimirs are dened as,
c 2 2k = µk2 + µk3. (3.44)
In both the extreme ases (i) no oupling (χ = 0) and (ii) high oupling (χ → ∞),
the olletive optimal problem simplies to studying the single agent problem as
onsidered in Setion 3.2. The χ = 0 ase is immediate. The details of the high
oupling limit χ → ∞ is worked out in [Justh and Krishnaprasad, 2015b℄ (see
setion 3()). Dening the quantities
N N N
1 ∑ 1 ∑ 1 ∑
α1 = µj1, α2 = µj2, α3 = µj3, (3.45)
N N N
j=1 j=1 j=1
we obtain the following dierential equations
α̇1 = −α2α3,
α̇ = α (3.46)2 1α3,
α̇3 = −α1α2.
63
These are the same equations we already obtained and analyzed in details for a
single agent ase (3.13).
3.4 Conluding Remarks
In this hapter, we studied an optimal ontrol problem of a uniyle on the plane
in detail. First order neessary onditions are obtained by using the Pontryagin's
Maximum Priniple and the Lie-Poisson redution tehnique. All possible mo-
tion types are properly ategorized by the relative values of the Casimir and the
Hamiltonian. In the later part, we presented a framework for studying a lass
of optimal ontrol problems involving many agents. These agents interat with
eah other by a pre-determined interation graph. The interation enters into the
optimal ontrols of the agents through the additive `ontrol-mismath' term in
ost funtional. This type of ost has been used in literature [Justh and Krish-
naprasad, 2015a,b℄ to apture the `allelomimeti behavior' in natural oks. The
single agent ase, onsidered in this hapter, emerges naturally in a synhroniza-
tion limit of the olletive model. This olletive framework gives us the starting
point to oneptualize a ontinuum ok where we study the limiting ase of
N → ∞ under a spei interation graph. These topis are desribed in detail
in later hapters of this thesis.
64
Part II
Analysis of Colletive Motion:
Top-down Approahes
65
Chapter 4
Continuum Floking and Control
4.1 Motivation
It is a ommon pratie in lassial mehanis to onsider a ontinuous desription
of a physial system. Appliations inlude vibrating rods, vibrating membranes,
uid mehanis et. [Goldstein et al., 2001; Chorin et al., 1990℄ The transition
from disrete-partile system to a ontinuum enables a ompat desription of
the system, often leading to partial dierential equations that reveal deep insights
into the system, e.g. wave-like phenomena, whih may be too obsure or inelegant
in the disrete ounterpart. In the same spirit, we attempt to oneptualize a
ontinuum ok and address its optimal maneuvering properties. Biologial oks
are known to show remarkable response to predator attaks. In the ase of an
attak, the whole ok seems to divert away from the predator. They an perform
these tasks my means of propagating information (in this ase, threat) through
the ok at a muh higher speed than the oking speed [Attanasi et al., 2014℄.
66
A key goal of our approah is to apture this phenomenon, i.e. to unover the
wave-like aspets of oking. An optimal ontrol problem for a ok of nite
agents is presented in Setion 3.3. We take the same framework and formulate a
general ontinuum version of the problem.
Consider an ensemble of N idential, self-steering partiles, eah obeying a
drift free left invariant system on a matrix Lie group G,
ġi = giξi, gi ∈ G, ξi ∈ g, i = 1, ..., N, (4.1)
where g is the Lie algebra of G. These partiles interat with eah other by
a pre-dened graph of interation. The olletive behavior of suh a system was
extensively studied in [Justh and Krishnaprasad, 2015b℄ by solving an appropriate
optimal ontrol problem. Borrowing the notations of [Justh and Krishnaprasad,
2015b℄, we write the ost funtional as self-energy term oupled with a mismath
in steering term
∫ T
J = L(ξ1(t), ..., ξN(t))dt, (4.2)
0
where,
( )
N N N
1 ∑ ∑∑
L(ξ1, ..., ξ ) = ‖ξ ‖2N i + χ aij ‖ξi − ξ 2j‖ , (4.3)
2
i=1 i=1 j=1
where binary valued aij 's populate the adjaeny matrix that denes the graph
of interation and χ ≥ 0 is a oupling onstant. Note that the inner produt
√
〈ξ, η〉 = (ξTtr η), and the orresponding trae norm, ‖ξ‖ = 〈ξ, ξ〉 are in eet,
where ξ, η ∈ g.
67
In an attempt to extend this view, we onsider an innite number of partiles,
i.e. the limiting ase of N → ∞. Here we onsider a one-dimensional ontinuum
1
of partiles, i.e. eah partile is labeled by a point on a irle S . This way, the
agents are thought as a virtual lament. Moreover, we onsider a yli interation
graph, i.e. eah partile is thought to be interating with the `next' partile on
the irle. We introdue the maps, g : R × S1 → G and ξ : R × S1 → g. The
mismath in steering term an then be written as the gradient of ξ in the limiting
ase. In other words, in ontinuum limit of N → ∞, the Lagrangian in (4.3) take
the form,
( )
∥ ∥
∫ 2
1 2π ∥∂ξ(t, θ)∥
L(ξ) = ‖ξ(t, θ)‖2 + χ ∥ ∥ dθ. (4.4)
2 ∥ ∂θ ∥0
Note that the summations over the number of partiles in (4.3) have been replaed
by integral over the irle in the ontinuum setting in (4.4).
Let n be the dimension of the Lie algebra g and {A1, A2, ..., An} denote an or-
1
thonormal basis of g. We introdue the ontrols ui : R × S → R, i = 1, ..., m,
so that ξ an be written as,
m
∑
ξ(t, θ) = ui(t, θ)Ai, (4.5)
i=1
where m < n. With this substitution, the Lagrangian in (4.4) an be rewritten
as,
( )
m ( )∫ 2
1 2π ∑ ∂ui
L(u1, ..., um) = u
2 + χ dθ. (4.6)
2 i0 ∂θi=1
68
Finally, we attempt to minimize the ost funtional,
∫ T
J = L(u1, ..., um)dt, (4.7)
0
( )
m
∂g(t,θ) ∑
subjet to the group dynamis, = g(t, θ)ξ(t, θ) = g(t, θ) ui(t, θ)Ai and∂t
i=1
the xed end-point onstraints, g(0, θ) = g0(θ) and g(T, θ) = gT (θ).
We will note that this optimal ontrol problem an be ast in a more onvenient
setting of loop groups, the group of smooth funtions from the irle to the Lie
group G. In Setion 4.2, we will develop a general framework for suh optimal
ontrol problems in loop group setting. Controllability results will be disussed
in Setion 4.3. This helps us to desribe optimal ontrol solutions in Setion 4.4.
Neessary onditions will be derived by both alulus of variations and Pontrya-
gin's Maximum Priniple approah. An example of ontinuum of nonholonomi
integrators will be studied in detail in Setion 4.5. Setion 4.6 will present deriva-
tion of optimal ontrol equations in the SE(2) ase along with their numerial
treatment. This is a joint work with Dr. E. Justh [Halder et al., 2019a℄.
4.2 A Control System on a Loop Group
Let G be a nite dimensional matrix Lie group and g be its Lie algebra of di-
1
mension n. We will study spaes of smooth maps from the irle S to G and
g,
∞ 1
G = C (S ;G), L = C∞(S1; g).
69
We an onstrut Sobolev ompletions of G and L as done in [Krishnaprasad
et al., 1983℄. We an always view the Lie algebra g as a subalgebra of the general
linear algebra gl(r,R) ∞ 1for some r > n. Dening the spae R = C (S ; gl(r,R)),
we have
G ⊂ R, L ⊂ R.
For any f ∈ R, use the Sobolev k-norm (k ≥ 1)
k ∣ ∣∫ l 2
∑
∣ ∣
‖ ‖ df ∣ ∣k = f(θ) dθ, (4.8)
∣dθl ∣
1 l=0S
where
|f |2 = tr(fTf).
Let the ompletions of G , L , and R in this norm be denoted as Gk, Lk, and Rk,
respetively. By Proposition 3.1 of [Krishnaprasad et al., 1983℄, Gk is atually a Lie
group under pointwise multipliation operation (g1g2)(θ) = g1(θ) · g2(θ), g1, g2 ∈
G 1k, θ ∈ S . Moreover, Lk is the Lie algebra of Gk under pointwise Lie braket
dened as [f1, f2](θ) = [f1(θ), f2(θ)], f1, f2 ∈ Lk, θ ∈ S1. The spaes Gk and Lk
are alled loop groups and loop algebras [Pressley and Segal, 1986℄.
Similar to the nite dimensional Lie groups, we introdue (pointwise) left ation
by Lg : Gk → Gk, h 7→ gh, the left translation by g ∈ Gk. The tangent map of Lg
is then given as ThLg : ThGk → TghGk. We now dene a left invariant vetor eld
on Gk as follows. A vetor eld X : Gk → TGk, h 7→ X(h) is alled left invariant if
ThLg(X(h)) = X(gh), ∀h ∈ Gk. (4.9)
70
Reognizing the Lie algebra Lk as the tangent spae at identity e of Gk (e = {f ∈
Gk|f(θ) ≡ eG}, where eG is the identity element of G), i.e. Lk = TeGk, we an
dene a left invariant ontrol system as,
dg(t)
= TeLg(t)(ξ(t)) = g(t) · ξu(t), (4.10)
dt
where a given ontrol input u(t) determines a ontrolled vetor ξu(t) in the Lie
algebra Lk. Note that the loop algebra Lk an be identied with the tensor
∞ 1
produt spae g ⊗ F , where F is the ring of real valued C funtions on S .
Choose a basis of g as {A1, A2, . . . An}. Then, any ξ ∈ Lk an be written as,
ξ(θ) = ξ1(θ)A1 + · · ·+ ξn(θ)A 1n, θ ∈ S ,
where eah of ξk's (k = 1, . . . , n) are smooth funtions on the irle. We will now
limit ourselves to the study of ontrol vetors ξu of the form,
ξu(t) = u1(t)A1 + · · ·+ um(t)Am, (4.11)
where m < n and the ontrol input u(·) = (u1(·), . . . , um(·)) belongs to the set U
of pieewise ontinuous U valued funtions, where U is vetor spae of Rm valued
smooth funtions on the irle, i.e. U := {u(·) : u is pieewise ontinuous in t, u(t) ∈
U = C∞(S1;Rm)}.
4.3 Controllability
Having onstruted the ontrol system on the loop group Gk, it is natural to ask
the question of ontrollability or aessibility, i.e. given any two points g1 and
71
g2 in Gk, if they an be onneted by a pieewise dierentiable urve, onsisting
of possibly nitely many piees, eah piee being an integral urve of a left in-
variant vetor eld dened by hoosing a ontrol u(·) ∈ U . In nite dimensional
1
analogue of this question, i.e. where we shrink the irle S down to a point, the
ontrollability question is answered by the well known Chow-Rashevsky theorem
[Wei-Liang, 1939; Rashevsky, 1938℄. In innite dimensional ases, however, it is
not immediate if the Chow-Rashevsky theorem remains valid. There is a body of
literature [Heintze and Liu, 1999; Salehani and Markina, 2014℄ that attempts to
attak this problem. It is the result of [Heintze and Liu, 1999℄ that we use in this
setion. This result addresses the ontrollability question in a weaker sense whih
we will make expliit.
Let M be a omplete onneted Hilbert manifold and let X(M) denote the
set of all smooth vetor elds dened on M. Let F ⊂ X(M) be a given family of
smooth vetor elds onM. LetRF(x) be the set of points inM that an be joined
from x ∈ M by means of a pieewise dierentiable urve, eah piee of whih is
an integral urve of a vetor eld in F . Let Lie F be the Lie subalgebra of X(M)
generated by F , and Liex F = {X(x) : X ∈ Lie F} - the evaluation of Lie F at
x ∈ M. If M is nite dimensional, the lassial Chow-Rashevsky theorem holds:
if Liex F = TxM for eah x ∈ M, then RF (x) = M, for every x ∈ M. In a
general Hilbert manifold M, the following generalized Chow-Rashevsky theorem
holds:
72
Theorem 4.3.1 ([Heintze and Liu, 1999℄). Let M be a omplete onneted Hilbert
manifold and F a family of smooth vetor elds dened on M. If Liex F is dense
in TxM for all x ∈ M, then RF (x) is dense in M for all x ∈ M.
Theorem 4.3.1 is a weaker statement than the one for nite dimensional ase.
Here we make preise the strong and weak notions of ontrollability. Consider
the ontrol system onstruted in (4.10)(4.11). Note that the loop group Gk
an be given a struture of a smooth Hilbert manifold [Eells Jr, 1966; Ebin and
m
Marsden, 1970℄. In this ase, the family F ∈ X(Gk) is given by {Xi}i=1 , where
Xi(g(t)) = g(t) · (ui(t)Ai), for g(t) ∈ Gk.
Denition 4.1. (Strong Controllability) The ontrol system (4.10)(4.11) is
said to be strongly ontrollable if RF = Gk, i.e. given any two points g1, g2 ∈ Gk,
we an nd a ontrol input that will transfer the system from g1 to g2.
Denition 4.2. (Weak Controllability) The ontrol system (4.10)(4.11) is
said to be weakly ontrollable if RF is dense in Gk, i.e. given any two points
g1, g2 ∈ Gk, we an nd a ontrol input that will transfer the system from g1 to a
state that is arbitrarily lose to g2.
The set {A1, ..., Am} is said to be braket generating if the iterated brakets
of its elements span the Lie algebra g. In the nite dimensional analogue of
the ontrol system dened in (4.10)(4.11), the (strong) ontrollability ondition
m
aording to Chow-Rashevsky theorem is equivalent to having the set {Al}l=1
braket generating in g. We will now try to establish (weak) ontrollability of
73
the innite dimensional loop group ase by assuming that {A1, ..., Am} is braket
generating in g.
Theorem 4.3.2. Consider the ontrol system (4.10)(4.11) on the loop group Gk.
m
Assume that the set {Al}l=1 is braket generating in g. Then the system is weakly
ontrollable.
Proof. The denition of left invariant vetor elds on Gk (4.9) an also be made
expliit by means of smooth funtions on Gk. Let D be the set of smooth real val-
ued funtions on Gk. Then given an element ξ ∈ Lk, we an dene a dierentiable
vetor eld Xξ : D → D as,
(Xξf)(g) = (Df)g · gξ, f ∈ D, (4.12)
whereD denotes the dierential operator. Given two vetor eldsXξ, Xη ∈ X(Gk),
we an alulate their Jaobi-Lie braket dened as,
[Xξ, Xη]f = Xξ(Xηf)−Xη(Xξf), f ∈ D.
We ompute,
Xξ(Xηf)(g) = (Xξ((Df)g · gη))(g)
= D((Df)g · gη)g · gξ
= (D2f)g · (gη, gξ) + (Df)g · (D(gη)g · gξ)
= (D2f)g · (gη, gξ) + (Df)g · (gξη).
Similarly,
Xη(X
2
ξf)(g) = (D f)g · (gξ, gη) + (Df)g · (gηξ).
74
2
The symmetry of the seond dierential operator D yields
[Xξ, Xη]f(g) = (Df)g · (g(ξη − ηξ))
= X[ξ,η]f(g),
where [ξ, η] is the usual (pointwise) Lie braket on Lk. This leads us to a detailed
study of the Lie braket of the loop algebra Lk. It is immediate that Lk is
mr
generated by the generators {Pr }m ∈Z, r∈{1,...,n} dened as,r
Pmr imrθr = e Ar. (4.13)
r
Let the struture onstants of g be denoted by Γpq. Then,
n
∑
[Pmp, Pmq ] = Γr mp+mqp q pqPr . (4.14)
r=1
With this notation, a ontrolled vetor ξ(t) = ξu(t) ∈ Lk an be expressed as
( )
m
∑ ∑
ξ mr mru(t) = ur (t)Pr , (4.15)
r=1 mr∈Z
mr
where for eah r = 1, ..., m, ur (t) ∈ C's are the Fourier oeients of the ontrol
ur(t) and umrr (t) = u
−mr
r (t) (sine the ontrols are real). We now dene a family
mr
of vetor elds on Gk as F = {Xr }m ∈Z, r∈{1,...,m}, wherer
Xmrr f(g) = (Df) · gPmrg r , f ∈ D.
Taking braket of any two vetor elds from the family F yields another vetor
led whih is governed by the ommutator relationship in (4.14). Note that sine
the set {A1, ..., Am} is braket generating in g, for eah l ∈ {m+ 1, ..., n}, we are
guaranteed to get the item
∑
a m
Pml := P r r rl l , a ∈ +r Z ∪ {0}, mr ∈ Z, r ∈ {1, ..., m},
75
at some depth of iterated brakets from the family F . By hoosing mr ∈ Z, r ∈
{1, . . . , m}, we an then ahieve any ml ∈ Z. We have thus proved that if the set
{Ar}r∈{1,...,m} is braket generating in g, Lieg F is dense in the tangent spae TgGk
at eah g ∈ Gk. The generalized Chow's theorem 4.3.1 then provides the required
(weak) ontrollability result. 
4.4 Optimal Control Problems
We start with the left invariant ontrol system on the loop group Gk as in (4.10)-
(4.11). Now for a given T > 0, onsider the following xed end-point optimal
ontrol problem,
∫ T
min J(u) = L(g(t), u(t))dt
u∈U
( 0PG) (4.16)
subjet to: ġ = g · ξ, g(0) = g0, g(T ) = gT ; g0, gT ∈ Gk.
We are interested in deriving neessary onditions for optimality for suh optimal
ontrol problems. Speial are needs to be taken sine the problem is posed in an
innite dimensional setting. We provide two dierent approahes for doing that.
4.4.1 Calulus of Variations Approah
∞ 1 r
Let x(t) = (x1(t), ..., xr(t)) ∈ C (S ,R ) =: X denote a vetor that an be used
to represent the omponents of g(t) ∈ G, for some r ≥ n. The group dynamis
ġ = g · ∑ξ(u) = g · ( i uiAi) eetively lets us write the ontrol u(t) as a funtion
of (x(t), ẋ(t)). The xed endpoint onstraints in g an be translated to some
76
nonholonomi onstraints of the form
∫ T
Φ(x, ẋ) = φ(x(t), ẋ(t))dt = 0, (4.17)
0
where φ(x(t), ẋ(t)) ∈ C∞(S1,Rl) =: Z, for some l < r. Then the problem (PG)
an be written as,
∫ T
min J = L(x(t), ẋ(t))dt, (4.18)
0
subjet to the nonholonomi onstraints (4.17). Here we have to keep in mind that
the variations in x are to be both in t and θ. This is a well known problem alled
the `Lagrange problem' in alulus of variations. The one-dimensional Lagrange
problem is well studied [Gelfand and Fomin, 1963; Elsgol, 2012℄. However, the
theory behind multidimensional problem is more ompliated and less omplete
[Giaquinta and Hildebrandt, 1996; Bliss, 1946℄. The diulty arises sine not all
the ẋ are freely variable. Aording to [Giaquinta and Hildebrandt, 1996, p. 112℄,
there exist a Lagrange multiplier λ ∈ C∞(R× S1;Rl), suh that we an nd the
free extremals of the augmented Lagrangian in an usual way. Moreover, sine the
onstraints (4.17) are of isoperimetri type, λ does not depend on t [Rund, 1966,
p. 349℄. We an write the augmented Lagrangian as
L̃ = L+ 〈λ, φ〉Z , (4.19)
∞ 1 l
where λ = (λ1, λ2, ..., λl) ∈ C (S ;R ).
Furthermore, we are interested in a speial struture of the Lagrangian, namely
∫ 2π
L(x, ẋ) = L(x, ẋ) dθ, (4.20)
0
77
where the funtional L is alled the Lagrangian density. The augmented ost
funtion take the form,
∫ T ∫ 2π l ∫ 2π
∑
J̃(x, ẋ) = L(x, ẋ) dθdt+ λj(θ)φj(x, ẋ) dθ
0 0 j=1 0
∫ T ∫ 2π
= L̃(x, ẋ) dθdt,
0 0
where the redened Lagrangian density is given by,
l
∑
L̃ = L+ λjφj. (4.21)
j=1
We will relabel L̃ by L in subsequent analysis for onveniene. By invoking
∂x ∂x ∂2x
notations xt = , xθ = , xtθ = = x et., we an write L = L(x, x , x∂t ∂θ ∂t∂θ θt θ t, xtθ).
In order to optimize this ost funtional, the variational priniple requires that,
∫ T ∫ 2π r ( )
∑ ∂L ∂L ∂L ∂L
δJ̃ = δxi + δxi,θ + δxi,t + δxi,tθ dθdt = 0,
0 0 ∂xi=1 i ∂xi,θ ∂xi,t ∂xi,tθ
where δy denotes variation of the quantity y that vanishes at the endpoints of t
and θ. Using integration by parts, for eah i, we may write,
∣2π
∫ 2π ∂L ( )∫∂L ∣ 2π− ∂ ∂L∣δxi,θdθ = δxi δxidθ∣
0 ∂xi,θ ∂xi,θ ∣ 0 ∂θ ∂xi,θ
0
∫ 2π ( )
− ∂ ∂L= δxidθ,
0 ∂θ ∂xi,θ
Similarly,
∫ T ∂L ∫ T (∂ ∂L )
δxi,tdt = − δxidt,
0 ∂xi,t 0 ∂t ∂xi,t
∫ T ∫ 2π L ∫ T ∫ 2π ( L )∂ ∂2 ∂
δxi,tθdθdt = δxi dθdt.
0 0 ∂xi,tθ 0 0 ∂t∂θ ∂xi,tθ
78
Hene,
T 2π ( L ( L ) (∫ ∫ ∑ ∂ − ∂ ∂ − ∂ ∂L
) ( ))
∂2 ∂L
δJ̃ = + δxi dθdt
0 0 ∂xi ∂θ ∂xi,θ ∂t ∂xi,t ∂t∂θ ∂xi i,tθ
The Euler-Lagrange equations an then be expressed as,
( ) ( ) ( )
∂L ∂ ∂L ∂ ∂L ∂2− − ∂L+ = 0, i = 1, 2, ..., r. (4.22)
∂xi ∂θ ∂xi,θ ∂t ∂xi,t ∂t∂θ ∂xi,tθ
4.4.2 Maximum Priniple Approah
In this setion we provide a brief exposure to Pontryagin's Maximum Priniple
(PMP) type argument in innite dimensional spaes. It is to be noted that PMP
does not automatially hold in general innite dimensional optimal ontrol prob-
lems, one requires some more assumptions for it to work. A detailed study on
this subjet is done in Appendix A. Here we only dene some notations and
assumptions to state the neessary theorem.
We onsider an abstrat dierential equation in a Hilbert spae X ,
dx(t)
= f(t, x(t), u(t)) a.e. in [0, T ], (4.23)
dt
where x(t) ∈ X , u(·) ∈ U , and T > 0. Here X is alled the state spae and U
is the set of all measurable funtions u(·) : [0, T ] → U , where U is a separable
metri spae alled the ontrol spae. With this setup, we formulate the following
optimal ontrol problem (P),
∫ T
min J(u) = L(t, x(t), u(t))dt
u∈U 0
(P) (4.24)
subjet to: ẋ = f(t, x, u) a.e. in [0, T ], x(0) = x0, x(T ) = xT .
79
We assume that both the funtions f(·, ·, ·) and L(·, ·, ·) are Bohner integrable in
t ∈ [0, T ] and Lipshitz ontinuous in x(t) ∈ X , with onstant K. Furthermore,
′
we require the existene and ontinuity of the Fréhet derivatives fx(t, x, u) and
L′x(t, x, u)
′ ′
. We also assume the funtions f, L and their derivatives fx, Lx to be
bounded, i.e. there exists an M > 0, suh that
‖f(t, x, u)‖ ≤ M, ‖f ′x(t, x, u)‖ ≤ M,
‖L(t, x, u)‖ ≤ M, ‖L′x(t, x, u)‖ ≤ M,
for all (t, x(t), u(t)) ∈ [0, T ]×X×U . Note that these hypotheses ensure a ontin-
uous and unique solution of (4.23) to exist [Avez, 1986℄. The following tehnial
details is one of the key ingredients in the proof of the PMP.
Denition 4.3. (Finite Codimensionality) [Fattorini, 1987℄ A subset S of a
Hilbert spae Z is alled to be nite odimensional in Z, if there exists a losed
subspae Zc ⊆ Z of nite odimension suh that Sc = Πc(o(S)), has nonempty
interior in Zc, where Πc denotes the orthogonal projetion from Z onto Zc and o
means losed onvex hull.
We will now make a key assumption to derive a nontrivial maximum priniple.
∗
Let a solution of problem (P) exist and the optimal ontrol be denoted by u ∈ U
∗
and let the orresponding optimal trajetory be denoted as x (·). Then dene the
`reahable set' as,
{
∫ t
′ ∗ ∗
R := z(T ) ∈ X | z(t) = fx(s, x (s), u (s)) · z(s)ds
0
}
∫ t
+ (f(s, x∗(s), v(s))− f(s, x∗(s), u∗(s))) ds, for some v(·) ∈ U (4.25)
0
80
(A1) The set R is nite odimensional in X .
Remark 4.1. In general, it is not lear whether there exists a relationship between
ontrollability (strong or weak) of the system and the nite odimensionality
assumption of the `reahable set' R. We will, however, prove that in a speial
ase of G = H(3), the Heisenberg group, the strong ontrollability implies nite
odimensionality of R. It is of future onsideration to address this question in a
general ase.
Using usual formalism, we invoke the pre-Hamiltonian funtion H : R×X ×
U ×R×X∗ → R as,
H(t, x(t), u(t), p0, p(t)) = p0L(t, x(t), u(t)) + 〈p(t), f(t, x(t), u(t))〉 , (4.26)
p(t) ∈ X∗where is alled the ostate variable. Then the PMP an be written as,
∗
Theorem 4.4.1. (Maximum Priniple) Let u ∈ U be an optimal ontrol for
∗
problem (P) and x (t) be the orresponding optimal trajetory. Then, there exist
(p∗a pair 0, p
∗(t)) ∈ R × X∗, t ∈ [0, T ] ∗ ∗ ∗ ∗, suh that (p0, p ) ≡6 (0, 0), p0 ≤ 0, p (·)
satises the dierential equation,
ṗ∗(t) = − (f ′ (t, x∗ ⋆x (t), u∗(t)) p∗(t)− p∗0L′x(t, x∗(t), u∗(t)), (4.27)
⋆
where by A we denote the adjoint operator of the operator A. The pointwise
maximization of the pre-Hamiltonian holds,
H(t, x∗(t), u∗(t), p∗0, p
∗(t)) = maxH(t, x∗(t), v, p∗, p∗(t)), (4.28)
v∈ 0U
81
for a.e. t ∈ [0, T ] ∗ ∗. Moreover, x and p satisfy Hamilton's anonial equations,
i.e.
dx∗ δH
= ∗ (t, x
∗, u∗, p∗ ∗
dt δp 0
, p ) (4.29)
dp∗ δH
= − ∗ (t, x
∗, u∗, p∗0, p
∗). (4.30)
dt δx
A proof of this theorem is rather ompliated and is given in Appendix A. We
now use this result to state a maximum priniple for the loop group ase.
∗
Theorem 4.4.2. (Maximum Priniple in loop group setting) Let u ∈ U
∗
be an optimal ontrol for problem (PG) and g (t) be the orresponding optimal
trajetory. Assume the nite odimensionality ondition (A1). Denote Rk, the
∞ 1
Hilbert spae of k-Sobolev ompletion of the spae R = C (S , gl(r,R)), for some
r > n (p∗, p∗(t)) ∈ R×R , t ∈ [0, T ] (p∗, p∗. Then, there exist a pair 0 k , suh that 0 ) 6≡
(0, 0), p∗0 ≤ 0, p∗(·) satises the dierential equation
ṗ∗(t) = −p∗(t) · ξ(u∗(t))T, (4.31)
and the pointwise maximization of the pre-Hamiltonian holds,
H(g∗(t), u∗(t), p∗, p∗0 (t)) = maxH(g
∗(t), v, p∗, p∗(t)), (4.32)
v∈ 0U
∗ ∗
for a.e. t ∈ [0, T ]. Moreover, g and p satisfy Hamilton's equations, i.e.
dg∗ δH
= ∗ (g
∗, u∗, p∗0, p
∗)
dt δp
(4.33)
dp∗ −δH= ∗ (g
∗, u∗, p∗, p∗).
dt δg 0
Proof. It is almost immediate that under a nite odimensionality assumption like
(A1), we an state a maximum priniple like Theorem 4.4.1 for problem (PG). The
82
only aveat is that the state spae Gk is not a Hilbert spae and hene Theorem
4.4.1 annot be applied diretly. However, adopting an `enlargement' tehnique
[Brokett, 1973; Justh and Krishnaprasad, 2015a℄, we an state an analogous
maximum priniple. We reognize that the loop group Gk is a subset of Rk. The
spae Rk an then be regarded as the `raised' state spae. The dynamis (4.10),
1
along with the initial ondition g(0, θ) = g0(θ) ∈ G for all θ ∈ S , ensures that
g(t, θ) remains in G for all (t, θ) ∈ [0, T ] × S1. Endow the spae gl(r,R) with
( )
the trae inner produt and an indued norm, i.e. 〈A,B〉 (r, ) = Ttr A B andgl R
√
‖A‖ Tgl(r) = tr (A A), for A,B ∈ gl(r,R). We now dene the pre-Hamiltonian
H : Rk × U ×R× Rk → R as,
H(g(t), u(t), p0, p(t)) = 〈p(t), g(t)ξ(u(t))〉 + p0L(u(t)), (4.34)Rk
where the duality pairing in the denition of H an be expliitly written as,
2π k 〈 〉∫ i
∑
〈 d d
i
p(t), g(t)ξ(u(t))〉 = p(t, θ), (g(t, θ)ξ(u(t, θ))) dθ
Rk dθi i0 dθi=0 gl(r,R)
∫ 2π k ( )i
∑ d di
= tr p(t, θ)T · (g(t, θ)ξ(u(t, θ))) dθ.
0 dθ
i dθi
i=0
∗
We are now all set to apply Theorem 4.4.1. If u ∈ U is an optimal ontrol, then
we have,
H(g∗(t), u∗(t), p∗, p∗0 (t)) = maxH(g
∗(t), v, p∗, p∗0 (t)). (4.35)
v∈U
δH (g∗, u∗, p∗, p∗) = g∗ξ(u∗) = ġ∗It is obvious that 0 . We an also derive for anyδp∗
g̃(t) ∈ Rk, (suppressing other arguments)
δH · H(g
∗ + αg̃)−H(g∗) 〈 〉
∗ g̃ = lim = 〈p
∗, g̃ξ(u∗)〉 = p∗ξ(u∗)T, g̃ ,
δg α→0 α Rk Rk
83
whih implies the adjoint equation to (4.10) is,
ṗ∗(t) = −p∗(t) · ξ(u∗(t))T. (4.36)

4.5 Speial Case : G = H(3)
The previous setion on optimal ontrol provides a onrete foundation in whih
we an state the maximum priniple for the onsidered optimal ontrol problem.
In this setion, we will explore a speial ase where we take the Lie group, G as
the Heisenberg group, H(3). Note that the nite number of partiles ase of this
problem has been onsidered in [Justh and Krishnaprasad, 2016℄ and hene this
work an be thought as a ontinuum ounterpart of it. In H(3), g(t, θ) ∈ H(3) has
the struture,
 
1 x1 x
x1x2
 3 + 2 
 
 
g = 0 1 x  ,
 2 
 
 
0 0 1
that satisfy the group evolution equation,
∂
g(t, θ) = g(t, θ) · (u1(t, θ)A1 + u2(t, θ)A2) , (4.37)
∂t
where,
   
0 1 0 0 0 0
   
   
   
A 1 = 0 0 0 , A2 = 0 0 1 ,
   
   
   
0 0 0 0 0 0
84
along with A3 = [A1, A2], form an orthonormal basis for the assoiated Lie algebra
h(3). We attempt to address the optimal ontrol problem formulated before, under
this Heisenberg group setting, i.e. we solve the following,
∫ T ∫ T ∫ 2π
min L(u)dt = L(u)dθdt
0 0 0
( ( ))
T 2π ( )2 ( )∫ ∫ 21 ( )2 2 ∂u1 ∂u2= u1 + u2 + χ + dθdt,2 0 0 ∂θ ∂θ
(4.38)
where L is alled the Lagrange density funtion.
4.5.1 Controllability
It is a diret exerise of the generalized Chow-Rashevsky theorem 4.3.1 to show
∞ 1
weak ontrollability in this ase. The loop algebra C (S , h(3)) is spanned by
{eim1θA1, eim2θA2, eim3θA3}. The family F of left invariant vetor elds that is
m1 m2
hosen by means of ontrol inputs is given by F = {X1 , X2 }, where
Xmrr f(g) = (Df) · geimrθg Ar, f ∈ D, r = 1, 2.
Sine the only non-vanishing brakets in h(3) are [A1, A2] = A3 = −[A2, A1], the
∞ 1
set Lieg F would span the tangent spae at every point g ∈ C (S ,H(3)).
However, we an provide an argument that establishes the strong ontrol-
lability in this ase. Here we desribe an approah to onstrut a andidate
0 T
smooth ontrol given any endpoint onditions xi (θ) and xi (θ), i = 1, 2, 3. With-
0
out loss of generality, we may assume xi (θ) = 0 for all i. Now sine, xi(t, θ) =
85
∫ t
ui(t, θ)dt, i = 1, 2, we an hoose smooth ontrols vi(t, θ), t ∈ [0, t̄ ], i = 1, 2, for0
some t̄ < T , suh that x1 and x2 reah their nal endpoints. At time t = t̄, let the
T
`error' in x3 variable be denoted as ∆x3(θ) = x3 (θ)− x3(t̄, θ). Note that without
loss of generality, we may assume that ∆x3(θ) > 0 for all θ. We know that in
a single nonholonomi integrator ase, if we omplete a loop in time for (x1, x2)
variables, the hange in x3 variable will be given by the area of the loop. We may
use the same idea in the ontinuum ase to onstrut smooth ontrols. We an
generate the following irular loops (in time) in (x1, x2) variables,
     
x1(t, θ) x
T
1 (θ)− r(θ) cos(ω(t− t̄))     
  =  + r(θ)  ,
     
x T2(t, θ) x2 (θ) sin(ω(t− t̄)) (4.39)
2π
ω = − , t ∈ (t̄, T ],T t̄
where r(·) is a smooth funtion in θ. The ontrols required to generate these loops
are given by
   
ṽ1(t, θ) − sin(ω(t− t̄))   
  = r(θ)ω  , (4.40)
   
ṽ2(t, θ) cos(ω(t− t̄))
We then ompute,
∫
1 T
x3(T, θ) = x3(t̄, θ) + (x1(t, θ)ṽ2(t, θ)− x2(t, θ)ṽ1(t, θ))dt
2 t̄
= xT3 (θ)−∆x3(θ) + πr2(θ).
Sine ∆x3(·) is smooth, we an always hoose smooth funtion r(·) suh that
86
∆x3(θ) = πr
2(θ), and thus the smooth ontrols



vi(t, θ), t ∈ [0, t̄ ]
ui(t, θ) = , i = 1, 2, (4.41)


ṽi(t, θ), t ∈ (t̄, T ]
make the required state transitions possible. This shows that in the Heisenberg
ase we have strong ontrollability.
Remark 4.2. We an assume ∆x3 > 0 for all θ beause if it was not the ase,
we ould add another piee of ontrols v̂i's before applying the ontrols ṽi's. The
purpose of the ontrols v̂i's would be to make the states (x1, x2) undergo a irular
2
loop of radius r̂, for all θ. This will produe a hange in the x3 variable by πr̂
2
for all θ. Hene the new error an be written as ∆x3(θ) = ∆x3(θ) + πr̂ . We an
always pik a r̂ so that ∆x3(θ) > 0 for all θ.
4.5.2 Equations of Optimal Control
Neessary onditions for optimality an be derived by various methods. We will
present two suh approahes to solve the optimal ontrol problem (4.38).
87
4.5.2.1 Calulus of Variations Approah
Note that the group dynamis in (4.37) an also be expressed by the dierential
equations,
∂x1
= u1, (4.42)
∂t
∂x2
= u2, (4.43)
∂t
∂x3 1
= (x1u2 − x2u1). (4.44)
∂t 2
∫ T
Sine we have an integral onstraint, namely, ∆x3(θ) =
1 (x u
2 0 1 2
− x2u1)dt, we
invoke a Lagrange multiplier, λ ∈ C∞(S1,R) to write the augmented Lagrangian
density funtion as,
1 [( ) ( ) ]L = x2 21,t + x2,t + χ x2 21,θt + x2,θt + λ (x1x2,t − x2x1,t) . (4.45)2
2
where we adopted the notation onventions x = ∂xii,θ , x
∂ xi
∂θ i,tθ
= = x
∂t∂θ i,θt
et.
Appliation of the Euler-Lagrange equations (4.22) gives the following equations,
∂2x1 ∂x2 ∂
4x1
= λ + χ
∂t2 ∂t ∂θ2∂t2
(4.46)
∂2x 42 − ∂x1 ∂ x2= λ + χ ,
∂t2 ∂t ∂θ2∂t2
whih yields the evolution equations for optimal ontrols,
A∂u1 = λu2
∂t
(4.47)
A∂u2 = −λu1,
∂t
where we denote A := (1− χ∆), ∆ being the Laplaian.
∞ 1
Remark 4.3. A is a positive denite self adjoint operator in C (S ,R), having
eigenvalues αn = 1 + χn
2 inθ
with assoiated eigenvetors en = e for n ∈ Z.
88
Lemma 4.5.1. The quantities,
( )
∫ 2π ∂nu ∂nu ∂n n1A 1 u2A∂ u2hn = + dθ, n = 0, 1, 2, ... (4.48)
∂tn ∂tn ∂tn ∂tn0
are onserved along any optimal trajetories satisfying (4.47).
dhn
Proof. It is easy to establish that for eah n, = 0, whih follows diretly from
dt
the way optimal ontrols behave in (4.47) and the fat that the operator A is self
adjoint. 
4.5.2.2 PMP Approah
∞ 1 3
We introdue the ostate variable p(t) = (p1(t), p2(t), p3(t)) ∈ C (S ;R ), t ∈
[0, T ]. The pre-Hamiltonian an be written as (onsidering only normal extremals,
i.e. where p0 6= 0 and an be normalized to −1),
H = 〈ẋ, p〉 − L
2π ( )∫ 1 1 ( ) χ ( )
= u1p1 + u2p2 + (x1u2 − x 2 2 2 22u1)p3 − u + u − u + u dθ.
0 2 2
1 2 2 1,θ 2,θ
(4.49)
Remark 4.4. Note that the nite odimensionality assumption is satised in this
ase. To see this, note that the omponents of members of the `reahability set'
(4.25) an be expressed as,
∫ T
z1(T ) = w1(t)dt
0
∫ T
z2(T ) = w2(t)dt
0
∫
1 T
∫
1 T
z3(T ) = (z1(t)u
∗
2(t)− z2(t)u∗1(t)) dt+ (x∗1(t)w ∗2(t)− x2(t)w1(t)) dt,2 0 2 0
89
where w(t) = (w1(t), w2(t)) is any arbitrary ontrol input. If we denote trajetories
u
orresponding to any input u(·) as x (t) = x(t, u(·)) z (T ) = xw, then, 1 1 (T ), z2(T ) =
xw2 (T ). Furthermore note that,
[∫ ∫
1 T T
z3(T ) = (x
∗u∗1 2 − x∗ ∗2u1)dt+ (xw w1 w2 − x2 w1)dt2 0 0
∫ T ]
− ((x∗ − xw)(u∗ − w )− (x∗ w1 1 2 2 2 − x2 )(u∗1 − w1))dt
0
= x∗3(T ) + x
w
3 (T )− x̃3(T ),
∗ w ∗
where x̃ = x − x . We may now hoose w = u . This makes x̃(T ) = 0
and (z1(T ), z2(T ), z3(T )) = (x
∗
1(T ), x
∗
2(T ), 2x
∗
3(T )). Sine the Heisenberg ase is
strongly ontrollable, the `reahability set' R spans the whole of the state spae
X , making it trivially nite odimensional in X .
We an now diretly apply Theorem 4.4.1 to derive neessary optimality on-
ditions. The maximum priniple would require us to maximize (4.49) pointwise
over the ontrols, i.e. we are attempting to nd the Hamiltonian as,
H(x, p) = sup H(x, {v }2i i=1, p)
v ∞ 11,v2∈C (S ;R)
(
∫ 2π 1
= sup v1p1 + v2p2 + (x1v2 − x2v1)p3
v ,v ∈C∞ 1 21 2 (S ;R) 0
)
1 ( ) ( )− 2 2 − χv1 + v 2 22 v1,θ + v2,θ dθ (4.50)2 2
This maximization results in two Euler-Lagrange equations,
( )
∂H − ∂ ∂H = 0, i = 1, 2, (4.51)
∂vi ∂θ ∂vi,θ
90
that yield the optimal ontrols,
Au1 = p1 −
1
x2p3,
2
(4.52)
A 1u2 = p2 + x1p3,
2
with the usual denition of A = (1− χ∆). The Hamiltonian an be read as,
∫
1 2πH = (u1Au1 + u2Au2) dθ
2 0
[( ) ( )
∫
1 2π − 1 1= p1 x −12p3 A p1 − x2p3
2 0 2 2
( ) ( )]
1
+ p2 + x1p3 A−1
1
p2 + x1p3 dθ (4.53)
2 2
The dynamis of the ostate variable p an be alulated from Hamilton's equa-
∂p = − δH δHtion, , where denotes the funtional derivative of H with respet
∂t δx δx
to x. Note that the Hamiltonian funtion H is smooth in x, so we an take this
〈 〉
δH
derivative. This an be dened as, DH(xi) · σ = , σ , i = 1, 2, 3, whereδxi
DH(xi) · σ is the Fréhet derivative of H at xi in the diretion of σ. This is
dened as,
∣
d ∣
DH(x ∣i) · σ = H(x1 + ǫσ) .
dǫ ∣ǫ=0
We may alulate, for i = 1,
∫
1 d 2π
[( ) ( )
H · − 1D (x ) σ = p −1 11 1 x2p3 A p1 − x2p3
2 dǫ 0 2 2
( ) ( )] ∣
1 1 ∣
+ p2 + (x1 + ǫσ)p A−13 p2 + (x1 + ǫσ)p ∣3 dθ
2 2 ∣ǫ=0
( ( ))
∫
1 2π 1
= σp A−13 p2 + (x1 + ǫσ)p3 dθ
2 0 2
2π ( ( ))∫1
= σp3A−1
1
p2 + x1p3 dθ
2 0 2
〈 ( ) 〉
p3
= A−1 1p2 + x1p3 , σ , (4.54)
2 2
91
so that we may write,
H ( )δ p3A−1 1= p2 + x1p3 . (4.55)
δx1 2 2
( )
δH
We an similarly alulate = −p3A−1 p − 11 x2p δH3 and = 0. The evolutionδx2 2 2 δx3
of p an then be expressed as,
( )
∂p1 −p3= A−1 1p2 + x1p3
∂t 2 2
( )
∂p2 p3
= A−1 1p1 − x2p (4.56)3
∂t 2 2
∂p3
= 0.
∂t
From (4.52) and (4.56), we notie that,
∂ Au1 = −p A−13 u2
∂t
(4.57)
∂ Au −12 = p3A u1
∂t
Hene, reognizing p3 as the negative of Lagrange multiplier λ in the previous
setion, we redisover (4.47).
4.5.3 Behavior of Optimal Control
Denote, z(t, θ) := u1(t, θ) + iu2(t, θ). Then, (4.47) an be expressed as,
A∂z = −iλz. (4.58)
∂t
Sine u1, u2 are periodi funtions in θ with period 2π, they have a Fourier series
∞
∑ (ν) (ν)
representation, uν(t, θ) = un (t)e
inθ, ν = 1, 2, where un 's are the Fourier
n=−∞
92
oeients of uν . Moreover, sine uν 's are real valued funtions of θ, we have,
(ν) (ν)
u−n = un . Transforming (4.58) in Fourier domain, we get,
dž
A = −iΛž, (4.59)
dt
( ∞ ) ∞ (1)
where we denote, A = diag {αn}n=−∞ , ž(t) = [zn(t)]−∞. zn(t) = un (t) +
(2)
iun (t), i.e. zn(t)'s are Fourier oeients of z(t, θ). Λ is the innite Toeplitz
matrix generated by the Fourier oeients of λ, i.e.
 
λ λ λ
 0 −1 −2 · · ·
 
 
 λ 
 1 λ0 λ−1 · · ·
Λ =   (4.60)
 
 λ λ λ · · ·
 2 1 0 
 
 
· · · · · · · · · · · ·
Sine λn's are Fourier oeients of real valued funtion λ, we have λ−n = λn.
∗
This leads to the observation that Λ is (innite) Hermitian matrix, i.e. Λ = Λ .
4.5.3.1 Trunation of Fourier Coeients
Here we will onsider rst N+1 Fourier oeients of z and provide an analysis of
T
(4.59) in the trunated nite dimensional ase. We write, žN = (z−N · · · z0 · · · zN ) ∈
2N+1
C . Then, the trunated version of (4.59) an be written as,
džN
AN = −iΛN žN , (4.61)
dt
where AN and ΛN 's are appropriately trunated matries from A and Λ, respe-
1/2 −1
tively. Note that AN ≻ 0. Let us denote, ẑ = AN žN and BN = AN , then we
get,
93
dẑ
= −iΛ̂N ẑ, (4.62)
dt
( )∗
1/2 1/2
where Λ̂N = BN ΛNBN = Λ̂N . Sine the matrix −iΛ̂N is skew Hermitian,
all its eigenvalues are on the imaginary axis. Denote them by −iσn, σn ∈ R, for
n = −N, ..., N . There exists a unitary matrix V that diagonalizes −iΛ̂N , i.e.
−iΛ̂N = V ∗DV, D = diag({−iσn}).
We perform another oordinate hange by,
1/2
z̃ = V ẑ = V AN žN , (4.63)
whih yields the deoupled equations,
dz̃
= Dz̃, (4.64)
dt
i.e.
dz̃n
= −iσnz̃n,
dt
=⇒ z̃ (t) = e−iσntn z̃n(0), n = −N, · · · , N. (4.65)
Performing the inverse Fourier operation, we see that,
N
∑
z̃(t, θ) = z̃n(t)e
inθ
n=−N
N
∑
=⇒ z̃(t, θ) = ei(nθ−σnt)z̃n(0). (4.66)
n=−N
94
This is an equation of superposition of 2N + 1 traveling waves with n being the
σn
wave number and vn = is the speed of propagation assoiated with n-th moden
of the wave.
4.5.3.2 Veloity of Propagation
We know that in the wave equation of (4.66), the veloity vn of propagation
σn
orresponding to n-th frequeny is determined as, vn = , where σn's are (real)n
1/2 1/2 1/2 1/2
eigenvalues of the Hermitian matrix Λ̂N = BN ΛNBN . Sine BN ΛNBN ∼
ΛNBN ∼ BNΛN (similar matries), the eigenvalues of Λ̂N an be haraterized by
those of BN and ΛN .
N 1 1
We have eig(BN ) = {βn}−N , where βn = = 2 . Now, ΛN is a Toeplitzαn 1+χn
2N
Hermitian matrix formed by the Fourier oeients {λn}−2N . Given those oe-
ients, it is in general not possible to write down losed form representation of its
eigenvalues. However, the bounds of eigenvalues of suh a matrix is well known.
We will make a little detour to state these results.
95
4.5.3.3 Toeplitz Matries and Eigenvalues
Let f be a periodi funtion over the interval [0, 2π) and {fn} are its Fourier
oeients. Let us denote Tn(f), the n×n Toeplitz Hermitian matrix dened as,
 
f0 f−1 f−2 · · · f −(n−1)
 
 
 f 
 1 f0 f−1 · · · f−(n−2)
 
 
T n(f) = f2 f1 f0 · · ·  . (4.67) 
 
 
. .
.
 
. . .
.
. .
 
 
 
f(n−1) f0
We also dene,
mf = ess inf(f), (4.68)
Mf = ess sup(f). (4.69)
Let the eigenvalues of Tn(f) be denoted by τn,k, k = 1, ..., n. Then,
mf ≤ τn,k ≤ Mf . (4.70)
Note that, max |τn,k| ≤ max(|mf | , |Mf |) ≤ M|f |. We reall another useful result
k
here.
Lemma 4.5.2. Let P,Q be Hermitian positive denite matries of same order.
If τ(X) denote eigenvalues of X, then,
τmax(PQ) ≤ τmax(P ) · τmax(Q) (4.71)
τmin(PQ) ≥ τmin(P ) · τmin(Q) (4.72)
96
4.5.3.4 Bounds on the Veloity
To get bounds on σn's, we need an useful assumption that will be apparent shortly.
A-1 mλ > 0, i.e. ΛN = T2N+1(λ) ≻ 0.
Sine µn = τ(BNT2N+1(λ)), Lemma 4.5.2 gives the following bound, for n =
−N, ..., N ,
τmin(BN)τmin(T2N+1(λ) ≤ σn ≤ τmax(BN)τmax(T2N+1(λ)),
(4.73)
⇒ mλ= ≤ σ ≤ M .
1 + χN2
n λ
Aordingly, the veloity is bounded by,
mλ ≤ ≤ Mλvn (4.74)
n(1 + χN2) n
Remark 4.5. The assumption A-1 an be extended to inlude the ase Mλ < 0,
i.e. ΛN ≺ 0 as well.
4.5.3.5 Speial Cases
1. Case - I: Constant λ :
λ(θ) = λ0 6= 0, a onstant. We may assume that λ0 > 0. In this ase, (4.59)
an be expliitly solved and the solution an be expressed as,
∞ ( )βnλ0
∑
in θ− t
z(t, θ) = e n z −jλ0tn(0) + e z0(0). (4.75)
n=−∞
n=6 0
λ0
Here the veloities vn an be written as, vn = , n ∈ Z \ {0}.n(1+χn2)
97
2. Case - II: Band-limited λ :
Here, we onsider only one frequeny omponent of λ, i.e. λ(θ) = λ0 +
λ eiθ1 + λ e
−iθ
−1 . Let A-1 hold, i.e. ΛN = T2N+1(λ) ≻ 0. In this ase, ΛN is
a tri-diagonal Toeplitz matrix. The eigenvalues of suh a matrix are known
to take the following form.
Lemma 4.5.3.
( )
kπ
τk(ΛN) = λ0 + 2 |λ1| cos , k = 1, ..., 2N + 1. (4.76)
2N + 2
This, ombined with Lemma 4.5.2, we get the following bound,
( )
λ0 + 2 |λ1| cos (2N+1)π2N+2 ( )≤ σ ≤ λ + 2 |λ | cos πn 0 1 (4.77)
1 + χN2 2N+2
4.5.4 Strong Coupling Limit, χ → ∞
It is interesting to note that in the limit χ → ∞, the equations (4.47) take simple
∞ 1
form. To see this, note that for some z ∈ C (S ;R), we may express
∞
∑
A−1z = βn 〈φn, z〉 φn,
n=−∞
1 −1
where βn = 2 are the eigenvalues of A , and φn = einθ. This implies,1+χn
( )
∑
A−1 1lim z = lim 〈z, φ 〉φ
χ→∞ χ→∞ n n1 + χn2
n
= 〈z, 1〉
∫ 2π
=⇒ lim A−1z = zdθ.
χ→∞ 0
98
t = 0 s t = 6.67 s t = 13.33 s t = 20 s
1 1 1 1
PSfrag replaements
0.5 0.5 0.5 0.5
0 0 0 0
-0.5 -0.5 -0.5 -0.5
-1 -1 -1 -1
0 5 0 5 0 5 0 5
θ θ θ θ
3 3 3 3
2 2 2 2
1 1 1 1
0 0 0 0
-1 -1 -1 -1
-2 -2 -2 -2
-3 -3 -3 -3
-2 0 2 -2 0 2 -2 0 2 -2 0 2
x1 x1 x1 x1
Figure 4.1: Numerial solution of (4.61) for experiment 1. Evolution of u1 (blue)
and u2 (red) is given in the rst row; the ontrols form a simple traveling wave
along the θ domain. The seond row shows evolution of x1, x2 variables.
With this realization, we may rewrite (4.47) in the strong oupling limit as,
∫
∂u 2π1
= λu2dθ
∂t 0
(4.78)
∫
∂u 2π2
= − λu1dθ.
∂t 0
It is lear from (4.52) that both u1 and u2 are independent of θ. Then the equations
(4.78) an be equivalently written as
u̇1 = λ̃u2
(4.79)
u̇2 = −λ̃u1,
99
u1, ux 22
x2
x2
x2
t = 0 s t = 4 s t = 8 s
1 0 0
0.5
-0.5 -1
0
-1 -2
-0.5
-1 -1.5 -3
1 1 1
1 1 1
0 0 0
PSfrag replaements 0 0 0
-1 -1 -1 -1 -1 -1
t = 12 s t = 16 s t = 20 s
0 0 0
-1 -2
-2
-2 -4
-4
-3 -6
-4 -6 -8
1 1 1
1 1 1
0 0 0 0 0 0
-1 -1 -1 -1 -1 -1
Figure 4.2: Evolution of x3 for experiment 1. The blue loop represents the irle
S1 while the height of eah point on the red loop is given by the value of x3 at
the orresponding point of θ.
∫ 2π
where λ̃ = λdθ. Equations (4.79) are optimal ontrol evolution equations for
0
a single agent [Justh and Krishnaprasad, 2016℄. Thus in the strong oupling limit
χ → ∞, the optimal ontrol solutions for the ontinuum of agents ollapses to
that of a single agent. This is alled the synhronization of the ok, where every
agent in the ok behaves the same way.
100
t = 0 s t = 6.67 s t = 13.33 s t = 20 s
1 1 1.5 1.5
PSfrag replaements
0.5 1 1
0.5
0.5
0 0.5
0 0
-0.5 0
-0.5
-0.5
-1 -0.5 -1
-1 -1.5 -1 -1.5
0 5 0 5 0 5 0 5
θ θ θ θ
3 3 33
2 2 22
1 1 11
0 0 00
-1 -1 -1-1
-2 -2 -2-2
-3 -3 -3-3
-2 0 2 -2 0 2 -2 0 2 -2 0 2
x1 x1 x1 x1
Figure 4.3: Numerial solution of (4.61) for experiment 2. Evolution of u1 (blue)
and u2 (red) is given in the rst row. The seond row shows evolution of x1, x2
variables.
4.5.5 Simulation Results
We simulate the evolution of optimal ontrols u1 and u2 governed by the linear
partial dierential equations (4.47) by means of Fourier analysis as presented in
the setion 4.5.3. In partiular, here we present the solutions of the trunated
ordinary dierential equations (4.61), where we only keep trak of rst N + 1
Fourier oeients of eah variable. Note that λ is assumed to have less than
N + 1 oeients. We will now present numerial solutions in dierent ases by
101
u1, ux 22
x2
x2
x2
t = 0 s t = 4 s t = 8 s
1 0 0
0.5 -1
-0.5
0 -2
-1
-0.5 -3
-1 -1.5 -4
1 1 1
1 1 1
0 0 0
PSfrag replaements 0 0 0
-1 -1 -1 -1 -1 -1
t = 12 s t = 16 s t = 20 s
0 0 0
-2
-2
-4 -5
-4
-6
-6 -8 -10
1 1 1
1 1 1
0 0 0 0 0 0
-1 -1 -1 -1 -1 -1
Figure 4.4: Evolution of x3 for experiment 2. The blue loop represents the irle
S1 while the height of eah point on the red loop is given by the value of x3 at
the orresponding point of θ.
varying the initial onditions and parameter values. For all the experiments pre-
sented here, N is taken to be 30 and the nal time T is set as 20 seonds and
four snapshots of the optimal ontrols are shown. The Hamiltonian is veried to
be staying a onstant (up to mahine preision) for all of the experiments. The
evolution of state variables is also reorded. For all the experiments presented
here, (x1, x2) is set to start from a unit irle and x3 is initially zero for all values
of θ.
102
PSfrag replaements
t = 0 s t = 6.67 s t = 13.33 s t = 20 s
4 1.5 1 1.5
1 0.5 1
3
0.5 0 0.5
2 0 -0.5 0
-0.5 -1 -0.5
1
-1 -1.5 -1
0 -1.5 -2 -1.5
0 5 0 5 0 5 0 5
θ θ θ θ
PSfrag replaements
t = 0 s t = 6.67 s t = 13.33 s t = 20 s
4 3 1.5 2
1
2 1.5
3 0.5
1
1 0
2 0.5
0 -0.5
0
1 -1
-1
-1.5 -0.5
0 -2 -2 -1
0 5 0 5 0 5 0 5
θ θ θ θ
PSfrag replaements
t = 0 s t = 6.67 s t = 13.33 s t = 20 s
4 4 3 4
3 3
3 2
2 2
2 1
1 1
1 0
0 0
0 -1 -1 -1
0 5 0 5 0 5 0 5
θ θ θ θ
Figure 4.5: Numerial solutions of (4.61) for experiment 3 for (a) rst row, χ =
0.1, (b) seond row, χ = 1, () third row, χ = 10. In all plots, u1 is shown in blue
and u2 is shown in red. The speed of information propagation dereases as the
oupling onstant inreases.
103
u1, u2 u1, u2 u1, u2
The simplest set of initial onditions u1(0, θ) = cos(θ), u2(0, θ) = sin(θ), λ(θ) =
5, χ = 1 generate a traveling wave solution, aording to equation (4.75). This
an be seen from Fig. 4.1. The states (x1, x2, x3) are also integrated from an
initial ondition for (x1, x2) forming a irle on the plane and x3 being identially
zero for all θ. The evolution of (x1, x2) an be seen from Fig. 4.1. The shape
of the (x1, x2) irle did not hange, although its size varied over time. The
evolution of x3 is shown in Fig. 4.2, whih appeared to derease steadily for all
θ. Next, for experiment 2, we onsider a band-limited λ, i.e. λ(θ) = 5 + cos(θ).
Keeping all other onditions same as in experiment 1, we get Fig. 4.3-4.4. Here
both size and shape of the (x1, x2) irle hanged over time. The value of x3
dereased in this ase as well but more asymmetrially than in experiment 1. In
experiment 3, we show how a loalized disturbane gets spread in the ontinuum.
For this experiment, we let the ontrol u1 is initially zero everywhere, u1(0, θ) = 0.
However, u2 has a loalized peak at a ertain spatial point. We took the example
of a Gaussian form,
1 − (θ−π)
2
u2(0, θ) = 2ρ2√ e , (4.80)
2πρ2
with ρ = 0.1. We then plot the solutions in three dierent settings of χ values,
χ = 0.1, 1, 10 in Fig. 4.5. λ is taken to be a onstant in all these ases, λ(θ) = 5.
It is disovered in the previous setion that the speed of traveling wave dereases
as χ inreases. This an be seen learly from Fig. 4.5 as the disturbane is seen
to be not well propagated for higher values of χ.
104
4.6 A Continuum of Agents on the Plane
As an extension of the problem onsidered in Setion 3.3, we will explore the ase
where we take the underlying group as the speial Eulidean group, SE(2). This
ase then an be seen as a ontinuum ounterpart of [Justh and Krishnaprasad,
2015b℄. Every g(t, θ) ∈ SE(2) an be represented as,
 
cosx3 − sin x 3 x1
 
 
g =  sin x3 cosx x
 ,
 3 2
 
 
0 0 1
with group evolution dynamis,
∂
g(t, θ) = g(t, θ) · (u1(t, θ)A1 + u2(t, θ)A2) , (4.81)
∂t
where,
   
0 −1 0 0 0 1
   
   
   
A1 = 1 0 0 , A =  2 0 0 0 ,
   
   
   
0 0 0 0 0 0
along with A3 = [A1, A2], form a basis for the assoiated Lie algebra se(2). Note
that sine {A1, A2} is braket generating in se(2), similar argument as in the
Heisenberg ase would provide (weak) ontrollability result in this ase as well.
We then seek to minimize the same ost funtional as in (4.38),
( ( ))
T 2π ( )2 ( )∫ ∫ 21
J = (u2 + u2
∂u1 ∂u2
) + χ + dθdt
2 1 20 0 ∂θ ∂θ
105
4.6.1 Equations of Optimal Control
4.6.1.1 Calulus of Variations Approah
The system dynamis an be equivalently expressed as,
∂x1
= u2 cosx3
∂t
∂x2
= u (4.82)2 sin x3
∂t
∂x3
= u1.
∂t
Note that the ontrols u1 and u2 an be expressed by the xj variables and their
derivatives as,
u1 = x3,t, u2 = x1,t cos x3 + x2,t sin x3.
Sine we have two integral onstraints in this ase, namely,
∫ T ∫ T
∆x1(θ) = u2 cosx3dt = (x1,t cosx3 + x2,t sin x3) · cosx3dt
0 0
∫ T ∫ T
∆x2(θ) = u2 sin x3dt = (x1,t cosx3 + x2,t sin x3) · sin x3dt,
0 0
∞ 1
we invoke Lagrange multipliers λ, µ ∈ C (S ;R) and the augmented Lagrangian
density an be read as,
L 1
[( ) ( )]
= x2 2 23,t + (x1,t cos x3 + x2,t sin x3) + χ x3,tθ + (x1,tθ cos x3 + x
2
2,tθ sin x3)
2
+ λ (x1,t cosx3 + x2,t sin x3) · cosx3 + µ (x1,t cosx3 + x2,t sin x3) · sin x3.
(4.83)
We reall the Euler-Lagrange equations from (4.22),
L ( L ) ( ) ( )∂ ∂ ∂ 2− − ∂ ∂L ∂ ∂L+ = 0,
∂xi ∂θ ∂xi,θ ∂t ∂xi,t ∂θ∂t ∂xi,θt
106
for i = 1, 2, 3. We note the following quantities,
u2 = x1,t cosx3 + x2,t sin x3,
u2,θ = x1,tθ cosx3 + x2,tθ sin x3,
(
∂ ∂L )
= u2,t cos x3 − u1 (u2 sin x3 + λ sin 2x3 − µ cos 2x3) ,
∂t ∂x1,t
(
∂2 ∂L )
= χ (u2,θθt cosx3 − x3,θ (u2,θt sin x3 + u2,θu1 cos x3)
∂θ∂t ∂x1,θt
− sin x3 (u2,θθu1 + u2,θu1,θ)) ,
( )
∂ ∂L
= u2,t sin x3 + u1 (u2 cosx3 + λ cos 2x3 + µ sin 2x3) ,
∂t ∂x2,t
( )
∂2 ∂L
= χ (u2,θθt sin x3 + x3,θ (u2,θt cosx3 − u2,θu1 sin x3)
∂θ∂t ∂x2,θt
+cosx3 (u2,θθu1 + u2,θu1,θ)) ,
( )
∂ ∂L
= u1,t,
∂t ∂x3,t
( L )∂2 ∂
= χu1,θθt.
∂θ∂t ∂x3,θt
Subsequently, the Euler-Lagrange equations for the Lagrangian (4.83) take the
form,
u2,t cos x3 − u1u2 sin x3 − u1 (λ sin 2x3 − µ cos 2x3)
= χ (u2,θθt cosx3 − sin x3 (u2,θθu1 + u2,θu1,θ)− x3,θ (u2,θt sin x3 + u2,θu1 cosx3)) ,
u2,t sin x3 + u1u2 cos x3 + u1 (λ cos 2x3 + µ sin 2x3)
= χ (u2,θθt sin x3 + cosx3 (u2,θθu1 + u2,θu1,θ) + x3,θ (u2,θt cosx3 − u2,θu1 sin x3)) ,
u1,t − χu1,θθt = −u2 (λ sin x3 − µ cosx3 − χu2,θx3,θ) ,
(4.84)
107
whih after some readjustments yield,
∂
(1− χ∆)u1 = −u2 (λ sin x3 − µ cosx3 − χx3,θu2,θ)
∂t
∂
(1− χ∆)u2 = u1 (λ sin x3 − µ cosx3 − χx3,θu2,θ) (4.85)
∂t
∂
(λ sin x3 − µ cosx3 − χu2,θx3,θ) = −u1(1− χ∆)u2.
∂t
Denoting A := (1 − χ∆) and µ3 := λ sin x3 − µ cosx3 − χu2,θx3,θ, we then an
express (4.85) as,
∂ Au1 = −µ3u2
∂t
∂ Au (4.86)1 = µ3u1
∂t
∂µ3
= −u1Au2.
∂t
4.6.1.2 PMP Approah
∞ 1 3
We introdue the ostate variable p(t) = (p1(t), p2(t), p3(t)) ∈ C (S ;R ), t ∈
[0, T ]. The pre-Hamiltonian (onsidering only normal extremals) an be written
as,
H = 〈ẋ, p〉 − L
2π ( )∫ 1 ( ) χ ( )
= u2 (p
2 2 2 2
1 cosx3 + p2 sin x3) + u1p3 − u + u − u
2 1 2 2 1,θ
+ u2,θ dθ.
0
(4.87)
Remark 4.6. Here we only onsider normal extremals, i.e. when p0 6= 0 and an
be normalized to −1. Note that in this ase, the emptiness (i.e. the full ostate
being identially zero) of the PMP would not our. It is of future eort to
investigate whether the nite odimensionality ondition is satised in the SE(2)
ase.
108
The maximum priniple (Theorem 4.4.1) would require us to maximize (4.87)
pointwise over the ontrols, i.e. we are attempting to nd the Hamiltonian as,
H = sup H
v1,v2∈C∞(S1;R)
∫ 2π
= sup (v2 (p1 cosx3 + p2 sin x3) + v1p3
v1,v ∈C∞2 (S1;R) 0
)
1 ( ) χ ( )− v21 + v22 − v21,θ + v22 2 2,θ dθ (4.88)
This maximization results in two Euler-Lagrange equations,
( )
∂H − ∂ ∂H = 0, i = 1, 2, (4.89)
∂vi ∂θ ∂vi,θ
that yields the optimal ontrols,
Au1 = p3,
(4.90)
Au2 = p1 cosx3 + p2 sin x3,
where we denote A = (1− χ∆), as usual. The Hamiltonian an be read as,
∫
H 1
2π
= (u1Au1 + u2Au2) dθ
2 0
∫
1 2π ( )
= p3A−1p3 + (p1 cosx3 + p2 sin x −13)A (p1 cosx3 + p2 sin x3) dθ
2 0
(4.91)
The dynamis of the ostate variable p an be alulated from Hamilton's equa-
∂p = − δH δHtion, , where denotes the funtional derivative of H with respet to
∂t δx δx
x. Expliit alulations yield,
∂p1
= 0
∂t
∂p2
= 0 (4.92)
∂t
∂p3
= −u2 (−p1 sin x3 + p2 cos x3) .
∂t
109
If we denote µ3 := −p1 sin x3 + p2 cosx3, then from (4.90) and (4.92), we an
redisover (4.86). Renaming µi = Aui, i = 1, 2, (4.86) an also be written as,
∂µ1
= −µ A−13 µ2
∂t
∂µ2
= µ A−1µ (4.93)3 1
∂t
∂µ3
= −µ A−12 µ1.
∂t
4.6.2 Strong Coupling Limit, χ → ∞
Similar to the Heisenberg ase, it an be shown that synhronization is ahieved in
∞ 1
strong oupling limit in this ase as well. We know that for some z ∈ C (S ;R),
− ∫1 2π
we have lim A z = zdθ. Thus, we may rewrite (4.93) in the strong oupling
χ→∞ 0
limit as,
( )
∫
∂µ 2π1
= −µ3 µ2dθ
∂t 0
( )
∫
∂µ 2π2
= µ3 µ1dθ (4.94)
∂t 0
( )
∫
∂µ 2π3
= −µ2 µ1dθ .
∂t 0
If we dene the variables,
∫ 2π
αi = µidθ, i = 1, 2, 3, (4.95)
0
then these variables will evolve aording to,
α̇1 = −α2α3
α̇2 = α α (4.96)1 3
α̇3 = −α1α2.
110
Equations (4.96) are the equations for a single agent senario and they are studied
in detail in Setion 3.2. This indiates the synhronization phenomenon in the
planar ontinuum ok. Note that α1 and α2 are essentially the optimal ontrols
u1 and u2, respetively. The Hamiltonian (4.88) simply beomes,
1 ( )
h∞ = α
2
1 + α
2
2 . (4.97)2
4.6.3 Simulation Results
While it has not been possible to haraterize general solutions of (4.93) analyti-
ally, here we demonstrate numerial solutions. To numerially solve the evolution
equations of µi variables, we used a nite dierene method. We partitioning the
spae domain [0, 2π] uniformly in M points, 0 = θ1, . . . , θM = 2π, so that the
2π
dierene between two onseutive spae points beome δθ = . In this dis-
M
rete setting any z(t, θ) an be approximated as an M vetor, z(t, θ) ≈ z(t) =
[z1(t), z2(t), . . . , zM(t)]
T
with the onstraint z1(t) = zM(t) for all t to respet the
periodiity property. Note also that in a seond-order entral dierene sheme,
the double partial spae derivative is expressed as,
∂2z(t, θj) zj+1 − 2zj + zj−1
= ,
∂θ2 δθ2
for j = 1, . . . ,M with appropriate adjustments for the boundary points j = 1,M .
The linear operatorA = (1−χ∆) an then be expressed by theM×M nonsingular
111
PSfrag replaements
t = 0 s t = 6.67 s t = 13.34 s t = 20 s
1 10 10 10
0.5 5 5 5
0 0 0 0 µ1µ2
µ3
-0.5 -5 -5 -5
-1 -10 -10 -10
0 5 0 5 0 5 0 5
θ θ θ θ
1.5 1 0.4 0.2
0.2 0
1 0.5
u1
0 -0.2 u2
0.5 0
-0.2 -0.4
0 -0.5 -0.4 -0.6
0 5 0 5 0 5 0 5
θ θ θ θ
0.01 0.5 1 2
0 0 0 0
-0.01 -0.5 -1 -2
0 10 5.5 6 6.5 6 7 8 4 6 8
x1 10 -3 x1 x1 x× 1
Figure 4.6: Numerial solution of (4.99) (rst two rows) along with the state
evolution (third row) for experiment 1. The arrows represent the diretion of
movement of the partile at that point.
matrix AM ,
 
( ) ( ) ( )
1 + 2χ − χ 0 . . . − χ 0
 δθ2 δθ2 δθ2 
 
 
( ) ( ) ( )
 − χ 1 + 2χ − χ 0 . . . 0 
 δθ2 δθ2 δθ2 
 
 
( ) ( ) ( )
AM =  0 − χ 1 + 2χ χ  . (4.98)
 δθ2 δθ2
− 2 . . . 0δθ 
 
 
. . .
.
 
. . . .
.
. . .
 
 
 
( ) ( ) ( )
0 − χ2 . . . 0 − χ2 1 + 2χδθ δθ δθ2
112
x2
u µ , µ , µ1, u2 1 2 3
With this notation, the partial dierential equations (4.93) an be expressed as a
system of ordinary dierential equations (ODE),
µ̇ = −µ A−11 3 M µ2
µ̇ = µ A−1µ (4.99)2 3 M 1
µ̇3 = −µ −12AM µ1.
These ODEs (4.99) are then solved using a mid-point based ODE solver in MAT-
LAB. The optimal ontrols ui's an be derived from the µ variables by the relation,
ui = A
−1
M µi, i = 1, 2 whih are used in the quadrature of the state variables xj 's.
Here we present results of some experiments with varying initial onditions. The
nal time T and spae disretization fator M is kept xed at T = 20 seonds
and M = 128 for all the experiments. A high value of M is hosen for a faithful
alulation of the spatial derivatives. In the subsequent experiments we try to
investigate the behavior of a simple loop under the optimal ontrols generated
by (4.99), i.e. we take, x1(0, θ) = 0.01 cos(θ) and x2(0, θ) = 0.01 sin(θ) so that
initially the partiles start on a irle. The remaining initial onditions of x3 and
µi variables and the parameter χ is varied in the following experiments. It is to be
noted that the Hamiltonian and the Casimir variables are validated to be onstant
in eah of the experiments.
Experiment 1
We take a simple example where eah agent start moving in the positive x axis
and with unit speeds, i.e. x3(0, θ) = 0, µ2(0, θ) = 1. The initial urvature of eah
113
th th th th
0.02
10 10
0.01 10
0 0 0 0
PSfrag replaements
-0.01
-10
-0.02 -10 -10
0 -10 0 10 -10 0 10 -10 0 10
x1 x1 x1 x1
Figure 4.7: State evolution is only shown for experiment 2. The arrows represent
the diretion of movement of the partile at that point.
agent is taken to be zero, µ1(0, θ) = 0 whih means every agent starts harmo-
niously with same veloity and urvature. The value of χ is taken as 1. Four
snapshots of the µ, u and x variables are shown in Fig. 4.6. The urvature eld is
seen to be forming two peaks in the spatial domain whih gives rise to the twisted
form of the initial irle.
Experiment 2
We keep all the initial onditions same as in experiment 1 exept the initial dire-
tion of movement of the partiles. It is simulated that almost half the partiles
try to go in one diretion while the other half in the opposite diretion. To write
this as a ontinuous periodi funtion, we take


( ( ( )))
π

 1 + tanh 100 θ − π , if 0 ≤ θ < π
2 2
x3(0, θ) = . (4.100)

 ( ( ( )))
π
 1− tanh 100 θ − 3π , if π ≤ θ ≤ 2π
2 2
This denition of the initial diretion means that the partiles on the `east' part
of the irle initially go to the right and the partiles on the `west' part go to
the left. The simulation results are shown in Fig. 4.7. Comparing to Fig. 4.6,
114
x2
th th th th
0.02 1010 10
0.01 5
0 0 0 0
PSfrag replaements
-0.01 -5
-0.02 -10
-10
-10
0 -10 0 10 -10 0 10 -5 0 5
x1 x1 x1 x1
Figure 4.8: State evolution is only shown for experiment 3. The arrows represent
the diretion of movement of the partile at that point.
the rst two rows are idential sine the initial onditions of µ variables did not
hange and hene they are omitted. What is interesting is that the irle splits
into two loops onneted by very small number of partiles.
Experiment 3
Similar to experiment 2, we try to investigate the eet of hange in initial dire-
tion of movement of the partiles. Here, we set the partiles to go on a radially
outward path, i.e. x3(0, θ) = θ, with keeping all other onditions same as in ex-
periment 2. The simulation results are shown in Fig. 4.8. The rst two rows are
not shown sine they are idential with experiment 1.
Experiment 4
In this experiment, again we x all the initial onditions and parameters same as
in experiment 1, exept the initial urvature is given a loal intensity. In other
words we hoose this Gaussian funtion,
1 (θ 2−π)
µ1(0, θ) = √ e− 2σ2 , (4.101)
2πσ2
115
x2
PSfrag replaements
t = 0 s t = 6.67 s t = 13.34 s t = 20 s
10 10 20 10
5
5 5 10
0 µ1µ2
0 0 0 µ3
-5
-5 -5 -10 -10
0 5 0 5 0 5 0 5
θ θ θ θ
1.1 1.5 1.5 1.5
1 1 1
1.05
u1
0.5 0.5 0.5 u2
1
0 0 0
0.95 -0.5 -0.5 -0.5
0 5 0 5 0 5 0 5
θ θ θ θ
2
1
0.01 1 1
0 0.5
0 0
-0.01
0 -1
-1
0 10 3 3.5 4 5 6 4 4.5 5 5.5
x1 -3 x x×10 1 1 x1
Figure 4.9: Numerial solution of (4.99) (rst two rows) along with the state
evolution (third row) for experiment 4. The arrows represent the diretion of
movement of the partile at that point.
with σ = 0.05. The purpose of hoosing this initial ondition is to see whether
a loalized information gets spread aross the ontinuum or not. The simulation
results are shown in Fig. 4.9.
Experiment 5
Finally, we demonstrate the eet of strong oupling limit by taking a large value
116
x2 u1, u2
µ1, µ2, µ3
PSfrag replaements
t = 0 s t = 6.67 s t = 13.34 s t = 20 s
2 10 20 20
5 10 10
1
0 0 0 µ1µ2
0 µ3
-5 -10 -10
-1 -10 -20 -20
0 5 0 5 0 5 0 5
θ θ θ θ
1.5 1 0.75 0.5
1
0
u1
0.5 0.5 0.7 u2
-0.5
0
-0.5 0 0.65 -1
0 5 0 5 0 5 0 5
θ θ θ θ
0.02 1.67 3.9
0.01 4.491.66 3.89
4.48
0 1.65
3.88 4.47
-0.01 1.64
3.87 4.46
0 10 6.46 6.47 6.48 12.5 12.515 17.4 17.42
x1 10 -3 x1 x x× 1 1
Figure 4.10: Numerial solution of (4.99) (rst two rows) along with the state
evolution (third row) for experiment 5. The arrows represent the diretion of
movement of the partile at that point. The u1, u2 solutions are almost `at',
indiating single agent solution or synhronization.
of χ = 1000. We note that even in the ase, µ1(0, θ) = cos(θ), µ2(0, θ) = 1 +
0.2 cos(θ), µ3(0, θ) = sin(θ), x3(0, θ) = π/4, we essentially get the system derived
by the optimal ontrols that are spatially non varying.
117
x2 u1, u2
µ1, µ2, µ3
4.7 Disussion and Sope of Future Researh
In this hapter, we have presented a general framework for a ontinuum desrip-
tion of a ok. We are interested in solving optimal ontrol problems to explain
olletive movement of suh a ok. We reognize this is a hallenging problem
that naturally provides several open questions for further researh. We itemize
few suh possibilities.
• It is shown that under a ertain nite odimensionality of a reahable set,
the PMP remains valid in a general Hilbert spae setting. One might want to
disover its relationship with the ontrollability ondition. In partiular, if a
system on the loop group is strongly ontrollable, does the PMP ondition
satisfy automatially? We have been able to show this to be true in the
Heisenberg ase. Does this remain valid if we only have weak ontrollability?
• It is of interest to extrat meaningful features of the optimal solutions of
the SE(2) ase. While we have not been able to solve (4.93) analytially,
we want to answer few questions about it. For example, do these equations
possess a traveling wave like solution (like the Heisenberg ase)? If so, what
is the speed of those waves? The answer might give an insight toward the
information transfer in biologial swarms. It is the inherent nature of the
numerial study presented here that there exist many possibilities by varying
initial ongurations whih makes it partiularly diult. We explored only
a tiny fration of possible variations in the initial ongurations in this
118
doument. A future work ould perform more extensive numerial study of
these partial dierential equations (4.93).
• We have presented the results in this hapter under the ase of a xed yli
interation topology. A more general, possibly state dependent (hene time
dependent, too) interation sheme an be modeled and subjeted under
similar questions.
• It will be an interesting future work to establish ontinuum parallel of the
Lie-Poisson redution in the loop group ase.
119
Chapter 5
Cognitive Cost of Floking: A Geometri
and Hamiltonian Perspetive
5.1 Introdution
It has been an appealing question to researhers from several elds to address
how natural olletives funtion at a fundamental level. Many theories have been
proposed to desribe this phenomenon over the past few deades. The lak of
aurate motion apturing tehnologies had limited the study of natural olletives
for many years. However, as motion apture beame more sophistiated, more
movement data of these olletives were reorded. This enabled researhers to
unover several underlying mehanisms behind oking [Ballerini et al., 2008a,b;
Attanasi et al., 2014; Cavagna et al., 2018; Nagy et al., 2010℄. These studies
shed light on how individual agents interat with its neighboring agents or how
information may be propagated through the ok.
120
Continuing the spirit of the `top-down' view of the ok, we will present a novel
perspetive for analyzing olletive motion data. The ok movement results in a
time-series of its kineti energy, whih an be divided into several energy modes.
Normalized modes dene a urve in some appropriate dimensional simplex whih
we attempt to desribe by an evolutionary game dynamis. Individual energy
modes are onsidered as pure strategies of suh a game. An optimal ontrol
problem is proposed to best t the data on the simplex, where the ontrol inputs
modulate the tness assoiated with the strategies. This is in ontrast to the
optimal ontrol problem posed in Chapter 4, where the ontrols are `low-level'
i.e. individual ontrol inputs are determined post-solution of the optimal ontrol
problem. In the present ontext however, the ontrol inputs are `high-level'. The
olletive itself is thought to be deiding the optimal alloation of its energy
among several dierent modes during a ight event. A notion of ognitive ost is
introdued to denote the optimal ost for the olletive to perform this alloation.
This work brings together several key ingredients for this data-driven approah.
In setion 5.2, the motion data of European oks is detailed. This data is then
subjeted to a linear data smoothing tehnique [Dey and Krishnaprasad, 2012℄
that reonstruts smooth trajetory data of eah bird in the ok. A nonlinear
data smoothing tehnique [Dey and Krishnaprasad, 2014℄ is later used for the
optimal energy alloation problems. These smoothing tehniques are based on
optimal ontrol theory and are desribed in setion 5.3. Setion 5.4 ontains the
geometri theory developed in [Mishiati and Krishnaprasad, 2017℄ to ompute
121
dierent energy modes. Finally, in setion 5.5, a generative model based on [Raju
and Krishnaprasad, 2018℄ is desribed to onstrut an optimal ontrol problem on
a simplex. Numerial solution of this optimal ontrol problem, as well as the idea
of ognitive ost, are presented in setion 5.6. This is a joint work with V. Raju
[Halder et al., 2019b℄.
5.2 Floking Data
We are provided with ight data of European starlings that were taken by Dr.
Andrea Cavagna and his ollaborators from the Colletive Behaviour in Biologial
Systems (COBBS) group at the Institute for Complex Systems (ISC-CNR), Uni-
versity of Rome La Sapienza". Starlings gather around urban areas during the
winter months in order to get extra warmth from the ities. Floks of these kind of
birds are well known to perform remarkable maneuvers, the purpose and meha-
nisms of whih still elude researhers. Equipped with modern imaging tehniques
and sophistiated algorithms for stereo reonstrution, these group of researhers
managed to apture a series of ight events with dierent ok sizes in the winter
months of 2011. See [Attanasi et al., 2014℄ for more details about the proess. We
will study eight partiular oking events, the details of whih are given in table
5.1.
122
Floking Flok Size Duration Data Capture Rate
Event (n) (seonds) (frames/seond)
1 175 5.4875 80
2 123 1.8176 170
3 46 5.6118 170
4 485 2.3471 170
5 104 3.8824 170
6 122 4.1588 170
7 380 5.7353 170
8 194 1.7588 170
Table 5.1: Details of aptured oking events
5.3 Data Smoothing
Given a time-indexed sequene of sampled observations on a manifold, genera-
tive models provide a meaningful way of apturing them through the use of an
underlying dynamial system omplete with ontrol inputs having useful interpre-
tations. The ontrol inputs are determined by solving an optimal ontrol problem,
where the ost funtion onsists of a tness term that penalizes mismath between
the generated trajetory and sampled data, and a smoothing term weighted by a
parameter λ that aets the smoothness of the generated trajetory. We disuss
two generative models that have been proposed to solve this problem.
123
5.3.1 A linear generative model
A rst approah to solving the data smoothing problem, presented in [Dey and
Krishnaprasad, 2012℄, is to formulate an optimal ontrol problem to minimize the
jerk path integral, with intermediary state osts determining the t error. Suppose
N 3
that {ri}i=0 denote the positions of the birds at eah sampling time, with ri ∈ R .
3
In order to reover a trajetory t r(t) : [t0, tN ] → R , one an use the jerk-driven
linear generative model,
ṙ(t) = v(t)
v̇(t) = a(t)
ȧ(t) = u(t) (5.1)
where v(t), a(t),u(t) denote the veloity, aeleration and jerk (input) of the tra-
jetory. The ost funtional to be minimized is
N tN∫
∑
J = || 2 2l r(ti)− r(t)|| + λ ||u(t)|| dt (5.2)
i=0 t0
where the minimization is over initial onditions r(t0),v(t0), a(t0) and the input
u(t). Dening the state and output as
 
r(t)
 
 
 
x(t) = 9 3 (t)  ∈ R ,y(t) = x(t) ∈ Rv
 
 
 
a(t)
124
we obtain the linear state equations
ẋ(t) = Ax(t) +Bu(t)
y(t) = Cx(t), where
   
0 13 0 0   
   
   
A =  0 0 1 , B =  0  , C = [13 0 0] (5.3)
 3   
   
   
0 0 0 13
Therefore, the problem of minimizing Jl subjet to (5.3) is a linear, quadrati
optimal ontrol problem, whih an be solved by a ompletion of squares of terms
in the ost by invoking a path independene lemma, or by applying the Pontryagin
Maximum Priniple as shown in [Dey and Krishnaprasad, 2012℄. This approah
has been used to smooth the starling ok data for all the events listed in table
5.1.
5.3.2 Data smoothing in the Eulidean setting
In this setion, we present a general result on the Pontyagin Maximum Priniple
n
based approah for data smoothing on the Eulidean spae R . Suppose that
{ }
xd
N
i denote the sampled data. For a generative model given by the dynamisi=0
ẋ = f(x, u) n mon R , with the ontrol u ∈ R , the optimal ontrol problem an be
formulated as:
∫ N
λ tN ∑
min J(x(t0), u) = ‖u‖2 dt+ Fi(x(ti), xdi ),
x(t0), u∈Rm 2 t0 i=0
(5.4)
subjet to: ẋ = f(x, u),
125
where parameter λ > 0 is a regularization parameter, and Fi's are suitably dened
t errors of the reonstruted trajetories and sampled data at the sampling
times. Using Pontryagin's Maximum Priniple, the optimal ontrol values an be
alulated as a funtion of the state and a o-state variable. The following result
from [Dey and Krishnaprasad, 2014℄ states this preisely.
Theorem 5.3.1. (PMP for data smoothing [Dey and Krishnaprasad, 2014℄ ) Let
u∗(·) ∗be an optimal ontrol input for (5.28), and let x (·) denote the orresponding
n
state trajetory. Then there exist a ostate trajetory p : [t0, tN ] → R , p 6= 0 ,
suh that
∗ ∂Hẋ = (t, x∗, p, u∗)
∂p
(5.5)
∂H
ṗ = − (t, x∗, p, u∗)
∂x
during t ∈ (ti, ti+1), i = 0, 1, ..., N − 1, and the Hamiltonian is given as
H(t, x∗, p, u∗) = max H(t, x∗, p, v), (5.6)
v∈Rm
for t ∈ [t0, tN ]\{t0, t1, ..., tN}, where the pre-Hamiltonian is dened asH(t, x, p, u) =
〈p, f(x, u)〉 − λ ‖u‖2. Moreover, jump disontinuities of the ostate variable an
2
be written as
p(t−0 ) = 0,
p(t+)− p(t− ∂Fi(x(ti))i i ) = , i = 0, 1, ..., N, (5.7)∂x(ti)
p(t+N) = 0.
The pieewise ontinuous nature of the o-state trajetory due to jump on-
ditions arising from mismath between the sampled data points and the reon-
126
struted state must be noted here. The initial ondition x(t0) is identied by using
the terminal ondition for the o-state, while the optimal value of λ is typially
obtained through leave-one-out or ordinary ross validation. The reonstruted
trajetory is then obtained as the projetion onto the state spae of the solution
of Hamilton's equations derived from the (maximized pre-) Hamiltonian. We refer
the reader to [Dey, 2015℄ for a detailed treatment of these problems. This is the
result that we will use in our data tting problem on a simplex.
5.4 Energy Modes
Avian oks display a variety of ight behaviors that may be haraterized as
olletive strategies suh as steady direted translation of enter of mass (whih
we denote by om), oherent rotation about enter of mass (rot), hange of form
(ens), internal re-shuing of relative positions (dem), rapid expansion or on-
tration of volume (vol) et. A oking event may display all of the mentioned
strategies to varying degrees as governed by the time-dependent alloation of ki-
neti energy to eah strategy. We take the viewpoint presented in [Mishiati and
Krishnaprasad, 2017℄ and study the frations of the total kineti energy of a ok
alloated to several `kinemati modes'  rigid translations, rigid rotations, inertia
tensor transformations, expansion and ompression, in order to desribe olletive
behavior.
If the positions of the birds in a ok are denoted by {r1, r2, ..., rn}, the enter
127
of mass an be written as,
n
1 ∑
r = ri, (5.8)om
n
i=1
where we treat every bird alike, i.e. their masses are taken to be equal. The
ensemble inertia tensor is dened by
n
∑
K = (ri − Tr ) (r − r ) . (5.9)om i om
i=1
Let the veloities of the birds be denoted as, {vr1, ...,vrn}, then the total kineti
energy is,
n
1 ∑
E = ‖ 2vri‖ . (5.10)
2
i=1
We an dene the position and veloity vetor with respet to the enter of mass,
i.e. c , [c1, ..., c ] ∈ R3×n 3×nn , where ci = ri − r ; vc , [vc1,vc2, ...,v ] ∈ R ,om cn
where vci = vri − v . Then,om
n
n 1 ∑
E = ‖ 2v ‖ , E ‖ 2, vci‖ . (5.11)om om rel
2 2
i=1
We thus have the splitting, E = E + E . As presented in [Mishiati and
om rel
Krishnaprasad, 2017℄, instantaneous relative energy alloations an be expressed
4 1
on a probability simplex (∆ ) by exploiting the ber bundle strutures of the
ok's total onguration spae to split the total kineti energy using (i) ensemble
bration or (ii) shape bration.
(i) Ensemble Fibration: We note that the ensemble inertia tensor K (5.9) is
a symmetri positive denite matrix. Hene its eigendeomposition an be
1 n
Note that in this hapter we will use ∆ to denote the n-dimensional simplex.
128
written as, K = QΛQT, with Λ = diag(λ1, λ2, λ3), where λ1 ≥ λ2 ≥ λ3 > 0.
Dene, F := cvT+v cTc and F̃ = [F̃
T
ij ] = Q FQ. Then the following energyc
modes an be alulated,
( )
1 F̃ 2 F̃ 2 F̃ 2
E , 12 + 13 + 23
ens.rot
2 λ1 + λ2 λ1 + λ3 λ2 + λ3
( )
(5.12)
1 F̃ 211 F̃
2 2
E , + 22
F̃
+ 33 .
ens.def
8 λ1 λ2 λ3
Furthermore,
( )
1 2 Ttr cv
E c, , (5.13)
vol
2 tr(K)
so that, E = E − E . We may also alulate E = E −
ens.res ens.def vol dem rel
E − E . Hene, in this bration we have the following splitting of
ens.rot ens.def
the kineti energy,
{ }
E E E E E
om dem ens.rot vol ens.res
, , , , ∈ ∆4 (5.14)
E E E E E
(i) Shape Fibration: Dene
n
∑
J = (ci × vci) ,
i=1
(5.15)
n
∑
( )
Ic = ‖ci‖2 T1− cici .
i=1
Then the rotational energy E an then be alulated as,
rot
1
E , JTI−1J, (5.16)
rot
2 c
The shape residual energy is given by E = E −E −E , whih
shp.res rel rot end.def
provides the splitting in this bration as below
{ }
E E E E E
om rot shp.res vol ens.res
, , , , ∈ ∆4 (5.17)
E E E E E
129
While we an split the kineti energy in 5 dierent modes (5.14),(5.17), many
oking events show a predominant alloation of nearly onstant energy of rigid
translation (E ). We exlude this omponent from the total E in our analysis,
om
and onsider the alloation of the remaining energy E to obtain a time dependent
rel
trae of eah event on a lower dimensional simplex. In partiular, we apture the
trae generated by the following deomposition of E using ensemble bration
rel
on the 1-simplex by two dierent methods,
{ }
E E
dem ens
( 1ENS-I) , ∈ ∆ , (5.18)
E E
rel rel
{ }
E E − E
ens.rot rel ens.rot
(ENS-II) , ∈ ∆1, (5.19)
E E
rel rel
where E = E−E , and E = E −E = E +E +E . Similarly,
rel om ens rel dem ens.rot vol ens.res
a one dimensional simplex desription using shape bration may be given by two
ways,
{
E E − }E
shp.res rel shp.res
( 1SHP-I) , ∈ ∆ , (5.20)
E E
rel rel
{ }
E E
rot shp
(SHP-II) , ∈ ∆1, (5.21)
E E
rel rel
where E = E − E = E + E + E .
shp rel rot shp.res vol ens.res
In this way, moment-to-moment deisions made by individuals in a ok, tak-
ing aount of the deisions of their neighbors, ontribute to ok-sale strategies
as aptured by suh time dependent traes on the probability simplex. Treating
the strategy prevalene as being given by the respetive energy frations, we resort
to a generative evolutionary game dynamis to model the ompetition between
the ok-sale strategies.
130
5.5 Generative model on the 1-simplex and the
data-smoothing problem
Sine we are interested in desribing the evolution of two ok strategies as in
eqs. (5.18) and (5.19) for ensemble bration or eqs. (5.20) and (5.21) for shape
bration, we apture the trae of oking events via a generative model on the
1-simplex. We onsider an evolutionary game model, namely repliator dynamis
equipped with a multipliative ontrol, in order to desribe their evolution in the
interior (0, 1) of the one-dimensional simplex. The hoie of repliator dynamis
is inuened by its universality in desribing simplex-preserving dynamis, and
by virtue of being an extremal for a variational problem [Svirezhev, 1972; Raju
and Krishnaprasad, 2018℄. Presently, with the inlusion of a ontrol variable,
we onsider a dierent variational problem that aims to perform data smoothing
using regularization as in [Dey and Krishnaprasad, 2014℄. To see this, let x =
[x x ]T ∈ ∆11 2 where xi, i = 1, 2 denote the prevalene of strategies i (to be
speied) on the simplex with the natural onstraint x1 + x2 = 1. xi = 1, i = 1, 2
orrespond to alloation of Erel entirely to one of the two pure strategies. Suppose
that the frequenies assoiated with the strategies are updated aording to the
rule
f i(x)
xi(t+ 1) = xi(t) (5.22)
f̄(x)
131
i
where the tness f (x) = Ax and f̄ = x1f
1(x) + x f 22 (x)
2
. Here, A = [aij ] ∈ R
th th
denes a payo matrix with aij denoting the payo of the i strategy against j
strategy. In the ase that the payos do not depend on the strategy j of against
whih it is mathed up, the olumns of A are idential. In the ode limit of (5.22),
after an inhomogeneous time-sale hange, we get the mean eld equations:
ẋi(t) = xi(t)(f
i(x)− f̄(x)), i = 1, 2 (5.23)
It an be readily veried that (5.23) is simplex-preserving, leaving the pure strate-
gies invariant. Sine addition of the same term to eah omponent of the tness
keeps the dynamis (5.23) unhanged, by subtrating a21 and a12 from the rst
and seond olumn elements of A respetively, we get the equivalent payo matrix
 
a
 11 − a21 0 
Ã =   (5.24)
 
0 a22 − a12
We introdue a ontrol input ũ that sales the tness, and hoose the parameters
of the matrix suh that a11 − a21 = −(a22 − a12) = 1 so that the tness an be
rewritten as:
 
1 0
 
f(x) = ũ  x (5.25)
 
0 −1
Due to the simplex onstraint, (5.23) is ompletely desribed using x = x1:
ẋ(t) = ũ(t)x(t)(1 − x(t))(f 1(x)− f 2(x)) (5.26)
with x = 0, 1 orresponding to the pure strategies 2 and 1 respetively. Due to
1 2
our hoie of the payo matrix parameters, f − f is a onstant. This allows
132
1 2
us to adopt a time-sale hange by the fator f − f to arrive at our generative
model:
ẋ(t) = u(t)x(t)(1− x(t)) (5.27)
This dynamis results in asymptoti onvergene to the pure strategy x = 1 in
the absene of ontrol, that is, when u(t) ≡ 1. However, the time-varying on-
trol variable u serves to model hanging preferenes for the ok strategies by
appropriate hanges in its sign and magnitude. Suh a temporal modulation of
the tness ensures feasibility of apturing arbitrary traes in the interior of the
simplex.
{xd d d dGiven a set of data points 0, x1, ..., xN} with eah xk ∈ (0, 1), k = 0, 1, ..., N ,
at time instants {t0, t1, ..., tN}, we formulate the optimal ontrol problem,
∫ N
λ tN ∑
min J(x(t0), u) = u
2dt+ Fi(x(ti)),
x(t0), u∈R 2 t0 i=0
(5.28)
subjet to: ẋ = ux(1− x),
where the t errors Fi's are given by the Kullbak-Leibler divergene measure of
mismath between the data and the state,
( ) ( )
xd 1− xd
Fi(x) = x
d i
i log + (1− xd) log ii − , i = 0, 1, ..., N. (5.29)x 1 x
We an diretly appeal to Pontryagin's Maximum Priniple (PMP) and theo-
rem (5.3.1) to write neessary onditions for optimality. We an write the pre-
Hamiltonian as,
λ
H(x, p, u) = upx(1− x)− u2. (5.30)
2
133
The Hamiltonian maximization ondition (5.6) yields an optimal ontrol in eah
time interval t ∈ (ti, tt+1), i = 0, 1, ..., N − 1,
1
u = px(1− x), (5.31)
λ
with Hamiltonian given by,
H 1(x, p) = p2x2(1− x)2. (5.32)
2λ
Hamilton's equations (5.5) read,
1
ẋ = px2(1− x)2
λ
(5.33)
1
ṗ = − p2x(1− x)(1− 2x).
λ
The jump onditions for p (5.7) an be written as,
p(t−0 ) = 0,
x(t di)− x
p(t+i )− p(t−) = i , i = 0, 1, ..., N, (5.34)i x(ti)(1− x(ti))
p(t+N) = 0.
du
Remark 5.1. Note that the optimal ontrol is pieewise onstant sine = 0
dt
for eah of these time intervals t ∈ (ti, ti+1), i = 0, 1, ..., N − 1.
Therefore, denoting xk = x(tk), k = 0, 1, ..., N , any optimal ontrol an be
desribed by a vetor (u0, u1, ..., uN) with the onditions
1
u0 = (x0 − xd
λ 0
),
1
uk − uk−1 = (x dk − xk), k = 1, 2, ..., N (5.35)λ
uN = 0.
134
Pieewise onstany of the ontrol input allows us to write the solution to the state
equation (5.26) expliitly. Suppose the sampling time of the trae is uniform, i.e.
∆t := tk+1 − tk, ∀k ∈ {0, ..., N − 1}, integrating the state equation (5.26) in
(tk, tk+1), we an write
x euk∆tk
xk+1 = − , k = 0, 1, ..., N − 1. (5.36)1 + xk (euk∆t 1)
By iteration, we an in turn write every xk as a funtion of x0 and u0, u1, ..., uk−1,
x0e
(u0+u1+···+uk )∆t−1
xk = xk(x0) = , k = 1, 2, ..., N. (5.37)
1 + x (e(u0+u1+···+uk 1)∆t0 − − 1)
The endpoint ondition (uN = 0) an then be written as,
x + x + · · ·+ x = xd0 1 N 0 + xd1 + · · ·+ xdN , (5.38)
where the left hand side of (5.38) is a funtion of x0. Solving the optimal ontrol
problem (5.28) thus boils down to solving (5.38) for x0 ∈ (0, 1).
Remark 5.2. The value of the regularization parameter λ is usually hosen
through ross validation tehnique. We do not employ any suh tehniques here.
The value of λ is hosen suh that the root nding algorithm for solving (5.38)
onverges for all events. For λ = 0.2, the roots were found with reasonably good
−5
auray with value of the funtion at the root being of the order of 1 × 10 or
lower for all events. For lower λ however, the problem beomes stier and left
hand side of (5.38) demonstrates `eetive disontinuity' in x0. This poses seri-
ous problem in solving (5.38). It is to be noted that the original aptured ight
data was subjeted to data-smoothing to obtain smooth trajetories [Dey, 2015℄.
135
0.5 1
PSfrag replaements PSfrag replaements
ENS-I
ENS-I
SHP-I
0.4 0.8 SHP-I
ENS-II
ENS-II
SHP-II
SHP-II
0.3 0.6
0.2 0.4
0.1 0.2
0 0
0 2 4 6 8 0 2 4 6 8
Event Event
(a) Time Averaged Hamiltonian Integral (b) Time Averaged Total Cost
Figure 5.1: Hamiltonian signatures
The data-smoothing problem in [Dey, 2015℄ onsidered a linear generative model
as in Setion 5.3.1 and used ordinary ross validation for trajetory of eah bird
to determine the appropriate weight to the regularization term. This generated
smooth trajetories with suppressed level of noise ompared to the original data.
d d
We then take the sampled data {x0, · · · , xN} from these smooth trajetories. This
an justify taking same value of λ aross all the events. As a future step, ross
validation ould be employed to arrive at a good value of λ in the range where
(5.38) an be solved.
5.6 Data Fitting Results
For all 8 events, we solve the optimal ontrol problem (5.28) and report the
results here. The value of the regularization weight λ is taken to be 0.2 and 100
data samples at regular time intervals are taken for all events. Given the data
136
∫
1 Hdt
T
J/T
vetor, we solve equation (5.38) for x0 ∈ (0, 1). In Table 5.2, we report time
averaged Hamiltonian integrals and time averaged total osts for all the dierent
games that we onsider in eqs. (5.18) to (5.21). These time averaged Hamiltonian
integrals are thought of as ognitive osts of the events. As seen from Table 5.2,
the trend of (ENS-I) losely follow the game (SHP-I), while the other two games
seem to follow eah other. This is graphially represented in Fig. 5.1. Optimal
ontrol solutions for the games ENS-I (5.18) and SHP-II (5.21) for individual
events are shown in Fig. 5.25.9. We note that more variation in the energy time
signal results in higher ognitive ost (in both measures). This is interpreted as
the olleting having to `think' more to properly alloate the modes, inurring
higher osts. These ognitive osts for a partiular game an thus indiate overall
physial behavior of the ok. For example, in the games (ENS-II) or (SHP-II)
where a rotational energy is onsidered as one of the pure strategies, relatively
higher ognitive osts for event 2, 5 indiate that the oks went through more
rotations than the other events during the ight periods. On the other hand,
low ost for event 4 is justied by almost retilinear overall motion. Similar
onlusions an be drawn for the other set of games (ENS-I) and (SHP-II), where
the respetive ognitive ost will stipulate nature of variation of the demorati
(reshuing within the ok) energy. The higher the ost is, more aggressively the
relative positions of the birds within the oks are hanged, leading to a more
omplex ight event.
137
∫
∫Hdt J(∫x0,u)
Duration
dt dt
(seonds) (ENS-I) (SHP-I) (ENS-II) (SHP-II) (ENS-I) (SHP-I) (ENS-II) (SHP-II)
5.4875 0.1232 0.1263 0.0976 0.1077 0.1981 0.1975 0.1454 0.1499
1.8176 0.1432 0.1018 0.2210 0.1760 0.2227 0.1619 0.3769 0.3118
5.6118 0.2735 0.2392 0.0613 0.1073 0.4595 0.4092 0.1557 0.2495
2.3471 0.1021 0.1270 0.0107 0.0190 0.2440 0.2702 0.0594 0.0610
3.8824 0.0779 0.2699 0.1587 0.1383 0.0896 0.3692 0.3001 0.3041
4.1588 0.1809 0.1634 0.0846 0.1105 0.2799 0.2706 0.2063 0.2090
5.7353 0.0804 0.1293 0.0576 0.0619 0.1127 0.2079 0.1087 0.1221
1.7588 0.4569 0.4069 0.0731 0.1090 0.8037 0.8361 0.2074 0.3810
Table 5.2: Hamiltonian Signature
5.7 Disussion
In this hapter, we have brought together several results from geometry, optimal
data-tting and evolutionary game theory to assoiate a ognitive aspet to ok-
ing. The ight data of Starling oks give rise to time-signals of energy mode
distributions. Here, the whole ok is oneptualized to making deisions about
how to optimally alloate its energy in several modes. The dierent energy modes
are thought as pure strategies of an evolutionary game and their tness is mod-
ulated by some deision or ontrol variables. These ontrols are then determined
by optimally tting this model to the observed energy mode distributions in the
data. The ost to this data-tting are referred to as ognitive ost for the ok.
In this work, we have only onsidered splitting energy into two modes. In this
138
setting, the optimal ontrol solutions present interesting harateristis. It will
be an important diretion to onsider energy splitting in several energy modes,
hene solving the tting problem in a higher dimensional simplex. It will also be
of interest to interpret the ognitive osts in suh senarios.
139
(a) Flok Trajetory
1 3
0.9 2.5
0.8 2
0.7 1.5
0.6 1
0.5 0.5
0.4 0
Trajectory fit
0.3 Data -0.5
0.2 -1
0.1 -1.5
0 -2
0 1 2 3 4 5 6 0 1 2 3 4 5 6
Time (s) Time (s)
(b) Trajetory Fit () Optimal Control
1 2
Trajectory fit
0.9 Data
1.5
0.8
1
0.7
0.5
0.6
0.5 0
0.4
-0.5
0.3
-1
0.2
-1.5
0.1
0 -2
0 1 2 3 4 5 6 0 1 2 3 4 5 6
Time (s) Time (s)
(d) Trajetory Fit (e) Optimal Control
E
dem
Figure 5.2: Event 1, λ = 0.2, Number of samples = 100, (b)-() x = (ENS-I),
E
rel
E
rot
(d)-(e) x = (SHP-II)
E
rel
140
x, x x, x
d d
u u
(a) Flok Trajetory
1 1.5
0.9 1
0.8 0.5
0.7
0
0.6
-0.5
0.5
-1
0.4
-1.5
0.3
0.2 -2
0.1 -2.5Trajectory fit
Data
0 -3
0 0.5 1 1.5 2 0 0.5 1 1.5 2
Time (s) Time (s)
(b) Trajetory Fit () Optimal Control
1 2
Trajectory fit
0.9 Data 1.5
0.8 1
0.7
0.5
0.6
0
0.5
-0.5
0.4
-1
0.3
0.2 -1.5
0.1 -2
0 -2.5
0 0.5 1 1.5 2 0 0.5 1 1.5 2
Time (s) Time (s)
(d) Trajetory Fit (e) Optimal Control
E
dem
Figure 5.3: Event 2, λ = 0.2, Number of samples = 100, (b)-() x = (ENS-I),
E
rel
E
rot
(d)-(e) x = (SHP-II)
E
rel
141
x, x x, x
d d
u u
(a) Flok Trajetory
1 3
0.9
2
0.8
1
0.7
0
0.6 Trajectory fit
Data
0.5 -1
0.4
-2
0.3
-3
0.2
-4
0.1
0 -5
0 1 2 3 4 5 6 0 1 2 3 4 5 6
Time (s) Time (s)
(b) Trajetory Fit () Optimal Control
1 2.5
Trajectory fit
0.9 Data 2
0.8 1.5
0.7 1
0.6 0.5
0.5 0
0.4 -0.5
0.3 -1
0.2 -1.5
0.1 -2
0 -2.5
0 1 2 3 4 5 6 0 1 2 3 4 5 6
Time (s) Time (s)
(d) Trajetory Fit (e) Optimal Control
E
dem
Figure 5.4: Event 3, λ = 0.2, Number of samples = 100, (b)-() x = (ENS-I),
E
rel
E
rot
(d)-(e) x = (SHP-II)
E
rel
142
x, x x, x
d d
u u
(a) Flok Trajetory
1 2.5
0.9
2
0.8
1.5
0.7
1
0.6
0.5 0.5
0.4
0
0.3
-0.5
0.2
-1
0.1 Trajectory fit
Data
0 -1.5
0 0.5 1 1.5 2 2.5 0 0.5 1 1.5 2 2.5
Time (s) Time (s)
(b) Trajetory Fit () Optimal Control
1 0.8
Trajectory fit
0.9 Data
0.6
0.8
0.4
0.7
0.2
0.6
0.5 0
0.4
-0.2
0.3
-0.4
0.2
-0.6
0.1
0 -0.8
0 0.5 1 1.5 2 2.5 0 0.5 1 1.5 2 2.5
Time (s) Time (s)
(d) Trajetory Fit (e) Optimal Control
E
dem
Figure 5.5: Event 4, λ = 0.2, Number of samples = 100, (b)-() x = (ENS-I),
E
rel
(d)-(e) x = Erot (SHP-II)
E
rel
143
x, x x, x
d d
u u
(a) Flok Trajetory
1 1.5
Trajectory fit
0.9 Data
1
0.8
0.7 0.5
0.6
0
0.5
-0.5
0.4
0.3 -1
0.2
-1.5
0.1
0 -2
0 0.5 1 1.5 2 2.5 3 3.5 4 0 0.5 1 1.5 2 2.5 3 3.5 4
Time (s) Time (s)
(b) Trajetory Fit () Optimal Control
1 2
Trajectory fit
0.9 Data 1.5
0.8 1
0.7
0.5
0.6
0
0.5
-0.5
0.4
-1
0.3
0.2 -1.5
0.1 -2
0 -2.5
0 0.5 1 1.5 2 2.5 3 3.5 4 0 0.5 1 1.5 2 2.5 3 3.5 4
Time (s) Time (s)
(d) Trajetory Fit (e) Optimal Control
E
dem
Figure 5.6: Event 5, λ = 0.2, Number of samples = 100, (b)-() x = (ENS-I),
E
rel
E
rot
(d)-(e) x = (SHP-II)
E
rel
144
x, x x, x
d d
u u
(a) Flok Trajetory
1 3
0.9
2
0.8
0.7
1
0.6
0.5 0
0.4
-1
0.3
0.2
-2
0.1 Trajectory fit
Data
0 -3
0 1 2 3 4 5 0 1 2 3 4 5
Time (s) Time (s)
(b) Trajetory Fit () Optimal Control
1 1
Trajectory fit
0.9 Data
0.5
0.8
0
0.7
-0.5
0.6
0.5 -1
0.4
-1.5
0.3
-2
0.2
-2.5
0.1
0 -3
0 1 2 3 4 5 0 1 2 3 4 5
Time (s) Time (s)
(d) Trajetory Fit (e) Optimal Control
E
dem
Figure 5.7: Event 6, λ = 0.2, Number of samples = 100, (b)-() x = (ENS-I),
E
rel
x = Erot(d)-(e) (SHP-II)
E
rel
145
x, x x, x
d d
u u
(a) Flok Trajetory
1 1.5
0.9
1
0.8
0.5
0.7
0
0.6
0.5 -0.5
0.4
-1
0.3
-1.5
0.2
-2
0.1 Trajectory fit
Data
0 -2.5
0 1 2 3 4 5 6 0 1 2 3 4 5 6
Time (s) Time (s)
(b) Trajetory Fit () Optimal Control
1 2
Trajectory fit
0.9 Data
1.5
0.8
0.7 1
0.6
0.5
0.5
0
0.4
0.3 -0.5
0.2
-1
0.1
0 -1.5
0 1 2 3 4 5 6 0 1 2 3 4 5 6
Time (s) Time (s)
(d) Trajetory Fit (e) Optimal Control
E
dem
Figure 5.8: Event 7, λ = 0.2, Number of samples = 100, (b)-() x = (ENS-I),
E
rel
E
rot
(d)-(e) x = (SHP-II)
E
rel
146
x, x x, x
d d
u u
(a) Flok Trajetory
1 5
0.9
4
0.8
3
0.7
2
0.6
0.5 1
0.4
0
0.3
-1
0.2
-2
0.1 Trajectory fit
Data
0 -3
0 0.5 1 1.5 2 0 0.5 1 1.5 2
Time (s) Time (s)
(b) Trajetory Fit () Optimal Control
1 2
Trajectory fit
0.9 Data
1.5
0.8
0.7
1
0.6
0.5 0.5
0.4
0
0.3
0.2
-0.5
0.1
0 -1
0 0.5 1 1.5 2 0 0.5 1 1.5 2
Time (s) Time (s)
(d) Trajetory Fit (e) Optimal Control
E
dem
Figure 5.9: Event 8, λ = 0.2, Number of samples = 100, (b)-() x = (ENS-I),
E
rel
E
rot
(d)-(e) x = (SHP-II)
E
rel
147
x, x x, x
d d
u u
Chapter 6
Conlusions and Diretions for Future Re-
searh
We have made an attempt to explain olletive behavior in natural oks in this
thesis. Possible roboti appliations in this ontext are also presented. The thesis
is distintively divided into two parts depending on the underlying approah 
either olletive behavior is viewed as an emergene of interations between small
number of agents in a `bottom-up' fashion or those interations are inferred in a
`top-down' way. We summarize below the ontributions of this dissertation along
with diretions in whih this line of researh an be ontinued.
In Chapter 2, we explored inter-agent interation strategies from both theoret-
ial and implementation perspetives. First, we onsider a two-agent senario in
whih one agent pursues the other using onstant-bearing (CB) pursuit law. The
pursued agent behaves like a moving beaon whose movement is independent to
the other. It is then shown that under partiular parameter setting of the CB law
148
and onstant urvature paths of the beaon, the ombined equations of motion of
the system resemble to that of a ouple of gravitating partiles. Periodi orbits
are shown to exist, eah orresponding to a xed energy. This result is used in a
roboti appliation subsequently. We have eetively utilized the results of this
problem in the problem of enirling a stati beaon that is sensed visually by a
mobile robot by means of a amera with limited eld-of-view (FOV). Proper feed-
bak law for the robot is advised to make a desired losed loop in the phase spae
asymptotially stable. This guarantees intermittent appearane of the beaon in
the amera's FOV. Laboratory demonstration of this problem inorporates online
estimation of the beaon's position when it falls out of the FOV. Seondly, labo-
ratory implementations of two biologially plausible feedbak laws are presented.
These laws inlude another dual-agent law alled Mutual Motion Camouage and
a multi-agent onsensus type law alled Topologial Veloity Alignment. In this
hapter we have shown how omplex olletive motion patterns an emerge from
simple interations among the agents in a ok.
We study the problem of optimal steering of a single agent in Chapter 3. The
agent is driven from initial to nal onguration on the plane while minimizing
the ontrol ost that penalizes both speed and urvature ontrol. Optimal ontrol
solution is obtained by using Pontryagin's Maximum Priniple (PMP) and Lie-
Poisson redution tehnique. Optimal trajetories are ategorized by the values
of the Hamiltonian and another onserved quantity alled Casimir. This problem
is then extended to apture the senario of a ok of agents moving on the plane.
149
These agents interat with eah other by a predened graph. The individual
ontrol osts are oupled with mismath in ontrol with the neighbors. This
hapter forms a bridge between the two parts of the thesis.
In Chapter 4, we take the problem of optimal steering of many agent ase and
onsider its ontinuum limit. A goal of suh an approah is to develop wave equa-
tions that may explain observed phenomenon of information transfer in natural
oks. We only onsider the yli graph of interation that enables us to present
the results in a ompat way. A general optimal ontrol problem in the loop
group ase is stated. General ontrollability result in innite dimensional setting
is shown to be helpful to onstrut suh optimal ontrol problem. The neessary
onditions for optimality, namely the Pontryagin's Maximum Priniple (PMP) in
Hilbert spae setting is only valid under a ondition of nite o-dimensionality of
a reahable set. Two speial ases of this problem are studied. The ase in whih
the underlying group is the Heisenberg group H(3), i.e. a ontinuum of nonholo-
nomi integrators is studied in detail. We have shown that the optimal ontrol
solutions possess traveling wave harater. Moreover, a synhronization result is
obtained in whih the innite oupling strength prohibits every agent in the ok
to behave dierently. The ase of planar ontinuum, i.e. agents moving in the
speial Eulidean group SE(2) is also onsidered. Optimal ontrol evolution equa-
tions are obtained by both alulus of variations and PMP approah. Similar to
the H(3) ase, synhronization result is obtained. Numerial simulations for both
these ases are presented. However, we have not been able to perform a thorough
150
analytial study of the partial dierential equations obtained in the SE(2) ase.
It is one of the diretions in whih future researh ould be onduted. A ritial
question will answer whether traveling wave solutions exist in this ase. Further-
more, dierent interation graphs an be onsidered to unover more interesting
details about this problem.
Chapter 5 presents a data-driven analysis of ight data of European Starling
oks, aptured in Rome. This data gives rise to temporal signals of the ok's
energy distribution in several energy modes. We use an optimal ontrol based
data-tting tehnique to explain this data as the outome of an evolutionary
game on a simplex. We all the data-tting ost funtionals of the underlying
optimal ontrol problem as `ognitive ost' that measures the ognitive eort of
the ok to alloate its energy in dierent modes. In our work, we have only
onsidered energy splitting into two modes so as to onsider a simple game on
the one-dimensional simplex. This an be extended to higher dimensions where
multiple energy modes are onsidered.
151
Appendix A
An Optimal Control Problem in an Innite
Dimensional Setting
A.1 Introdution
After Pontryagin provided his method for solving optimal ontrol problems in -
nite dimensional setting [Pontryagin et al., 1962℄, there have been many attempts
to try and prove similar priniples in innite dimensions. However, the ounterex-
ample of Egorov [Egorov, 1963℄ posed a serious hallenge to that proess. This
ounterexample showed that the Pontryagin's maximum priniple does not gen-
erally hold in innite dimensions. In partiular, the ostate variable an beome
identially zero, making the maximum priniple empty. The advanements in the
following deades [Ekeland, 1979; Fattorini, 1987; Li and Yong, 2012; Krastanov
et al., 2011℄ showed that it is possible to state PMP in some ases where some
additional assumptions are made. In this work, we adopt a similar path to prove
152
the maximum priniple set in a muh friendlier setting.
We onsider an abstrat dierential equation in a Hilbert spae,
dx(t)
= f(t, x(t), u(t)), a.e. in [0, T ], (A.1)
dt
where x(t) ∈ X , u(·) ∈ U , and T > 0. Let X be a Hilbert spae alled the state
spae and U be the set of all measurable funtions u(·) : [0, T ] → U , where U is
a separable metri spae alled the ontrol spae. With this setup, we formulate
the following optimal ontrol problem (P),
∫ T
min J(u) = L(t, x(t), u(t))dt
u∈U 0
(P) (A.2)
subjet to: ẋ = f(t, x, u), a.e. in [0, T ], x(0) = x0, x(T ) = xT .
We assume that both the funtions f(·, ·, ·) and L(·, ·, ·) are Bohner integrable
in t ∈ [0, T ] and Lipshitz ontinuous in x(t) ∈ X , with onstant K. Further-
′
more, we require the existene of the ontinuous Fréhet derivatives fx(t, x, u) and
L′x(t, x, u)
′ ′
. We also assume the funtions f, L and their derivatives fx, Lx to be
bounded, i.e. there exists an M > 0, suh that
‖f(t, x, u)‖ ≤ M, ‖f ′x(t, x, u)‖ ≤ M,
(A.3)
‖L(t, x, u)‖ ≤ M, ‖L′x(t, x, u)‖ ≤ M,
for all (t, x(t), u(t)) ∈ [0, T ]× X × U . Note that these hypotheses ensure a on-
tinuous and unique solution of (A.1) to exist [Avez, 1986℄. Let the spae U be
endowed with the distane funtion,
d(u, v) = meas{t : u(t) 6= v(t)}, (A.4)
153
where meas{·} denotes the usual Lebesgue measure in [0, T ]. Then, aording
to Theorem 5.3 of [Fattorini, 1987℄, the spae U is omplete with respet to the
distane d. A diret onsequene of the assumptions stated above leads to an
important result [Krastanov et al., 2011℄.
Lemma A.1.1. Let u1(·), u2(·) be any two arbitrary members of U . Denote the
state trajetories assoiated with these ontrols by, xi(·) = x(·, ui(·)), i = 1, 2.
Then there exist positive onstants C1, C2 suh that,
sup ‖x1(t)− x2(t)‖ ≤ C1d(u1, u2), (A.5)
t∈[0,T ]
|J(u1)− J(u2)| ≤ C2d(u1, u2). (A.6)
Proof. Let S ⊂ [0, T ] be the set where the ontrols u1 and u2 dier, i.e. d(u1, u2) =
meas{S}. We know that (A.1) an also be written as,
∫ t
x(t) = x0 + f(s, x(s), u(s))ds.
0
Then,
∫ t
x1(t)− x2(t) = (f(s, x1(s), u1(s))− f(s, x1(s), u1(s))) ds
0
∫ t
= (f(s, x1(s), u1(s))− f(s, x2(s), u1(s))) ds
0
∫ t
+ (f(s, x2(s), u1(s))− f(s, x2(s), u2(s))) ds
0
∫ t
= (f(s, x1(s), u1(s))− f(s, x2(s), u1(s))) ds
0
∫
+ (f(s, x2(s), u1(s))− f(s, x2(s), u2(s))) ds.
[0,t]∩S
154
Taking norm of both sides we get,
∫ t
‖x1(t)− x2(t)‖ ≤ K ‖x1(s)− x2(s)‖ ds+ 2Md(u1, u2),
0
where the Lipshitz property of f in x is used in the rst term and the uni-
form boundedness is used in the seond term. By the use of Gronwall-Bellman
inequality, we arrive at (A.5). The result (A.6) an be derived analogously. 
We now proeed to solve this problem by a maximum priniple based approah.
Before we state the maximum priniple, let us give some additional tehnial
details that are going to be essential for the proof of the maximum priniple.
Denition A.1. (Finite Codimensionality) [Fattorini, 1987℄ A subset S of a
Hilbert spae Z is alled to be nite odimensional in Z, if there exists a losed
subspae Zc ⊆ Z of nite odimension suh that Sc = Π(o(S)), has nonempty
interior in Zc, where Πc denotes the orthogonal projetion from Z onto Zc and o
means losed onvex hull.
We will now make a key assumption to derive a nontrivial maximum priniple.
∗
Let a solution of problem (P) exist and that optimal ontrol is denoted as u ∈ U
∗
and let the orresponding optimal trajetory be denoted as x (t). Then dene the
`reahable set' as,
{ ∫ t
R := z(T ) ∈ X | z(t) = f ′ (s, x∗x (s), u∗(s)) · z(s)ds
0
t }∫
+ (f(s, x∗(s), v(s))− f(s, x∗(s), u∗(s))) ds, for some v(·) ∈ U (A.7)
0
(A1) The set R is nite odimensional in X .
155
A.2 Maximum Priniple
Using usual formalism, we invoke the pre-Hamiltonian funtion H : R×X ×U ×
R×X∗ → R as,
H(t, x(t), u(t), p0, p(t)) = p0L(t, x(t), u(t)) + 〈p(t), f(t, x(t), u(t))〉 , (A.8)
where p(t) ∈ X∗ is alled the ostate variable. Intuitively, we want to make the
following statement of the maximum priniple that needs to be validated.
∗
Theorem A.2.1. (Maximum Priniple) Let u ∈ U be an optimal ontrol for
∗
problem (P) and x (t) be the orresponding optimal trajetory. Then, there exist
(p∗, p∗(t)) ∈ R × X∗a pair 0 , t ∈ [0, T ] ∗ ∗, suh that (p0, p ) ≡6 (0, 0), p∗0 ≤ 0 p∗, (·)
satises the dierential equation,
ṗ∗(t) = − (f ′ ⋆(t, x∗(t), u∗(t)) p∗(t)− p∗L′ (t, x∗(t), u∗x 0 x (t)), (A.9)
⋆
where by A we denote the adjoint operator of the operator A. The pointwise
maximization of the pre-Hamiltonian holds,
H(t, x∗(t), u∗(t), p∗0, p
∗(t)) = maxH(t, x∗(t), v, p∗0, p
∗(t)), (A.10)
v∈U
∗ ∗
Moreover, x and p satisfy Hamilton's anonial equations, i.e.
dx∗ δH
= ∗ (t, x
∗, u∗, p∗ ∗0, p ) (A.11)dt δp
dp∗ −δH= (t, x∗, u∗∗ , p
∗, p∗0 ). (A.12)dt δx
0
Proof. At the outset, we begin by introduing the variable, x (t) ∈ R that obeys
 
x0(t)
 
0 0
 
the dynamis, ẋ = L(t, x, u), x (0) = 0. Dene, y(t) = ∈ R × X, so
 
x(t)
156
that
 
dy(t) L(t, x, u) 
=   =: g(t, y, u), a.e. in [0, T ]. (A.13)
dt  
f(t, x, u)
The ore of the argument in proving the maximum priniple will follow re-
sults of Ekeland [Ekeland, 1979℄, and tehniques developed in [Fattorini, 1987;
Krastanov et al., 2011; Li and Yong, 2012℄. We will now state a result known as
Ekeland variational priniple [Ekeland, 1979℄.
Lemma A.2.1. (Ekeland Variational Priniple) Let V be a omplete metri
spae with respet to the distane funtion d(·, ·) and let F : V → R ∪ {+∞} be
lower semiontinuous and bounded below with F 6≡ +∞. Let ǫ > 0 and u ∈ V be
suh that
F (u) ≤ inf{F (w) : w ∈ V }+ ǫ. (A.14)
Then there exists v ∈ V suh that
√
d(u, v) ≤ ǫ (A.15)
√
F (w)− F (v) ≥ − ǫ d(w, v), ∀w ∈ V. (A.16)
Let us proeed by assuming that an optimal ontrol to the problem (P) exists
∗ ∗ 0,∗ ∗
and is denoted by u and let y = (x , x ) be the orresponding optimal traje-
∗
tory. We write the minimum ost by η0, i.e. η0 = J(u ). Now, for given ǫ > 0, we
onsider the funtion Jǫ : U → R,
√
Jǫ(u) = (J(u)− η0 + ǫ)2 + ‖x(T )− xT ‖2X . (A.17)
157
It is evident that Jǫ(u) > 0 for all u ∈ U and ǫ > 0. Moreover,
J (u∗ǫ ) = ǫ ≤ inf Jǫ(u) + ǫ,
u∈U
ǫ
whih by the Ekeland variational priniple yields the existene of u ∈ U suh
that
√
d(uǫ, u∗) ≤ ǫ, (A.18)
√
J (w)− J (uǫǫ ǫ ) ≥ − ǫ d(w, uǫ), ∀w ∈ U . (A.19)
ǫ
Next, we introdue a variation in ontrol u what is known as needle variations".
For any v(·) ∈ U , let h : [0, T ] → R×X ,
 
L(t, xǫ(t), v(t))− L(t, xǫ(t), uǫ(t))
 
h(t) = (g(t, xǫ(t), v(t))− g(t, xǫ(t), uǫ(t))) =   .
 
f(t, xǫ(t), v(t))− f(t, xǫ(t), uǫ(t))
(A.20)
Then, aording to Corollary 3.9 (p. 144) of [Li and Yong, 2012℄ , for any ρ ∈ (0, 1],
there is a measurable set Fρ ⊂ [0, T ] suh that meas{Fρ} = ρT and
∫ t ∫
ρ h(s)ds = h(s)ds+ o(ρ), (A.21)
0 Fρ∩[0,t]
o(ρ)
where → 0 as ρ ↓ 0, uniformly in t ∈ [0, T ]. The perturbed ontrol is then
ρ
dened as,



uǫ(t), t ∈/ Fρ,
uǫρ(t) = . (A.22)


v(t), t ∈ Fρ
ǫ
It is of interest to express the perturbation in trajetory when the ontrol uρ is
ǫ ǫ
applied, i.e. we want a Taylor like expansion of yρ(t) = y(t, uρ) with respet to ρ
158
at ρ = 0 ǫ, i.e. at y (t) = y(t, uǫ). Let us write,
yǫρ(t) = y
ǫ(t) + ρψǫ(t) + o(ρ).
Then,
yǫρ(t)− yǫ(t)
ψǫ(t) = lim
ρ↓0 ρ
( )
∫ t ∫1 t
= lim g(s, yǫ(s), uǫρ ρ(s))ds− g(s, yǫ(s), uǫ(s))ds
ρ↓0 ρ 0 0
∫ t g(s, yǫ(s), uǫρ ρ(s))− g(s, yǫ, uǫρ(s))
= lim ds
ρ↓0 0 ρ
∫
g(s, yǫ(s), v(s))− g(s, yǫ, uǫ(s))
+ lim ds
ρ↓0 ρ
[0,t]∩Fρ
∫ t ∫ t
= g′y(s, y
ǫ(s), uǫ(s)) · ψǫ(s)ds+ (g(s, yǫ(s), v(s))− g(s, yǫ, uǫ(s))) ds,
0 0
(A.23)
′
where the seond term follows from (A.21). gy is the Fréhet derivative of g with
respet to y and an be deomposed as,
 
L′x(t, x
ǫ(t), uǫ(t)) · q
g′
 
(t, yǫy (t), u
ǫ(t)) · q̄ =   , for any q̄ = (q0, q) ∈ R×X.
 
f ′x(t, x
ǫ(t), uǫ(t)) · q
∣
ǫ ǫ ǫ ǫ d ǫ
∣
Let's write ψ (t) = (z0(t), z (t)). In partiular, we have z0(T ) = J(u ) anddρ ρ ρ=0
∣
zǫ(t) = d xǫ ∣ρ(t)) , eah of whih an be spelled out separately from (A.23),dρ ρ=0
∫ T ∫ T
zǫ(T ) = L′0 x(s, x
ǫ(s), uǫ(s)) · zǫ(s)ds+ (L(s, xǫ(s), v(s))− L(s, xǫ(s), uǫ(s))) ds,
0 0
(A.24)
∫ t ∫ t
zǫ(t) = f ′x(s, x
ǫ(s), uǫ(s)) · zǫ(s)ds+ (f(s, xǫ(s), v(s))− f(s, xǫ(s), uǫ(s))) ds.
0 0
(A.25)
159
ǫ ǫ
Now, as a next step, we derive neessary onditions for the pair (y (t), u (t)) to
be suboptimal. We do this by using the Ekeland variational priniple and letting
ρ tend to zero. In (A.19), we set w = uǫ ǫ ǫρ. Note that the ontrols u and uρ dier
only in the set Fρ, whih has a measure ρT . Then,
√
Jǫ(u
ǫ
ρ)− Jǫ(uǫ) ≥ − ǫd(uǫρ, uǫ)
≥ −√ǫρT.
i.e.
Jǫ(u
ǫ ǫ
ρ)− Jǫ(u ) ≥ − √T ǫ. (A.26)
ρ
Now, note that,
∣
Jǫ(u
ǫ
ρ)− J ǫǫ(u ) dJ ǫ ∣ǫ(uρ) ∣
lim =
∣
ρ↓0 ρ dρ ∣
ρ=0
[
1 dJ(uǫρ)
= 2(J(uǫ )− η0 + ǫ)
2Jǫ(uǫρ)
ρ dρ
∣
( ǫ )] ∣
∥ ∥ ∥ ∥ dxρ(T ) ∣
+2 ∥xǫ ǫ∥ ∥ ∥ρ(T )− xT xρ(T )− xT · ∣dρ ∣
ρ=0
〈 〉
(J(uǫ)− η0 + ǫ) xǫ(T )− xT
= zǫ0(T ) + , z
ǫ(T ) . (A.27)
J (uǫǫ ) Jǫ(uǫ)
Thus, taking the limit in (A.26), we an write,
√
ξǫzǫ0 0(T ) + 〈ξǫ, zǫ(T )〉 ≥ −T ǫ, (A.28)
ǫ ǫ
ξǫ = (J(u )−η0+ǫ) ξǫ = x (T )−xT ∈ X∗where 0 ǫ and ǫ . Note additionally that,Jǫ(u ) Jǫ(u )
(ξǫ)2 + ‖ξǫ‖20 = 1. (A.29)
160
(yǫEquation (A.28) an be regarded as the neessary onditions for (t), uǫ(t)).
Finally, we will let ǫ ∗ ∗tend to zero to obtain neessary onditions for (y (t), u (t))
to be optimal. Dene,
∫ T ∫ T
z0 := L
′
x(s, x
∗(s), u∗(s)) · z(s)ds+ (L(s, x∗(s), v(s))− L(s, x∗(s), u∗(s))) ds,
0 0
(A.30)
∫ t ∫ t
z(t) := f ′x(s, x
∗(s), u∗(s)) · z(s)ds + (f(s, x∗(s), v(s))− f(s, x∗(s), u∗(s))) ds.
0 0
(A.31)
Sine v(·) is any arbitrary element in U , z(T ) ∈ R, .f. (A.7).
Lemma A.2.2. The following results hold true.
lim |zǫ0(T )− z0| = 0,
ǫ↓0
(A.32)
lim sup ‖zǫ(t)− z(t)‖ = 0,
ǫ↓0 t∈[0,T ]
ǫ ∗
Proof. Let us denote Sǫ = {t ∈ [0, T ] : u (t) 6= u (t)}. Then, meas{Sǫ} =
√
d(uǫ, u∗) ≤ ǫ, by (A.18). From the denition of z(t), we nd,
∫ t
‖z(t)‖ ≤ M ‖z(s)‖ ds+ 2MT,
0
′
where boundedness of both f and fx have been used. Applying the Gronwall-
MT
Bellman inequality, we get ‖z(t)‖ ≤ 2MTe , for all t ∈ [0, T ]. We now write,
∫ t
zǫ(t)− z(t) = f ′ (s, xǫ(s), uǫ(s)) · (zǫx (s)− z(s))ds
0
∫ t
+ (f ′ (s, xǫx (s), u
ǫ(s))− f ′x(s, x∗(s), u∗(s))) · z(s)ds
0
∫ t
+ (f(s, xǫ(s), v(s))− f(s, x∗(s), v(s))) ds
0
∫ t
+ (f(s, xǫ(s), uǫ(s))− f(s, x∗(s), u∗(s))) ds.
0
161
′
Taking norm on both sides and utilizing boundedness of fx and z(t), we obtain,
∫ t
‖zǫ(t)− z(t)‖ ≤ M ‖zǫ(s)− z(s)‖ ds
0
∫ t
+ 2MTeMT ‖f ′x(s, xǫ(s), uǫ(s))− f ′ (s, x∗x (s), u∗(s))‖ ds
0
∫ t
+ ‖f(s, xǫ(s), v(s))− f(s, x∗(s), v(s))‖ ds
0
∫ t
+ ‖f(s, xǫ(s), uǫ(s))− f(s, x∗(s), u∗(s))‖ ds. (A.33)
0
The last term in (A.33) an be written as,
∫ t
‖f(s, xǫ(s), uǫ(s))− f(s, x∗(s), u∗(s))‖ ds
0
∫
= ‖f(s, xǫ(s), u∗(s))− f(s, x∗(s), u∗(s))‖ ds
[0,t]\Sǫ
∫
+ ‖f(s, xǫ(s), uǫ(s))− f(s, x∗(s), u∗(s))‖ ds
[0,t]∩Sǫ
∫
≤ K ‖xǫ(s)− x∗(s)‖ ds+ 2Md(uǫ, u∗)
[0,t]\Sǫ
≤ √ ǫ↓0(KC ǫ ∗1T + 2M)d(u , u ) ≤ (KC1T + 2M) ǫ −→ 0.
Note that we have used the Lipshitz property of f and result of Lemma A.1.1.
The seond and third term an be treated in a similar fashion to show they are of
o(1) whih goes to 0 as ǫ tends to 0. Note that, instead of Lipshitz ontinuity, we
′
would use ontinuity of fx in x in order to use appropriate upper bound. Hene,
(A.33) an be written as,
∫ t
‖zǫ(t)− z(t)‖ ≤ M ‖zǫ(s)− z(s)‖ ds+ o(1),
0
162
whih by Gronwall-Bellman inequality yields,
‖zǫ(t)− z(t)‖ ≤ eMT −ǫ↓0o(1) → 0,
uniformly in t ∈ [0, T ] ǫ. The onvergene of z0(T ) an be shown analogously. 
Using (A.28) and the Cauhy-Shwarz inequality, we may now write,
ξǫ0z0 + 〈ξǫ, z(T )〉 = ξǫ0zǫ0 + 〈ξǫ, zǫ(T )〉 − ξǫ ǫ0 (z0 − z0)− 〈ξǫ, zǫ(T )− z(T )〉
≥ − √T ǫ− |ξǫ0| |zǫ0 − z | − ‖ξǫ0 ‖ ‖zǫ(T )− z(T )‖
≥ − √T ǫ− |zǫ0 − z0| − ‖zǫ(T )− z(T )‖ . (A.34)
ǫ 2 ǫ 2 ǫ ǫ
The last inequality follows, sine (ξ0) + ‖ξ ‖ = 1, both |ξ0| ≤ 1 and ‖ξ ‖ ≤ 1.
√
ǫ ǫ ǫ
Denote, κ = −T ǫ − |z0 − z0| − ‖z (T )− z(T )‖ and by the onvergene results
ǫ
(A.32), we see that κ → 0 as ǫ ↓ 0. Thus, (A.34) an be expressed as,
ξǫz + 〈ξǫ0 0 , z〉 ≥ −κǫ, ∀z0 ∈ R, z ∈ R, (A.35)
ǫ ǫ↓0
where, κ −→ 0. Now the assumption (A1) that the set R is nite odimensional
ǫ ǫ
in X is going to be useful in proving nontriviality of the limit of the pair (ξ0, ξ )
as ǫ goes to 0. Here we state the following lemma from [Fattorini, 1987℄, as a
onsequene of nite odimensionality.
Lemma A.2.3. Let Q be a set of nite odimension in a Hilbert spae Z and let
{zn} be a sequene of vetors in Z suh that
0 < c ≤ ‖zn‖ ≤ C.
163
Assume additionally that,
〈zn, q〉 ≥ −θn,
for q ∈ Q and θn → 0 as n → ∞. Then there exists a subsequene of {zn} that
onverges weakly to z ∈ Z, and z 6= 0.
Now hoose a sequene {ǫ(n)} suh that ǫ(n) → 0 as n → ∞. Sine both
ǫ(n)
sequenes {ξ0 } {ξǫ(n) ǫ(nk)and } are bounded, there exist subsequenes {ξ } and
{ξǫ(nk)} ∗that onverge weakly to some ξ̄0 ∈ R and ξ̄ ∈ X . For simpliity, let the
subsequenes be denoted by themselves. Showing nontriviality of the pair (ξ̄0, ξ̄)
is a ruial step in proving maximum priniple in innite dimensional ase. Reall
ξǫ = (J(u
ǫ)−η0+ǫ)
that, 0 ǫ , so that ξ
ǫ > 0, ∀ǫ > 0. Hene, we may only have ξ̄ ≥ 0.
J 0ǫ(u ) 0
If ξ̄0 6= 0 ǫ, we are done proving that (ξ̄0, ξ̄) 6= (0, 0). Otherwise, let ξ0(n) → 0 as
( )
∥ ∥
→ ∞ ǫ(n) 2
2
ǫ(n)
n ∥ ∥. Then from the relation (A.29), we get, 1 ≥ ξ = 1 − ξ0 ≥
1 − δ > 0, for some δ > 0, for n large enough. Finally, by the lemma A.2.3, we
get ξ̄ 6= 0 in the ase ξ̄0 = 0. Hene we onlude that,
∗
(ξǫ ǫ0, ξ ) ⇀ (ξ̄0, ξ̄) 6= (0, 0), ξ̄0 ≥ 0. (A.36)
Then, nally taking the limit ǫ ↓ 0 in (A.35), we get for any z ∈ R and z0 as
speied in (A.30), there exists a pair R×X∗ ∋ (ξ̄0, ξ̄) 6= (0, 0), ξ̄0 ≥ 0, so that,
〈 〉
ξ̄0z0 + ξ̄, z ≥ 0. (A.37)
∗ ∗
Now, let us introdue the ostate variable p (t) ∈ X , that obeys the following
dierential equation,
ṗ∗
⋆
(t) = − (f ′x(t, x∗(t), u∗(t)) p∗(t)− p∗0L′x(t, x∗(t), u∗(t)), (A.38)
164
∗
with p (T ) = −ξ̄ and p∗0 = −ξ̄0 ≤ 0 ⋆. Here by A we denote the adjoint operator
of the operator A. Now sine z(0) = 0, we have,
〈p∗(T ), z(T )〉
= 〈p∗(T ), z(T )〉 − 〈p∗(0), z(0)〉
∫ T d
= 〈p∗(t), z(t)〉 dt
0 dt
∫ T
[〈
= − (f ′ ⋆
〉
x(t, x
∗(t), u∗(t)) p∗(t)− p∗L′0 x(t, x∗(t), u∗(t)), z(t)
0
+ 〈p∗(t), f ′ ∗ ∗x(t, x (t), u (t))z(t) + f(t, x∗(t), v(t))− f(t, x∗(t), u∗(t))〉] dt
∫ T
= −p∗0 〈L′x(t, x∗(t), u∗(t)), z(t)〉 dt
0
∫ T
+ (〈p∗(t), f(t, x∗(t), v(t))− f(t, x∗(t), u∗(t))〉) dt (A.39)
0
This, ombined with the denition of z0 (A.30) and equation (A.37) yields,
〈p∗(T ), z(T )〉+ p∗0z0
∫ T
= [H(t, x∗(t), v(t), p∗ ∗0, p (t))−H(t, x∗(t), u∗(t), p∗0, p∗(t))] dt ≤ 0. (A.40)
0
Sine the ontrol set U is separable, the similar argument as in [Krastanov et al.,
2011; Li and Yong, 2012℄ would give the pointwise maximization riterion of the
pre-Hamiltonian,
H(t, x∗(t), v(t), p∗0, p
∗(t)) ≤ H(t, x∗(t), u∗(t), p∗0, p∗(t)), a.e. in [0, T ], ∀v ∈ U .
(A.41)
From the denition of the Hamiltonian, we an nally ompute its derivatives. In
what follows, the appropriate arguments will be suppressed for simpliity and the
notation |∗ will imply the funtion has been evaluated at optimal parameters. We
165
then nd for any (x̃, p̃) ∈ X ×X∗,
∣ ∗ − ∗ ∣δH ∣ · H(p + σp̃) H(p )
∣
∣ p̃ = lim ∣
δp∗ ∣ ∣∗ σ→0 σ ∗
∣
d ∗ ∣= H(p + σp̃)∣
dσ ∣σ=0,∗
= 〈p̃, f |∗〉
and,
∣
δH H(x∗∣ + σx̃)− ∣H(x∗) ∣
∣ ∣
δx∗
· x̃ = lim
∣ σ→0 σ ∣∗ ∗
∣
d ∣
= H(x∗ + σx̃)∣
dσ ∣σ=0,∗
∣
d ∣
= (〈p∗, f(x∗ + σx̃)〉+ p∗0L(x∗ + σx̃)) ∣dσ ∣σ=0,∗
= 〈p∗, f ′x|∗ · x̃〉+ p∗0 〈L′x|∗, x̃〉
〈 ′ | ⋆
〉
= (fx ∗) p
∗ + p∗L′0 x|∗, x̃ .
Thus, we may write the anonial Hamilton's equations of motion,
dx∗ δH
= ∗ (t, x
∗, u∗, p∗ ∗
dt δp 0
, p ),
(A.42)
dp∗ −δH= ∗ (t, x
∗, u∗, p∗, p∗).
dt δx 0
This ompletes the proof of the maximum priniple.

166
Bibliography
Agrahev, A. and Sahkov, Y. 2004. Control Theory from the Geometri
Viewpoint. Springer, Berlin.
Attanasi, A., Cavagna, A., Del Castello, L., Giardina, I., Grigera,
T. S., Jeli¢, A., Melillo, S., Parisi, L., Pohl, O., Shen, E., and Viale,
M. 2014. Information transfer and behavioural inertia in starling oks. Nature
physis 10:691696.
Austin, M. A., Krishnaprasad, P. S., and Wang, L.-S. 1993. Almost
poisson integration of rigid body systems. Journal of Computational Physis
107:105117.
Avez, A. 1986. Dierential alulus. John Wiley & Sons In.
Ballerini, M., Cabibbo, N., Candelier, R., Cavagna, A., Cisbani, E.,
Giardina, I., Leomte, V., Orlandi, A., Parisi, G., Proaini, A.,
Viale, M., and Zdravkovi, V. 2008a. Interation ruling animal olletive
behavior depends on topologial rather than metri distane: Evidene from a
eld study. Proeedings of the National Aademy of Sienes 105:12321237.
Ballerini, M., Cabibbo, N., Candelier, R., Cavagna, A., Cisbani, E.,
Giardina, I., Orlandi, A., Parisi, G., Proaini, A., Viale, M.,
et al. 2008b. Empirial investigation of starling oks: a benhmark study in
olletive animal behaviour. Animal behaviour 76:201215.
Bellaïhe, A., Jean, F., and Risler, J. J. 1998. Geometry of nonholonomi
systems. In J.-P. Laumond (ed.), Robot Motion Planning and Control. Springer,
London.
Bialek, W., Cavagna, A., Giardina, I., Mora, T., Silvestri, E., Viale,
M., and Walzak, A. M. 2012. Statistial mehanis for natural oks of
birds. Proeedings of the National Aademy of Sienes 109:47864791.
Birkhoff, G. D. 1915. The restrited problem of three bodies. Rendionti del
Cirolo Matematio di Palermo 39:265334.
Bishop, R. L. 1975. There is more than one way to frame a urve. Amerian
Mathematial Monthly pp. 246251.
167
Bliss, G. A. 1946. Letures on the alulus of variations, phoenix siene series.
The University of Chiago Press, Chiago 37:16.
Brokett, R. 1973. Lie theory and ontrol systems dened on spheres. SIAM
Journal on Applied Mathematis 25:213225.
Byrd, P. F. and Friedman, M. D. 1971. Handbook of Ellipti Integrals for
Engineers and Sientists. Springer Berlin Heidelberg.
Cavagna, A., Cimarelli, A., Giardina, I., Parisi, G., Santagati, R.,
Stefanini, F., and Viale, M. 2010. Sale-free orrelations in starling oks.
Proeedings of the National Aademy of Sienes 107:1186511870.
Cavagna, A., Conti, D., Giardina, I., and Grigera, T. S. 2018. Propagat-
ing speed waves in oks: a mathematial model. Physial Review E 98:052404.
Chiu, C., Reddy, P. V., Xian, W., Krishnaprasad, P. S., and Moss, C. F.
2010. Eets of ompetitive prey apture on ight behavior and sonar beam
pattern in paired big brown bats, Eptesius Fusus. Journal of Experimental
Biology 213:33483356.
Chorin, A. J., Marsden, J. E., and Marsden, J. E. 1990. A mathematial
introdution to uid mehanis, volume 3. Springer.
Corbet, P. S. 1999. Dragonies: Behavior and Eology of Odonata. Comstok
Publishing Assoiates.
Cuker, F. and Smale, S. 2007. Emergent behavior in oks. IEEE Transa-
tions on automati ontrol 52:852862.
Davis, H. T. 1962. Introdution to nonlinear dierential and integral equations,
volume 971.Dover Publiation, In.
Dey, B. 2015. Reonstrution, Analysis and Synthesis of Colletive Motion. PhD
thesis, University of Maryland, College Park.
Dey, B. and Krishnaprasad, P. S. 2012. Trajetory smoothing as a linear
optimal ontrol problem. In Proeedings of Annual Allerton Conferene on
Communiation, Control, and Computing, pp. 14901497, Urbana-Champaign,
IL. IEEE.
Dey, B. and Krishnaprasad, P. S. 2014. Control-theoreti data smoothing.
In Proeedings of 53rd IEEE Conferene on Deision and Control, pp. 5064
5070. IEEE.
Dubins, L. E. 1957. On urves of minimal length with a onstraint on aver-
age urvature, and with presribed initial and terminal positions and tangents.
Amerian Journal of Mathematis 79:497516.
168
Ebin, D. G. and Marsden, J. 1970. Groups of dieomorphisms and the motion
of an inompressible uid. Annals of Mathematis pp. 102163.
Eells Jr, J. 1966. A setting for global analysis. Bulletin of the Amerian
Mathematial Soiety 72:751807.
Egorov, Y. V. 1963. The optimal ontrol in banah spae. Dokl. Akad. Nauk
SSSR 150:241244.
Ekeland, I. 1979. Nononvex minimization problems. Bulletin of the Amerian
Mathematial Soiety 1:443474.
Elsgol, L. D. 2012. Calulus of variations. Courier Corporation.
Euler, L. 1744. Methodus inveniendi lineas urvas maximi minimive proprietate
gaudentes sive solutio problematis isoperimetrii latissimo sensu aepti.
Fattorini, H. 1987. A unied theory of neessary onditions for nonlinear non-
onvex ontrol systems. Applied Mathematis and Optimization 15:141185.
Fraihard, T. and Sheuer, A. 2004. From Reeds and Shepp's to ontinuous-
urvature paths. IEEE Transations on Robotis 20:10251035.
Galloway, K., Justh, E., and Krishnaprasad, P. 2009. Geometry of yli
pursuit. In Proeedings of 48th IEEE Conferene on Deision and Control,
IEEE, pp. 36293643.
Galloway, K. S., Justh, E. W., and Krishnaprasad, P. S. 2013. Sym-
metry and redution in olletives: yli pursuit strategies. Proeedings of the
Royal Soiety A 469:20130264.
Gelfand, I. and Fomin, S. 1963. Calulus of variations. Revised English edition
translated and edited by Rihard A. Silverman. Prentie-Hall In., Englewood
Clis, NJ.
Ghose, K., Horiuhi, T. K., Krishnaprasad, P. S., and Moss, C. F.
2006. Eholoating bats use a nearly time-optimal strategy to interept prey.
PLoS Biology 4:865873.
Giaquinta, M. and Hildebrandt, S. 1996. Calulus of Variations: Calulus
of Variation. I. Springer-Verlag.
Goldstein, H., Poole, C., and Safko, J. 2001. Classial mehanis.Addison-
Wesley.
Halder, U. and Dey, B. 2015. Biomimeti algorithms for oordinated motion:
theory and implementation. In Proeedings of International Conferene on
Robotis and Automation, pp. 54265432. IEEE.
169
Halder, U., Justh, E. W., and Krishnaprasad, P. S. 2019a. Continuum
oking and ontrol. In preparation for submission.
Halder, U. and Kalabi, U. 2017. Time-optimal solution for a uniyle path
on se2 with a penalty on urvature. IFAC-PapersOnLine 50:63206325.
Halder, U., Raju, V., Krishnaprasad, P. S., and et al. 2019b. Cognitive
ost of oking: A geometri and hamiltonian perspetive. In preparation for
submission.
Halder, U., Shlotfeldt, B., and Krishnaprasad, P. S. 2016. Steer-
ing for beaon pursuit under limited sensing. In Proeedings of 55th IEEE
Conferene on Deision and Control, pp. 38483855. IEEE.
Heintze, E. and Liu, X. 1999. Homogeneity of innite dimensional isopara-
metri submanifolds. Annals of mathematis 149:149181.
Inada, Y. and Kawahi, K. 2002. Order and exibility in the motion of sh
shools. Journal of Theoretial Biology 214:371  387.
Justh, E. W. and Krishnaprasad, P. S. 2003. Steering laws and ontinuum
models for planar formations. In Proeedings of 42nd IEEE Conferene on
Deision and Control, pp. 36093615. IEEE.
Justh, E. W. and Krishnaprasad, P. S. 2004. Equilibria and steering laws
for planar formations. Systems & Control Letters 52:25  38.
Justh, E. W. and Krishnaprasad, P. S. 2005. Natural frames and interating
partiles in three dimensions. In Proeedings of 44th IEEE Conferene on
Deision and Control, pp. 28412846. IEEE.
Justh, E. W. and Krishnaprasad, P. S. 2006. Steering laws for motion
amouage. Proeedings of Royal Soiety A 462:36293643.
Justh, E. W. and Krishnaprasad, P. S. 2015a. Enlargement, geodesis,
and olletives. In International Conferene on Networked Geometri Siene
of Information, pp. 558565. Springer.
Justh, E. W. and Krishnaprasad, P. S. 2015b. Optimality, redution and
olletive motion. Proeedings of Royal Soiety A 471:20140606.
Justh, E. W. and Krishnaprasad, P. S. 2016. Subriemannian geodesis
for oupled nonholonomi integrators. In Proeedings of Amerian Control
Conferene, pp. 72817288. IEEE.
Kinet. Mirosoft kinet. http://www.xbox.om/en-US/xbox-one/
aessories/kinet-for-xbox-one. [Online℄.
170
Krastanov, M. I., Ribarska, N., and Tsahev, T. Y. 2011. A pontryagin
maximum priniple for innite-dimensional problems. SIAM Journal on Control
and Optimization 49:21552182.
Krishnaprasad, P., Marus, S. I., and Hazewinkel, M. 1983. Current
algebras and the identiation problem. Stohastis: An International Journal
of Probability and Stohasti Proesses 11:65101.
Krishnaprasad, P. S. 1993. Optimal ontrol and Poisson redution. Tehnial
Report 93-87, Institute for Systems Researh, University of Maryland.
Kudrolli, A., Lumay, G., Volfson, D., and Tsimring, L. S. 2008. Swarm-
ing and swirling in self-propelled polar granular rods. Physial review letters
100:058001.
LaValle, S. M. 2006. Planning Algorithms. Cambridge university press.
Li, X. and Yong, J. 2012. Optimal ontrol theory for innite dimensional
systems. Springer Siene & Business Media.
MATLAB. MATLAB Robotis Toolbox. http://www.mathworks.om/
produts/robotis/. [Online℄.
Mishiati, M. and Krishnaprasad, P. S. 2010. Motion amouage for ov-
erage. In Proeedings of Amerian Control Conferene, pp. 64296435. IEEE.
Mishiati, M. and Krishnaprasad, P. S. 2011. Mutual motion amouage
in 3D. In Proeedings of the 18th IFAC World Congress, pp. 44834488.
Mishiati, M. and Krishnaprasad, P. S. 2012. The dynamis of mutual
motion amouage. Systems & Control Letters 61:894903.
Mishiati, M. and Krishnaprasad, P. S. 2017. Geometri deompositions
of olletive motion. Proeedins of Royal Soiety A 473:20160571.
Mizutani, A., Chahl, J. S., and Srinivasan, M. V. 2003. Inset behaviour:
Motion amouage in dragonies. Nature 423:604604.
Mora, T. and Bialek, W. 2011. Are biologial systems poised at ritiality?
Journal of Statistial Physis 144:268302.
Nagy, M., Akos, Z., Biro, D., and Visek, T. 2010. Hierarhial group
dynamis in pigeon oks. Nature 464:890893.
Nghiem, A. T., Auvinet, E., and Meunier, J. 2012. Head detetion using
Kinet amera and its appliation to fall detetion. In Proeedings of 11th
International Conferene on Information Siene, Signal Proessing and their
Appliations, pp. 164169. IEEE.
171
Ohsawa, T. 2013. Symmetry redution of optimal ontrol systems and prinipal
onnetions. SIAM Journal on Control and Optimization 51:96120.
Olberg, R. M., Worthington, A. H., and Venator, K. R. 2000. Prey
pursuit and intereption in dragonies. Journal of Comparative Physiology A
186:155162.
Oliver, A., Kang, S., Wünshe, B. C., and MaDonald, B. 2012. Using
the Kinet as a navigation sensor for mobile robotis. In Proeedings of the
27th Conferene on Image and Vision Computing, pp. 509514. ACM.
OpenCV. http://openv.org/. [Online℄.
Pioneer. Pioneer P3-DX. http://www.mobilerobots.om/researhrobots/
pioneerp3dx.aspx. [Online℄.
Pontryagin, L. S., Boltyanskiy, V. G., Gamkrelidze, R. V., and
Mishhenko, E. F. 1962. Mathematial theory of optimal proesses. In-
tersiene.
Pressley, A. and Segal, G. B. 1986. Loop Groups. Clarendon Press.
Raju, V. and Krishnaprasad, P. S. 2018. A variational problem on the
probability simplex. In Proeedings of 57th IEEE Conferene on Deision and
Control, pp. 35223528. IEEE.
Rashevsky, P. 1938. About onneting two points of a ompletely nonholonomi
spae by admissible urve. Uh. Zapiski Ped. Inst. Libknehta 2:8394.
Reddy, P. V., Justh, E., and Krishnaprasad, P. S. 2006. Motion amou-
age in three dimensions. In Proeedings of 45th IEEE Conferene on Deision
and Control, pp. 33273332. IEEE.
Reeds, J. and Shepp, L. 1990. Optimal paths for a ar that goes both forwards
and bakwards. Pai Journal of Mathematis 145:367393.
ROS. Robot Operating System. http://www.ros.org. [Online℄.
ROS-ARIA. http://wiki.ros.org/ROSARIA. [Online℄.
Rund, H. 1966. The Hamilton-Jaobi theory in the alulus of variations: its
role in mathematis and physis. Krieger Pub Co.
Salehani, M. K. and Markina, I. 2014. Controllability on innite-dimensional
manifolds: a howrashevsky theorem. Ata Appliandae Mathematiae
134:229246.
Sussmann, H. J. and Tang, G. 1991. Shortest paths for the Reeds-Shepp ar:
A worked out example of the use of geometri tehniques in nonlinear optimal
ontrol. Tehnial Report 91-10, Rutgers Center for Systems and Control.
172
Svirezhev, Y. M. 1972. Optimum priniples in genetis. Studies on theoretial
genetis pp. 86102.
Thurrowgood, S., Moore, R. J., Sool, D., Knight, M., and Srini-
vasan, M. V. 2014. A biologially inspired, vision-based guidane system for
automati landing of a xed-wing airraft. Journal of Field Robotis 31:699727.
Topaz, C. M., Bertozzi, A. L., and Lewis, M. A. 2006. A nonloal ontin-
uum model for biologial aggregation. Bulletin of mathematial biology 68:1601.
Vásárhelyi, G., Virágh, C., Somorjai, G., Tarai, N., Szörényi, T.,
Nepusz, T., and Visek, T. 2014. Outdoor oking and formation ight
with autonomous aerial robots. In Proeedings of IEEE/RSJ International
Conferene on Intelligent Robots and Systems, pp. 38663873. IEEE.
Vion. Vion motion apture system. http://www.vion.om. [Online℄.
Visek, T., Czirók, A., Ben-Jaob, E., Cohen, I., and Shohet, O. 1995.
Novel type of phase transition in a system of self-driven partiles. Physial
Review Letters 75:12261229.
Wei-Liang, C. 1939. Uber systeme von linearen partiellen dierentialgleihungen
erster ordnung. Math. Ann 117:98105.
Young, G. F., Sardovi, L., Cavagna, A., Giardina, I., and Leonard,
N. E. 2013. Starling ok networks manage unertainty in onsensus at low
ost. PLoS omputational biology 9:e1002894.
YouTube. Implementation videos. http://ter.ps/ICRA2015. [Online℄.
Zhang, H.-P., Beâer, A., Florin, E.-L., and Swinney, H. L. 2010.
Colletive motion and density utuations in baterial olonies. Proeedings of
the National Aademy of Sienes 107:1362613630.
173