II. RELATED WORK
Numerical methods that estimate the RoA given a closed-
form expression of the system dynamics include maximal
Lyapunov functions (LFs) and linear matrix inequalities
(LMIs). Ellipsoidal RoA approximation via LMIs [4], [5] has
been used for mobile robots [6], [7], and LMI relaxations can
also approximate the RoA of polynomial systems [8]. LFs
constructed by restricting them to be sum-of-squares (SoS)
polynomials [9] have been used in building randomized trees
with LQR feedback [10], funnel libraries [11] and stability
certificates for rigid bodies [12].
Reachability analysis [3], i.e., computing a backward
reachable tube to obtain the RoA without shape imposition,
for computing RoAs of dynamical walkers [13], has been
combined with machine learning to maintain safety over
a given horizon [14]. GPs can learn barrier functions for
ensuring the safety of unknown dynamical systems [15].
Similarly, barrier certificates (BCs) can identify areas for
exploration to expand the safe set [16].
Machine learning can learn LFs by alternating between
a learner and a verifier [17], [18], or via stable data-driven
Koopman operators [19]. Rectified Linear Unit (ReLU) acti-
vated neural networks can learn robust LFs for approximated
dynamics [20]. The Lyapunov Neural Network [21] can
incrementally adapt the RoA’s shape given an initial safe set.
As an alternative, GPs can obtain a Lyapunov-like function
[22], or an LF can be synthesized to provide guarantee’s on
a controller’s stability while training [23].
GPs are a popular choice to reduce data requirements
while modeling dynamical systems [24]. Some of their
applications in robotics include model-based policy search
[25], modeling non-smooth dynamics of robots with contacts
[26], and stabilizing controllers for control-affine systems
[27]. For RoA estimation problems, given an initial safe set
computed using a Lyapunov function, a GP can approximate
the model uncertainties on a discrete set of sampling points
from the safe region while expanding it [28].
Topology has multiple applications in robotics, such as
deformable manipulation and others [29]–[34]. Morse theory
can help incrementally build local minima trees for multi-
robot planning [35] and finds paths to cover 2D or 3D
spaces [36]. In recent work [37], Morse graphs are shown
to be effective in compactly describing the global dynamics
of a control system without an analytical expression of
its dynamics. To the best of the authors’ knowledge, the
current work is the first to apply surrogate modeling with
uncertainty quantification in conjunction with topological
tools to identify the global dynamics of robot controllers.
III. PROBLEM SETUP
This work aims to provide a data-efficient framework for
the analysis of global dynamics of robot controllers based on
combinatorial dynamics and order theory [37]–[40]. Consider
a non-linear, continuous-time control system:
˙x=f(x, u),(1)
where x(t)∈X⊆RMis the state at time t,Xis a compact
set, u:X7→ U⊆RMis a Lipschitz-continuous control
as defined by a deterministic control policy u(x), and f:
X×U7→ RMis a Lipschitz-continuous function. Neither
f(·)nor u=u(x)are necessarily known analytically. For
a given time τ > 0, let φτ:X→Xdenote the function
obtained by solving Eq. (1) forward in time for duration τ
from everywhere in X. A trajectory (or an orbit) is defined
as a sequence of states obtained by integrating Eq. 1 forward
in time.
The analysis of the global dynamics can reveal the sys-
tem’s attractors, which include fixed points, such as a state
that the control law manages to bring the system to; or limit
cycles, such as a periodic behavior of the system. It will also
reveal a Region of Attraction (RoA) which is a subset of the
basin of attraction of an attractor A. The basin of attraction
is the largest set of points whose forward orbits converge to
A, or more formally, the maximal set Bthat has the property:
A=ω(B) := \
n∈Z+
cl ∞
[
k=n
φk
τ(B)!
where φk
τis the composition φτ◦ · · · ◦ φτ(ktimes) and cl
is topological closure.
Since fand uare Lipschitz-continuous, φτis too; fur-
thermore any RoA of Eq. (1) is an RoA under φτ. Hence, it
is possible to study Eq. (1) by analyzing the behavior of the
dynamics according to φτ, which is not assumed, however,
to be computable and available.
IV. TOPOLOGICAL FRAMEWORK AND UNCERTAINTY
QUANTIFICATION VIA GPS
There are two key components for capturing meaningful
conditions of the dynamics according to φτ. First, identifying
effective combinatorial representations of the attractors and
maximal RoAs. And second, to achieve data efficiency by
employing GP-based surrogate modeling with uncertainty
quantification.
Morse Graphs for Understanding Global Dynamics:
Fig. 2(left) is used as a running 1-dim. example. The function
φτis first approximated by decomposing the state space X
into a collection of regions X, for instance, by defining
agrid. Fig. 2(left) shows a grid on the interval [−3,3]
decomposed into sub-intervals athrough e. Given a region
ξ∈ X (a cell), the system is forward propagated for
multiple initial states within ξfor a time τto identify regions
reachable from ξ. Consider, for example, the sub-interval b
as such a cell ξin Fig. 2(left). The arrows from the boundary
of bdepict the forward propagation of the dynamics where
bmaps to itself given the underlying dynamics
Then, a directed graph representation Fstores each region
in ξ∈ X as a vertex and edges pointing from ξto each
region reachable from ξ. In Fig. 2(left), Fis the graph
containing nodes given by the grid cells ato e. Each edge
represents a pair given by a cell and its image according
to the dynamics. For instance, (b,b)and (b,c)are two
edges added since bboth maps to itself and also maps to c.
Condensing all the nodes belonging to a strongly connected
components (SCCs) of Finto a single node, results in the
condensation graph CG(F). Then, edges on CG(F)reflect