Stabilization of slow-fast systems at fold points
H. Jardón-Kojakhmetov a , Jacquelien M.A. Scherpen a , D. del Puerto-Flores b ,
a

arXiv:1704.05654v2 [cs.SY] 25 Apr 2017

b

University of Groningen

Universidad de Guadalajara

Abstract
In this document, we deal with the stabilization problem of slow-fast systems (or singularly perturbed Ordinary Differential
Equations) at a non-hyperbolic point. The class of systems studied here have the following properties: 1) they have one fast
variable and an arbitrary number of slow variables, 2) they have a non-hyperbolic singularity of the fold type at the origin. The
presence of the aforementioned singularity complicates the analysis and the controller design of such systems. In particular, the
classical theory of singular perturbations cannot be used. We show a novel design process based on geometric desingularization,
which allows the stabilization of a fold point of singularly perturbed control systems. Our results are exemplified on an electric
circuit.
Key words: Nonlinear control; Slow-fast systems; singular perturbations.

1

Introduction

Slow-fast systems are characterized by having more than
one timescale. There are many phenomena in nature that
behave in two or more timescales such as population dynamics, cell division, electrical circuits, power networks,
chemical reactions, neuronal activity, etc. [11, 23, 26,
30, 45]. A particular property of slow-fast systems (under certain hyperbolicity conditions) is that they have a
structure suitable for model order reduction. Simply put,
certain slow-fast systems can be decomposed into two
simpler subsystems, the slow and the fast. The analysis
of those two subsystems allows a complete understanding of the more complex and higher dimensional one. A
mathematical theory supporting the previous fact is Geometric Singular Perturbation Theory (GSPT) [9], see
also [25]. The above, however, relies on the strong assumption of global timescale separation. When that does
not hold, then the classical technique of model order reduction cannot be employed. Thus, new mathematical
tools need to be introduced in order to deal with problems without global timescale separation.
Regarding the latter situation, many interesting phenomena are characterized by not having a global
Email addresses: h.jardon.kojakhmetov@rug.nl (H.
Jardón-Kojakhmetov), j.m.a.scherpen@rug.nl
(Jacquelien M.A. Scherpen),
dunstano.delpuerto@academicos.udg.mx (D. del
Puerto-Flores).

Preprint submitted to Automatica

timescale separation. This means that the variables
of the system do not always have the same timescale
relation throughout the phase-space. Mathematically
speaking, this phenomenon is characterized by singularities of the critical manifold, see Section 2, where the
timescale separation does not hold; and in a qualitative
sense, one usually observes jumps in the phase portrait
of the slow-fast system. Such an effect is also called loss
of normal hyperbolicity. Prototypical examples of twotimescale systems without global timescale separation
are the van der Pol oscillator [53], neuronal models [45],
and electrical circuits with impasse points [5, 6, 39].
Since loss of normal hyperbolicity is present in many
models, there is an increased need of their accurate understanding. The distinction between hyperbolic 1 and
non-hyperbolic 2 slow-fast systems is much more than
qualitative. From an analysis point of view, the available techniques are quite different; while the hyperbolic
case is well established, the non-hyperbolic scenario still
presents many challenges.
In the context of control systems, hyperbolic slow-fast
systems and the related model order reduction are nowadays well understood and have been used in many applications, e.g. [12, 23, 42, 43, 47]. The main and powerful
idea in the classical hyperbolic context is to design controllers for the reduced subsystems, which later guaran1
2

Global time scale separation.
No global timescale separation.

26 April 2017

tee stability of the overall slow-fast system. On the other
hand, non-hyperbolic slow-fast system are far from being well understood, and to the best knowledge of the
authors, there is not much progress yet. From a dynamical systems perspective, the technique called geometric
desingularization [8, 27] has been used to understand
the complex behavior of slow-fast systems around nonhyperbolic points. In this regard, the authors have made
preliminary progress in bringing such technique to the
control systems community [19, 20] in the planar case.

The parameter ε is always assumed 0 < ε  1. Let N be
an n-dimensional manifold,
a vector field X : N → T N
Pn
∂
is written as X = i=1 Xi ∂x
, where (x1 , . . . , xn ) is a
i
coordinate system in N .

The main contribution of this document is the development of a control design method based on geometric
desingularization. This provides a solution to the problem of the stabilization of non-hyperbolic singular perturbation problems. Although the technique presented
here is completely different from the classical one [25],
the idea remains the same: to obtain simpler subsystems where the control design becomes more accessible.
Furthermore, this document also generalizes our preliminary results of [19, 20] to systems with an arbitrary
number of slow variables and one fast variable.

where x ∈ Rns (slow variable), z ∈ Rnf (fast variable),
f and g are sufficiently smooth functions, and the independent variable is the slow time t. One can also define
a new time parameter τ = εt called the fast time, and
then (1) is rewritten as

A slow-fast system (SFS) is a singularly perturbed ordinary differential equation of the form
ẋ = f (x, z, ε)
εż = g(x, z, ε),

x0 = εf (x, z, ε)
z 0 = g(x, z, ε),

(2)

where the prime 0 denotes derivative with respect to τ .
Note that (1) and (2) are equivalent as long as ε > 0. For
convenience of notation we refer to (2) as Xε , and write
∂
(2) the ε-family of vector fields Xε = εf (x, z, ε) ∂x
+
∂
g(x, z, ε) ∂z
.

The rest of this document is organized as follows: in
Section 2 we provide preliminary information regarding
slow-fast systems followed by a description of the geometric desingularization method in Section 3. In Section 4 the specific problem under study is stated. Next,
Section 5, we apply the geometric desingularization to
control systems. Afterwards, in Section 6, we develop a
controller based on the method previously introduced.
Interestingly, we show that it is possible to inject a hyperbolic behavior to a non-hyperbolic slow-fast control
system even though the fast variable is not actuated.
Later, in Section 7 we exemplify our results for an electrical circuit. We finish in Section 8 with some concluding
remarks and a digression on open problems and future
work.
2

(1)

A first step towards understanding the dynamics of (1)
or (2) is to consider the limit equations when ε → 0. In
such a limit (1) becomes a Differential Algebraic Equation (DAE) of the form
ẋ = f (x, z, 0)
0 = g(x, z, 0).

(3)

On the other hand Xε becomes
X0 = 0

Preliminaries

∂
∂
+ g(x, z, 0) ,
∂x
∂z

(4)

which is called the layer equation. In principle, (3) and (4)
are not related anymore, however, the critical manifold
draws a bridge between them.

Abbreviations: SFS stands for Slow-Fast System,
SFCS for Slow-Fast Control System, and NH for Normally Hyperbolic.

Definition 1 The critical manifold is defined as

Notation: R, Z and N denote the fields of real, integer, and natural numbers respectively. Given a field F,
Fn denotes the n-cross-product F × · · · × F. The symbols R+ and R≥0 are respectively used to denote the
positive and the non-negative real numbers. In denotes
the n-dimensional identity matrix. The dimension of
the slow and fast variables is ns and nf respectively,
and N = ns + nf . Given x = (x1 , . . . , xn ) ∈ Rn and
αn
1
α = (α1 , . . . , αn ) ∈ Zn , we denote xα = (xα
1 , . . . , xn ).
n
The symbol S denotes the n-sphere. Given a matrix
A ∈ Rn×m , [A]i denotes its i-th row. Whenever a matrix A is square, |A| denotes its determinant. The vector
ê1 denotes the canonical vector ê1 = (1, 0, . . . , 0) ∈ Rn .

S = {(x, z) ∈ Rns × Rnf | g(x, z, 0) = 0} .

(5)

Note that S is ns -dimensional and serves as the phasespace of the DAE (3) and as the set of equilibrium points
of the layer equation X0 . A very important property of
critical manifolds is normal hyperbolicity.
Definition 2 A point s ∈ S is called hyperbolic if it is
a hyperbolic equilibrium point of X0 . The manifold S is
called normally hyperbolic (NH) if each point s ∈ S is a
hyperbolic equilibrium point of X0 .
2

flow along S is given by

The importance of NH critical manifolds is explained by
Geometric Singular Perturbation Theory (GSPT) [9, 21,
22, 30]. Briefly speaking, the following hold.

ẋ = f (x, H(x), us (x)).

• Let S0 ⊆ S be a compact, NH invariant manifold of
a SFS. Then, for ε > 0 sufficiently small, there exists
a locally invariant manifold Sε which is diffeomorphic
to S0 and lies within distance of order O(ε) from S0 .
• The flow of the SFS Xε along Sε converges to the flow
of the DAE along S0 as ε → 0.
• Sε has the same stability properties as S0 .
• Sε is in general not unique. However, any two perturbations of S0 lie within exponentially small distance
O(exp(−c/ε)), c > 0, from each other. Any of the representatives of Sε is called the slow manifold

Note that the reduced system is, as expected, independent of the fast variable z and of the fast controller
uf . Thus, the slow controller us is designed to make
{x = 0} an asymptotically equilibrium point of (7). After us has been designed, one studies the layer problem z 0 = g(x, z, us (x) + uf (x, z)), where x is taken as a
fixed parameter, and where the fast controller is designed
so that z = h(x, us (x)) is a set of asymptotically stable equilibrium points. The previous strategy plus some
extra (technical) interconnection conditions guarantee
that the origin (x, z) = (0, 0) ∈ Rns × Rnf is an asymptotically stable equilibrium point of the closed-loop system (6) for ε > 0 but sufficiently small [14, 24, 25, 33].
This feature of normal hyperbolicity has been exploited
in many applications, a few examples are [10, 34, 38, 43,
51, 52, 54, 55].

On the other hand, we say that a SFS is non-hyperbolic
if there exists a point p ∈ S such that the matrix ∂g
∂z (p)
has at least one eigenvalue with zero real part. The loss
of normal hyperbolicity can be related to jumps or rapid
transitions in, e.g., biological systems, climate models,
chemical reactions, nonlinear electric circuits, or neuron
models [7, 26, 28, 40, 41, 53].

3

The contribution of this paper is the stabilization of
a non-hyperbolic point of a control system with two
timescales. Before presenting the main tool to be used,
let us recall the classical (hyperbolic) setting of control
of singularly perturbed systems.
2.1

The Geometric Desingularization method

It is evident that the composite control strategy described in Section 2.1 is not applicable around nonhyperbolic points. There have been already some efforts
to study the stabilization problem of non-hyperbolic
points of nonlinear systems, see [32, 35]. However, these
do not address nonlinear two timescales systems, where
many open problems remain to be solved. We propose
to use the geometric desingularization method to design
controllers in a rather simple and standard way for multiple timescale systems around non-hyperbolic points.

Composite control of slow-fast control systems

Let us define a slow-fast control systems (SFCS’s) as
ẋ = f (x, z, ε, u)
εż = g(x, z, ε, u),

(7)

(6)

Before going into details, note that Xε (2) is an εparameter family of N -dimensional 3 vector fields. For
their analysis it is more convenient to lift such family
up and consider instead a single (N + 1)-dimensional
vector field defined as

where x ∈ Rns , z ∈ Rnf , and u ∈ Rm is a control input.
Let us now briefly recall a classical method to design
controllers for hyperbolic SFCS’s [25]. Assume that the
critical manifold of (6) is normally hyperbolic in a compact set (x, z, u) ∈ Ux × Uz × Uu ⊂ Rns × Rnf × Rm .
The goal is to design a control u that stabilizes the origin
(x, z) = (0, 0) ∈ Rns × Rnf for ε > 0 sufficiently small.
The hyperbolicity assumption implies that the critical
manifold S = {g(x, z, u) = 0} can be locally expressed
as a graph z = h(x, u) for (x, u) ∈ Ux × Uu . The idea
of composite control is to design u as a sum of two simpler controllers, namely u = us + uf where us = us (x) is
“the slow controller” and uf = uf (x, z) is “the fast controller”. When designing u, uf must be chosen so that it
does not destroy the normal hyperbolicity of the system,
meaning that g(x, z, us (x)+uf (x, z)) = 0 must have (locally) a unique root z = H(x). Thus, it is often required
that the effects of the fast controller uf disappear along
the S, that is uf |S = 0. In this way z = h(x, us (x)) is
(locally) a unique root of g(x, z, u) = 0 and the reduced

 0
x
X : z0
 0
ε

= εf (x, z, ε)
= g(x, z, ε)
=0

(8)

The geometric desingularization method, also known as
blow up, is a geometric tool introduced in [8] for the
analysis of SFSs around non-hyperbolic points, see also
[15, 18, 26, 28–30, 48, 49]. In an intuitive way, the blow
up method transforms non-hyperbolic points of SFSs to
(partially) hyperbolic ones. Below we provide a brief introduction to the technique, further details can be found
in [8, 30] .
3

3

Recall that in this document N = ns + nf .

Definition 3 A (quasi-homogeneous) 4 blow up transformation is a map
Φ : SN × R → RN +1
Φ(x̄, z̄, ε̄, r̄) 7→ (r̄α x̄, r̄β z̄, r̄γ ε̄),

with their vector fields, have to be “patched” together
via transition maps, defined as follows.
Definition 6 Let κv and κw be charts such that κv ∩
κw 6= ∅. A transition map φvw is a strictly positive
smooth function, well defined in κv ∩ κw , which makes
the following diagram commute.

(9)

where (x̄, z̄, ε̄) ∈ SN , α ∈ Zns , β ∈ Znf and γ ∈ Z.goo
In many applications, and in particular in this paper, it
is enough to consider r̄ ∈ R≥0 , which implies ε̄ ∈ R+ .
Thus, let B = SN × R≥0 , and Z = SN × {0}. A blow up
induces a vector field on B as follows.

κv

Definition 4 Let X : RN +1 → T RN +1 ' RN +1 be a
smooth vector field, and Φ : B → RN +1 a blow up map.
The blow up of X is a vector field X̄ : B → T B induced
by Φ in the following sense
X̄ = DΦ−1 ◦ X ◦ Φ.

B

κw

φvw

RN +1

RN +1

A schematic of the above definitions is shown in Figure
1.

(10)
In practice, we do not need to study all the charts covering B but rather a few of them. A good way to choose
the relevant charts is to first perform a qualitative analysis of the dynamics of a SFS to realize which charts
provide useful local information. After one has studied
the vector fields in the charts and has connected the
information in “adjacent” charts, one performs a blow
down, the inverse of the blow up. In summary, the steps
to follow when studying the dynamics of a SFS near a
non-hyperbolic point and via the geometric desingularization method are as follows.

It may happen that the vector field X̄ degenerates along
Z. In such a case one defines the desingularized vector
field X̃ as
1
X̃ = m X̄ ,
(11)
r̄
for some well suited m ∈ N so that X̃ is not degenerate,
and is well defined along Z. Note that the vector fields
X̄ and X̃ are equivalent on SN × {r̄ > 0}. Moreover, if
the weights (α, β, γ) are well chosen, the singularities of
X̃ |r̄=0 are partially hyperbolic or even hyperbolic, making the analysis of X̃ simpler than that of X . Due to the
equivalence between X and X̃ , one obtains all the local
information of X around 0 ∈ RN +1 from the analysis of
X̃ around Z.

(1) Find a suitable blow up. This is, find the proper
weights of the blow up map and rescaling so that
the non-hyperbolic point is desingularized.
(2) Find which charts are relevant and express the blow
up in the respective local coordinates.
(3) Investigate the dynamics of the local vector fields in
the charts.
(4) Connect the local results via transition maps and
blow down.

While doing computations, it is more convenient
to study the vector field X̃ in charts. A chart is a
parametrization of a hemisphere of B.
Definition 5 A chart is obtained by setting one of the
coordinates (x̄, z̄, ε̄) ∈ SN to ±1 in the definition of Φ.
If such coordinate is, for example, x̄i , then the chart is
denoted as κ±x̄i .

r3
x3
ε3
z̄

For example, the chart κε̄ is defined by the coordinates (r̄α x̄, r̄β z̄, r̄γ ), while the chart κ−x̄1 by
(−r̄α1 , r̄α2 x̄2 , . . . , r̄αn x̄n , r̄β z̄, r̄γ ε̄). Similarly to the blow
up map, a chart κ±w induces a vector field to be denoted
by X̃ ±w . Note that each chart covers only part of B but
all of them define an open cover of B. Therefore, one
uses several charts in order to have a full understanding
of the dynamics of X̃ around B. These charts, together

z2
x̄

z1

ε̄

r2
ε2

x1
r1

4

Fig. 1. Illustration of a blow up and three charts with corresponding local coordinates.

A homogeneous blow up (or simply blow up) refers to all
the exponents α, β, γ set to 1.

4

4

A1. The matrix B is invertible. Therefore m = ns . Such
an assumption implies that the slow-dynamics are
fully actuated.
A2. f (x, z, ε) = f0 (x, z) + εf1 (x, z, ε) with f1 (0, 0, ε) =
0. This ensures that ε does not introduce any constant drift near the origin.

Setting of the problem

The forthcoming sections are dedicated to the stabilization of a particular class of SFCS. We consider systems
with one fast variable (nf = 1) and an arbitrary number of slow variables (ns ). Preliminary results for the
case (ns , nf ) = (1, 1) are presented in [19, 20]. Here we
extend such results by allowing an arbitrary number of
slow variables. To start, we consider that the critical
manifold has a particular geometric structure.

Remark 8 The fully actuated case is in fact solvable by
a minor modification of the composite control method: 1)
propose the controller u = Uh (z) + Uc (x, z), where the
purpose of Uh is to make the closed-loop critical manifold
NH. 2) Once the critical manifold is NH, one proceeds by
designing Uc as in the composite control (recall Section
2.1).

Definition 7 (Fold point) Consider the SFS (8) and
let p = (x0 , z0 ) ∈ S. If
g(x0 , z0 , 0) = 0,
∂2g
(x0 , z0 , 0) 6= 0,
∂z 2

∂g
(x0 , z0 , 0) = 0,
∂z

5

(12)

In this section we perform the geometric desingularization of the system

then p is called a fold point of the SFS.
We shall study a SFCS whose critical manifold (in open
loop) has a fold point at the origin. We do this by defining
the system
ẋ = f (x, z, ε, u)
εż = g(x, z, ε, u),

 0
x
X : z0
 0
ε

(13)

= ε (f (x, z, ε) + B(x, z, ε)u)
= −(z 2 + x1 )
= 0,

(15)

which is equivalent to (14). Notice that X is a vector field
on Rns +2 . Without loss of generality we can write f =
A + L(x, z) + F (x, z, ε), where A = f (0, 0, 0), L(x, z) is a
linear map, i.e., L(x, z) = L1 x + L2 z with L1 ∈ Rns ×ns ,
L2 ∈ Rns , and F (x, z, ε) stands for all the higher order
terms and satisfies, due to Assumption A2, F (0, 0, ε) =
0. Thus, (15) is rewritten as

where x = (x1 , . . . , xns ) ∈ Rns , z ∈ R, u ∈ Rm , and g
satisfies 5 g(x, z, 0, 0) = −(z 2 +x1 ). The choice of (13) is
due to: 1) Fold singularities are generic in one-parameter
families of smooth functions [2, 3], therefore, they are
expected to appear in nonlinear SFSs with at least one
slow variable. 2) In many applications, jumps and rapid
transitions occur through fold points, the best-known
example are relaxation oscillations in the van der Pol
oscillator [53], although other examples can be found
in e.g. quantum electronics [37], neuroscience [7, 45] or
biochemistry [49].

 0
x
X : z0
 0
ε

= ε (A + L(x, z) + B(x, z, ε)u + F (x, z, ε))
= −(z 2 + x1 )
= 0.
(16)
The blow up map Φ : Sns +1 × [0, ∞) → Rns +2 is defined
by
x = r̄α x̄, z = r̄γ z̄, ε = r̄ρ ε̄,
(17)
Pns
where x̄ = (x̄1 , . . . , x̄ns ), i=1 x̄i + z̄ + ε̄ = 1, α ∈ Zns ,
(γ, ρ) ∈ Z2 and r̄ ∈ R≥0 . In principle, we could set
all the exponents in the blow up map to 1, but this
would require more than one coordinate transformation
to completely desingularize (16) [4, 15, 30]. A good choice
of weights is 6 α1 = · · · = αn = 2, γ = 1, and ρ = 3.
As described in Section 3, it is more convenient to work
with charts, which in this case are defined as κ±x̄i =
{x̄i = ±1}, κ±z̄ = {z̄ = ±1}, and κε̄ = {ε̄ = 1}, and its
induced local vector fields. To simplify the exposition,
let us adopt the following notation.

We shall design a feedback controller u = u(x, z, ε) that
stabilizes the origin of (13). Moreover, we consider the
slow-actuated (see Remark 8 below) and affine case,
namely
ẋ = f (x, z, ε) + B(x, z, ε)u(x, z, ε)
εż = g(x, z, ε) = −(z 2 + x1 ),

Geometric desingularization of a folded slowfast control system

(14)

where B = B(x, z, ε) ∈ Rns ×m depends smoothly on its
arguments. We make the following two assumptions.
5

In a more singularity theory language, g = 0 coincides with
the critical set of the Fold catastrophe [2, 17, 50]. Moreover,
the choice of the negative sign is just for convenience and a
completely equivalent analysis follows otherwise.

6

The choice of the weights depends on the quasihomogenity
type of the vector field, for details see [1, 15, 16, 31].

5

Notation: Let κ±w be a chart. If H(x, z, ε) is a function
in the original space ((x, z, ε) ∈ Rns +2 ) its blow up is
denoted by H̄, i.e., H̄ = H ◦ κ±w .

To avoid confusion of notation, the local coordinates in
κε̄ are denoted as (χ, ζ, ρ) = (χ1 , . . . , χn , ζ, ρ) ∈ Rns +2 .
Therefore the blow up map reads as
(x1 , . . . , xns , z, ε) = (ρ2 χ1 , . . . , ρ2 χns , ρζ, ρ3 ),

Remark 9 In the analysis to be performed below, the
most important chart is κε̄ . This is because it is precisely
in such a chart where the singular dependence on ε is
resolved. Next we provide the blown up vector field in the
chart κε̄ .

and then, from Proposition 10 X̃ ε̄ reads as
 0
ρ
ε̄
X̃ : χ0
 0
ζ

Proposition 10 The vector field X̃ ε̄ , which corresponds
to the blow up of (16) in the chart κε̄ , reads as
 0
r̄
X̃ ε̄ : x̄0
 0
z̄

=0
= A + r̄2 L1 x̄ + r̄L2 z̄ + B̄ ū + F̄
= −(z̄ 2 + x̄1 ),

=0
= A + ρ2 L1 χ + ρL2 ζ + B̄ ū + F̄
= −(ζ 2 + χ1 ),

(22)

where B̄ = B̄(ρ, χ, ζ) = B(ρ2 χ, ρζ, ρ3 ), and similarly
for F̄ and ū. Moreover F̄ |{ρ=0} = 0. Observe that Assumption A1 readily implies local controllability of χ0 =
A + ρ2 L1 χ + ρL2 ζ + B̄ ū + F̄ for all ρ ≥ 0 sufficiently
small.

(18)

where ū = ū(r̄, x̄, z̄) and F̄ = F̄ (r̄, x̄ z̄).

Remark 12 Note that in (22) ρ is a regular perturbation
parameter. Thus, the idea is to solve first the stabilization problem for ρ = 0 and then use regular perturbation
arguments [36] to guarantee the stability of the origin of
(22) for ρ ≥ 0 sufficiently small.

PROOF. The coordinates in the chart κε̄ are defined
by (x1 , . . . , xns , z, ε) = (r̄2 x̄1 , . . . , r̄2 x̄ns , r̄z̄, r̄3 ). From
ε = r̄3 , it follows that r̄0 = 0. Then, x0 = 2r̄r̄0 x̄ + r̄2 x̄0 ,
which implies

x̄0 = r̄ A + r̄2 L1 x̄ + r̄L2 z̄ + B̄ ū + F̄ ,
(19)

The main argument to relate the stability of the blown
up vector field X̃ ε̄ and X is the following.

where B̄ = B̄(r̄, x̄, z̄) = B(r̄2 x̄, r̄z̄, r̄3 ), and similarly for
ū and F̄ . In a similar way one obtains

z̄ 0 = −r̄ z̄ 2 + x̄1 .
(20)

Proposition 13 Consider a SFS X given by (16) and
its blow up version X̃ ε̄ in the chart κε̄ given by (18). If
for each ρ ∈ (0, ρ0 ] the point (χ, ζ) = (0, 0) is G.A.S. for
X̃ ε̄ , then for each ε ∈ (0, ρ1/3 ] the point (x, z) = (0, 0) is
G.A.S. for X .

Note, however, that the blown up vector fields (19)-(20)
degenerate at {r̄ = 0}. Thus, one defines the desingularized vector field by rescaling by a factor 1r̄ , which leads
to the result.
2

PROOF. Let ρ = ρ0 > 0 be fixed. Then the change of
coordinates is defined as

Remark 11 The local vector fields in the other charts
are obtained similarly to Proposition 10.

φ(χ1 , . . . , χns , ζ) = (ρ20 χ1 , . . . , ρ20 χns , ρ0 ζ)
= (x1 , . . . , xns , z).

The stability properties of the blown up vector field X̃ ε̄ of
(18) are carried over similar properties into the original
SFS. The main argument is that “stability should be
invariant under changes of coordinates”, as is analyzed
in the next section.
6

(21)

(23)

Note that φ is a diffeomorphism with positive definite Jacobian. On the other hand, the hypothesis that (χ, ζ) =
(0, 0) is G.A.S. for X̃ ε̄ implies that there exists a ρ-family
of Lyapunov functions V̄ρ (χ, ζ) satisfying
• V̄ρ > 0 for all (χ, ζ) ∈ Rns +1 \ {0},
• V̄ρ0 < 0 for all (χ, ζ) ∈ Rns +1 \ {0},
• V̄ρ is radially unbounded.

Controller design via geometric desingularization

In this section we show how to use geometric desingularization to stabilize a non-hyperbolic point of a SFCS.
The method is to design a controller in the blown up
space and then to blow such controller down. We do this
mainly in the central chart κε̄ as argued above. However, the directional charts may be used for other purposes, for example, to increase performance or to shape
the transitory behavior of the trajectories of a SFS, etc.

Define the ε-family of Lyapunov candidate functions
Vε = V̄ρ ◦ φ−1 , where ε = ρ3 . From the properties of φ
(23), namely that φ is a diffeomorphism with positive
definite Jacobian, it follows that Vε > 0 and Vε0 < 0 for
all (x, z) ∈ Rns +1 \ {0}. Therefore, Vε is an ε-family of
Lyapunov functions. Finally, let ||(χ, ζ)|| → ∞, which
clearly implies that ||(x, z)|| = ||φ(χ, ζ)|| → ∞. From
6

the definition of Vε we have Vε (x, z) = V̄ρ ◦ φ−1 (x, z) =
V̄ρ (ε−2/3 x, ε−1/3 z), which in fact shows that Vε is also
radially unbounded.
2

Then (26) is rewritten as

Proposition 13 means that we can design controllers to
stabilize a SFCS by designing it in the blown up space.
The way the controller ū is actually designed depends on
the specific context of the problem. Below we present a
particularly interesting case where even though the origin is non-hyperbolic and the fast variable is not actuated, one is able to inject a hyperbolic behavior through
actuating the slow variables.

Let us now propose the Lyapunov function V1 (Y, ζ) =
1 T
1
2
2 γζ + 2 Y Y , γ > 0, it follows that

6.1

with c1 > 0, and then V10 = −c0 ζ 2 − c1 Y T Y . Therefore, such choice of controller ensures that the origin
(Y, ζ) = (0, 0) ∈ Rn+1 is a GAS equilibrium point of
(28). Moreover, note that (ξ, ζ) → (0, 0) imply χ1 → 0,
which shows the claim. Finally, for the above algorithm
to work, it is convenient to ensure that ξ = χ1 − α → 0
fast enough. This is because it is precisely at ξ = 0 when
the dynamics of ζ become hyperbolic. For the latter to
happen, let us consider the closed-loop matrix associated to (28) 7 , which reads as

Y 0 = A − W + B̄0 ū0
ζ 0 = −(c0 ζ + Y1 ).


V10 = −c0 ζ 2 + Y T A − W − γζ ê1 + B̄0 ū0 .

B̄0 ū0 = −c1 Y − A + W + γζ ê1

In this section we design a controller which induces a
hyperbolic behavior in (22) around the origin. We do
this via the backstepping algorithm [44].
Proposition 14 Consider (22). If ū is designed such
that
(24)

where B̄0 = B̄|{ρ=0} , and

Φ=

β̄(χ, ζ) = (2ζ − c0 )(ζ 2 + χ1 ) − c1 ζ 2 ,

(29)

Thus we design ū0 so that

Hyperbolicity injection

B̄0 ū = −A − c1 χ + (c0 c1 + γ)ζ ê1 + β(χ, ζ)ê1 ,

(28)

(25)

"
#
−c1 In γê1
−êT1 −c0

∈ R(ns +1)×(ns +1) .

(30)

(31)

It is straightforward to show that the eigenvalues of Φ
are




(32)
−c1 , . . . , −c1 , λ1,2 ,
{z
}

|

with 0 < γ  c1 and 0 < c0  c1 , then the origin of
(22) is rendered globally asymptotically stable for ρ ≥ 0
sufficiently small.

ns −1

where λ1,2 6= −c1 are the solutions of the quadratic polynomial s2 +(c0 +c1 )s+γ+c0 c1 = 0. Both eigenvalues λ1,2
have negative real part, so the closed-loop system (28)
is GAS, however, as argued above, it is convenient that
ξ = χ1 − α → 0 fast enough so that the backstepping algorithm quickly induces the desired hyperbolic behavior
in ζ. To ensure this, note that
eigenval if γ is small, the 




ues of Φ are approximately −c1 , . . . , −c1 , −c0 , which
{z
}


|


PROOF. First, we restrict the analysis to {ρ = 0}.
Therefore we have that (22) reads as
X̃ ε̄ |ρ=0 = A + B̄0 ū0

 ∂
∂
− (ζ 2 + χ1 ) ,
∂χ
∂ζ

(26)

where ū0 denotes the restriction ū|ρ=0 . Next, consider
the system
ζ 0 = −(ζ 2 + χ1 ).
(27)

ns times

The idea is to first treat χ1 in (27) as a virtual controller
which injects a hyperbolic behavior on ζ in (27). For
this let χ1 = α(ζ) and consider the Lyapunov function
V0 (ζ) = 12 ζ 2 . It follows that V00 = ζζ 0 = −ζ(ζ 2 + α).
Define α(ζ) = −ζ 2 + c0 ζ with c0 > 0. Then, by Lyapunov arguments it follows that ζ = 0 is a global asymptotic equilibrium point of the (virtual) closed-loop system ζ 0 = −(ζ 2 + α) = −c0 ζ. Note that, indeed, ζ = 0 is
a hyperbolic equilibrium point of ζ 0 = −(ζ 2 + α).

motivates the condition 0 < c0  c1 . The expressions
(24)-(25) are obtained by returning from (30) to the
coordinates (χ, ζ). Since (22) is a regular perturbation
problem in ρ, and the origin is still an isolated equilibrium point of (22) for ρ > 0, all the above arguments
hold for ρ > 0 sufficiently small [36].
2
The final step is to obtain the controller induced by
Proposition 14.

To simplify notation let ξ = χ1 −α, Y = (Y1 , . . . , Yns ) =
(ξ, χ2 , . . . , χns ) ∈ Rns and W = (2ζ − c0 )(ξ + c0 ζ)ê1 .

7

7

In fact, the closed-loop system (28) is linear.

proposed to regularize Σ1 by adding a parasitic capacitance in parallel to R as depicted in Figure 2.(b) (electric circuit Σ2 ), see Example 6 of [46] for more details,
and [13]. The capacitance Cε is assumed to be small,
e.g. Cε = ε. Thus, as ε → 0, the behavior of of circuit
Σ2 approaches that of Σ1 . Let (x1 , x2 , z) = (IL , VC , Vε ),
where Vε denotes the voltage at the capacitor Cε . The
equations describing the behavior of the circuit Σ2 read
as

Theorem 15 Consider the SFCS
 0
x
X : z0
 0
ε

= ε (f (x, z, ε) + Bu(x, z, ε))
= −(z 2 + x1 )
= 0,

(33)

If the control law u is designed such that
B0 u = −A −

(c0 c1 + γ)
c1
x+
zê1 + β(x, z)ê1 (34)
2/3
ε
ε1/3

Lẋ1 = −z − x2
C ẋ2 = x1
1
εż = − z 3 + z + x1 .
3

where A = f (0, 0, 0), B0 = B(0, 0, 0) and


β = ε−2/3 −c1 z 2 + (2ε−1/3 z − c0 )(z 2 + x1 ) ,

(37)

(35)
We immediately note that the corresponding critical
manifold is actually given by the characteristic equation
of the nonlinear resistor, namely

with 0 < c0  c1 and γ > 0, then the origin (x, z) =
(0, 0) ∈ Rn+1 is rendered globally asymptotically stable
for ε > 0 sufficiently small.


S=
PROOF. The proof follows from Proposition 14, the
corresponding blow down of the controller ū, and Proposition 13.
2

(x1 , x2 , z) ∈ R3 | x1 =


1 3
z −z .
3

Moreover, it is straightforward to see that the re-

Remark 16 In Theorem 15 we only require the knowledge of A = f (0, 0, 0), while the values of L and F are
not important as long as they are bounded. Theorem 15
can be extended to an adaptive, but local, version, where
the value of A may be assumed unknown, see [20].

Vε

R

Iε

IL
C

L

R
IR

IL

7

(38)

(a)

Example: Stabilization at a fold point of a nonlinear electric circuit

L
I2

Cε

Let us consider the electric circuit Σ1 as shown in Figure
2.(a), where the Capacitor C and Inductor L are usual
elements, but the Resistor R is assumed to be nonlinear,
that is IR = f (VR ) where f is some nonlinear function
and IR and VR denote the current and voltage of R respectively, see Example 5. of [46].

I1
Iε

u1

R2

R1

C

R

IC
IR

IL

+
−

C

L

For the purpose of this example we shall assume that
f (VR ) =

1 3
V − VR ,
3 R

(b)

(36)

(c)

Fig. 2. (a) Electric circuit Σ1 with a nonlinear resistive load.
(b) The regularization of circuit (a), denoted by Σ2 . (c) The
controlled circuit Σ3 .

and then the characteristic curve of the nonlinear Resistor R is as depicted in Figure 3. We remark that the chosen type of nonlinear behavior is qualitatively the same
as the one in [46], and also appears in several other nonlinear elements, such as tunnel diodes [39].

f (VR )
p2

The characteristic curve depicted in Figure 3 has two
fold points p1,2 located at (VR , IR ) = (±1, ∓ 23 ). At these
points the differential equation describing the behavior
of the circuit becomes singular. To overcome this it is

VR
p1

Fig. 3. The characteristic curve of the nonlinear Resistor R.

8

u2

Let us choose, for example, the operating point P =
(x1 , x2 , z) = (− 23 , 0, 1). We remark that since P is a nonhyperbolic point, the classical theory e.g. [25] cannot be
used, and instead we shall employ the theory developed
in this paper. The first step is to define new coordinates
(X1 , X2 , Z) = (−x1 − 23 , x2 , z −1), in this way (39) reads
as

gion of S between p1 and p2 is repelling while the
rest of S is attracting. Furthermore, one can easily show that (37) has three
equilibrium
points
√ √
√
√
{q1 , q2 , q3 } =
(0, 0, 0), (0, − 3, 3), (0, 3, − 3) ,
where q1 is unstable, while q2 and q3 are stable. With
this information we can qualitatively describe the dynamics of (37) as follows: for initial conditions away
from S, the trajectories of (37) quickly approach a stable region of S, and then evolve around it. Here two
things may happen, trajectories may converge to an
equilibrium point q2 or q3 , or they can approach a fold
point. When a trajectory reaches a fold point, then the
trajectory jumps towards a stable region of S and then
follows the same behavior as as described before. A
sample of trajectories of (37) is provided in Figure 4.

2
LẊ1 = 1 − R1 − R1 X1 + X2 + Z − R1 u1
3
2
X2
u2
C Ẋ2 = − − X1 +
−
3
R2
R2
1
εŻ = −(Z 2 + X1 ) − Z 3 .
3

Note that up to leading order terms, (40) is of the form
studied above. After the change of coordinates the operating point is located at the origin of (40). According to
Theorem 15 the controller that stabilizes P is given by


2
c1
c0 c1 + γ
1
L
R1 − 2/3 X1 +
Z +β
− +
u1 = −
R1
L 3L
ε
ε1/3


c1
2
u2 = −CR2
− 2/3 X2 ,
3C
ε
(41)

S
q2
q1
q3

z

(40)

x1

x2


where β = ε−2/3 −c1 Z 2 + (2ε−1/3 Z − c0 )(Z 2 + X1 ) .
To witness the effects of the controller, let us choose
parameters: C = 1F, L = 1H, R1 = R2 = 1Ω, ε = 0.05,
and controller gains c0 = γ = 1 and c1 = 10. Next we
choose initial conditions near the limit cycle, to allow
the system oscillate. For simulation purposes we let the
system evolve 10 seconds in open-loop. Then, at t =
10s we activate the controller. The results are shown in
Figure 5, where we see that when the controller takes
action, the trajectories quickly converge to the operating
point P .

Fig. 4. Trajectories of the open-loop system (37). The surface
is the critical manifold S while dashed lines represent lines
of folds, and the dots stand for the equilibrium points. Note
the existence of a limit cycle, which can be conjectured from
the stability of q2 and q3 . However such analysis falls off the
scope of this document and shall not be discussed further.

Our goal, however, is to set as operating point a fold
point. In order to do this we introduce controllers as depicted in Figure 2.(c) leading to circuit Σ3 , where we add
a current source u1 in parallel with the resistor R1 , all
these in series with the inductor L. The capacitor voltage x2 is modulated by placing a voltage source u2 and
a resistor R2 in series, all these in parallel with the capacitor C. Both ”virtual” resistors R1 and R2 are understood to be part of the controller, and in physical terms
are for proper coupling of the inputs and may also serve
as tuning parameters in a real implementation.

8

In this paper we have introduced the geometric desingularization technique to control systems. The main contribution is a controller design method for non-hyperbolic
fold points of slow-fast systems. With this we are able
to deal with singular perturbation problems for which
the classical theory does not apply. An essential feature
of our contribution is that the controller only actuates
on the slow variables, making it more suitable for applications. As a case study, we have provided a controller
based on the backstepping algorithm that renders the
origin of a SFCS globally asymptotically stable.

The behavior of the circuit Σ3 is described by
Lẋ1 = −z − x2 − R1 x1 + R1 u1
u2
x2
−
C ẋ2 = x1 +
R2
R2
1 3
εż = − z + z + x1 .
3

Conclusions

Further research directions in view of the potential applications include: the assumption that the slow system
is under-actuated, output feedback control, trajectory
and path-following along sets of non-hyperbolic points,

(39)

9

X1 [A]

2
0

[8]

−2

[9]
X2 [V ]

0

[10]
−1

Z [V ]

3

[11]

0

−4

0

5

10

15

[12]

20

t [s]

Fig. 5. Simulation results of the controlled circuit Σ3 depicted
in Figure 4.(c) and with the controller of Theorem 15.

[13]

and the extension of the theory to a general class of
SFCS with one fast variable but with arbitrarily degenerate non-hyperbolic points. The latter class of systems
are of the form

[14]

ẋ = f (x, z, ε) + B(x, z, ε)u
!
k−1
X
k
i−1
εż = − z +
xi z
+ O(ε),

[15]

(42)

[16]

i=1

where k ≥ 2. Another challenging framework is to consider non-hyperbolic slow-fast systems over networks.

[17]

References
[18]
[1]

[2]
[3]
[4]
[5]
[6]
[7]

V. I. Arnold, V. V. Goryunov, O. V. Lyashko, A. Iacob, and V. A. Vasil’ev. Singularity Theory: I. Number VI in Encyclopaedia of mathematical sciences.
Springer Berlin Heidelberg, 1998.
V. I. Arnold, S. M. Gusein-Zade, and A. N.
Varchenko. Singularities of Differentiable Maps,
Volume I, volume 17. Birkhäuser, 1985.
T. Bröcker. Differentiable Germs and Catastrophes,
volume 17 of Lecture Note Series. Cambridge University Press, 1975.
H. W. Broer, F. Dumortier, S. J. van Strien, and
F. Takens. Structures in Dynamics. North-Holland,
1991.
L. O. Chua and A. Deng. Impasse points. Part I:
Numerical aspects. International Journal of Circuit
Theory and Applications, 17(2):213–235, 1989.
L. O. Chua and A. Deng. Impasse points. Part II:
Analytical aspects. International Journal of Circuit
Theory and Applications, 17(3):271–282, 1989.
M. Desroches, J. Guckenheimer, K. Krauskopf,
C. Kuehn, H. M. Osinga, and M. Wechselberger.

[19]

[20]

[21]
[22]

[23]

10

Mixed-mode oscillations with multiple time scales.
SIAM Review, 54(2):211–288, 2012.
F. Dumortier and R. Roussarie. Canard Cycles and
Center Manifolds, volume 121. American Mathematical Society, 1996.
N. Fenichel. Geometric singular perturbation theory for ordinary differential equations. Journal of
Differential Equations, 31(1):53–98, 1 1979.
Z. Gajic and M. Lelic. Improvement of system order
reduction via balancing using the method of singular perturbations. Automatica, 37(11):1859 – 1865,
2001.
M. Galtier and G. Wainrib. Multiscale analysis of
slow-fast neuronal learning models with noise. The
Journal of Mathematical Neuroscience, 2(1):1–64,
2012.
X. Han, E. Fridman, and S. K. Spurgeon. Sliding mode control in the presence of input delay:
A singular perturbation approach. Automatica,
48(8):1904 – 1912, 2012.
E. Ihrig. The regularization of nonlinear electrical
circuits. Proceedings of the American Mathematical
Society, pages 179–183, 1975.
A. Isidori. Nonlinear Control Systems. Communications and Control Engineering. Springer-Verlag
London, 1995.
H. Jardón-Kojakhmetov. Classification of constrained differential equations embedded in the theory of slow-fast systems. PhD Thesis, University of
Groningen, 2015.
H. Jardón-Kojakhmetov. Formal normal form of Ak
slow–fast systems. Comptes Rendus Mathematique,
353(9):795–800, 2015.
H. Jardón-Kojakhmetov and H. W. Broer. Polynomial normal forms of constrained differential equations with three parameters. Journal of Differential
Equations, 257(4):1012–1055, 2014.
H. Jardón-Kojakhmetov, H. W. Broer, and
R. Roussarie. Analysis of a slow-fast system near a
cusp singularity. Journal of Differential Equations,
2016.
H. Jardón-Kojakhmetov and J. M. A. Scherpen.
Stabilization of a planar slow-fast system at a nonhyperbolic point. In Proceedings of the 22nd International Symposium on Mathematical Theory of
Networks and Systems, July 2016.
H. Jardón-Kojakhmetov, J. M. A. Scherpen, and
D. del Puerto-Flores. Nonlinear adaptive stabilization of a class of planar slow-fast systems at a nonhyperbolic point. In Proceedings of the American
Control Conference, 2017.
C. K. R. T. Jones. Geometric singular perturbation
theory. In Dynamical Systems, LNM 1609, pages
44–120. Springer-Verlag, 1995.
T. J. Kaper. An Introduction to Geometric Methods and Dynamical Systems Theory for Singular
Peturbation Problems. In Symposia in Applied
Mathematics, volume 56, pages 85–131. AMS, 1999.
P. V. Kokotovic. Applications of Singular Pertur-

[24]

[25]

[26]

[27]

[28]
[29]
[30]
[31]

[32]

[33]
[34]

[35]
[36]
[37]

[38]
[39]

[40] A. Roberts, J. Guckenheimer, E. Widiasih, A. Timmermann, and C. K. R. T. Jones. Mixed-mode oscillations of el niño–southern oscillation. Journal of
the Atmospheric Sciences, 73(4):1755–1766, 2016.
[41] H. G. Rotstein. Mixed-Mode Oscillations in Single
Neurons, pages 1–9. Springer New York, New York,
NY, 2013.
[42] V. R. Saksena, J. O’Reilly, and P. V. Kokotovic.
Singular perturbations and time-scale methods in
control theory: Survey 1976-1983. Automatica,
20(3):273–293, 1984.
[43] R. G. Sanfelice and A. R. Teel. On singular perturbations due to fast actuators in hybrid control
systems. Automatica, 47(4):692 – 701, 2011.
[44] S. Sastry. Nonlinear Systems: Analysis, Stability,
and Control. Springer-Verlag New York, 1999.
[45] A. Shilnikov. Complete dynamical analysis of a
neuron model. Nonlinear Dynamics, 68(3):305–328,
2012.
[46] S. Smale. On the mathematical foundations of electrical circuit theory. Journal of Differential Geometry, 7(1-2):193–210, 1972.
[47] M. W. Spong. Modeling and Control of Elastic
Joint Robots. Journal of Dynamic Systems, Measurement, and Control, 109(4):310, 1987.
[48] P. Szmolyan and M. Wechselberger. Canards in R3 .
Journal of Differential Equations, 177(2):419–453,
2001.
[49] P. Szmolyan and M. Wechselberger. Relaxation oscillations in R3 . Journal of Differential Equations,
200(1):69–104, 6 2004.
[50] F. Takens. Constrained equations; a study of implicit differential equations and their discontinuous
solutions. In Structural stability, the theory of catastrophes, and applications in the sciences, pages 143–
234. Springer, 1976.
[51] Y. Tang, C. Prieur, and A. Girard. Singular perturbation approximation of linear hyperbolic systems
of balance laws. IEEE Transactions on Automatic
Control, 61(10):3031–3037, Oct 2016.
[52] A. R. Teel, L. Moreau, and D. Nesic. A unified
framework for input-to-state stability in systems
with two time scales. IEEE Transactions on Automatic Control, 48(9):1526–1544, Sept 2003.
[53] B. van der Pol and J. van der Mark. The heartbeat
considered as a relaxation oscillation, and an electrical model of the heart. The London, Edinburgh,
and Dublin Philosophical Magazine and Journal of
Science, Ser.7,6:763–775, 1928.
[54] D. Del Vecchio and J. J. E. Slotine. A contraction theory approach to singularly perturbed systems. IEEE Transactions on Automatic Control,
58(3):752–757, March 2013.
[55] L. Zhou, C. Yang, and W. Zhang. Control for a class
of non-linear singularly perturbed systems subject
to actuator saturation. IET Control Theory Applications, 7(10):1415–1421, July 2013.

bation Techniques to Control Problems. SIAM Review, 26(4):501–550, 1984.
P. V. Kokotovic, R. E. O’Malley, and P. Sannuti.
Singular perturbations and order reduction in control theory An overview. Automatica, 12(2):123–
132, 1976.
P. V. Kokotovic, J. O’Reilly, and H. K. Khalil.
Singular Perturbation Methods in Control: Analysis and Design. Academic Press, Inc., Orlando, FL,
USA, 1986.
I. Kosiuk and P. Szmolyan. Scaling in Singular Perturbation Problems: Blowing Up a Relaxation Oscillator. J. Applied Dynamical Systems, 10(4):1307–
1343, 2011.
M. Krupa and P. Szmolyan. Extending geometric singular perturbation theory to non hyperbolic
points: fold and canard points in two dimensions.
SIAM J. Math. Anal., 33:286–314, 2001.
M. Krupa and P. Szmolyan. Relaxation oscillation
and canard explosion. J. Diff. Eqns., 174:312–368,
2001.
M. Krupa and M. Wechselberger. Local analysis
near a folded saddle-node singularity. Journal of
Differential Equations, 248(12):2841–2888, 2010.
C. Kuehn. Multiple Time Scale Dynamics. Springer
International Publishing, 2015.
E. Lombardi and L. Stolovitch. Normal forms of
analytic perturbations of quasihomogeneous vector
fields: rigidity, invariant analytic sets and exponentially small approximation. Ann. Sci. Éc. Norm.
Supér., 43(4), 2010.
L. Marconi, L. Praly, and A. Isidori. Robust asymptotic stabilization of nonlinear systems with nonhyperbolic zero dynamics. IEEE Transactions on
Automatic Control, 55(4):907–921, April 2010.
R. Marino and P. V. Kokotovic. A geometric approach to nonlinear singularly perturbed control
systems. Automatica, 24(1):31–41, 1 1988.
W. Marszalek and Z. W. Trzaska. Singularityinduced bifurcations in electrical power systems.
IEEE Transactions on Power Systems, 20(1):312–
320, Feb 2005.
L. Menini and A. Tornamb. Stability analysis of planar systems with nilpotent (non-zero) linear part.
Automatica, 46(3):537 – 542, 2010.
J. A. Murdock. Perturbations: theory and methods,
volume 27. SIAM, 1999.
E Neumann and A Pikovsky. Slow-fast dynamics in josephson junctions. The European Physical Journal B-Condensed Matter and Complex Systems, 34(3):293–303, 2003.
Y. Pan and H. Yu. Dynamic surface control via
singular perturbation analysis. Automatica, 57:29–
33, 7 2015.
G. Reissig. Differential-algebraic equations and impasse points. IEEE Transactions on Circuits and
Systems I: Fundamental Theory and Applications,
43(2):122–133, Feb 1996.

11

