Transverse Contraction Criteria for Stability of Nonlinear Hybrid Limit Cycles  Justin Z. Tang and Ian R. Manchester  Abstract — In this paper, we derive differential conditions guaranteeing   the   orbital   stability   of   nonlinear   hybrid   limit cycles. These conditions are represented as a series of point- wise linear matrix inequalities (LMI), enabling the search for stability certificates via convex optimization tools such as sum- of-squares programming. Unlike traditional Lyapunov-based methods, the transverse contraction framework developed in this paper enables proof of stability for hybrid systems, without prior knowledge of the exact location of the stable limit cycle in state space. This methodology is illustrated on a dynamic walking example.  I. I NTRODUCTION  Nonlinear hybrid dynamical systems with periodic solu- tions are widely found in diverse engineering and scientific fields   such   as   electronics   and   mechanics.   These   hybrid systems contain continuous-time and discrete-time dynamics which interact with each other. Stability of these dynamical systems is often a fundamental requirement for their practical value in applications. In this paper, we address the question: do all solutions of a hybrid nonlinear system starting in a particular set   K  converge to a stable unique limit cycle? A major motivation of this work is the study of underac- tuated bipedal locomotion [1], which can be represented as limit cycles in the state space [2]. The control design and sta- bility analysis of these “dynamic walkers” are difficult since their dynamics are inherently hybrid and highly nonlinear [3], [4]. The   most   well-known   stability   analysis   tool   for   limit cycles is the Poincar ́ e map [5], which describes the re- peated passes of the system through a single transversal hypersurface. However, for nonlinear systems, the Poincar ́ e map generally cannot be found explicitly. Further, since the system’s evolution is only analyzed on a single surface, regions of stability in the full state space are difficult to evaluate. In practice, stability in the full state space is often esti- mated using exhaustive simulation, such as via cell-to-cell mapping [6], which has been applied to analysis of walking robots [7]. However, computational costs of these methods are exponential in the dimension of the system. In recent years, convex optimization methods have been widely applied in search for a “stability certificate” based on Lyapunov theory [8]. To characterize regions of stability  *This work was supported in part by the Australian Research Council. The   authors   are   with   the   Australian   Centre   for   Field   Robotics (ACFR),   Department   of   Aerospace,   Mechanical   and   Mechatronic   En- gineering,   University   of   Sydney,   NSW   2006,   Australia   { j.tang i.manchester } @acfr.usyd.edu.au  for limit cycles, [9] and [10] introduced the notion of the Surface Lyapunov function, which verifies stability based on the “impact map” between one switching surfaces to the next switching surface. The method is limited to Piecewise Linear Systems. In [11], nonlinear limit cycle stability analysis was performed by constructing Lyapunov functions in the transverse dynamics. However, these Lyapunov based methods require knowl- edge of the exact location of the limit cycle in state space, and hence are not applicable when the system dynamics are uncertain, since uncertainty will generally change the location of the limit cycle. An alternative approach to Lyapunov methods is to search for a contraction metric [12], [13]. By defining stability incrementally between two arbitrary nearby trajectories, con- traction analysis answers the question of whether the limiting behaviour of a given dynamical system is independent of its initial conditions. For analysis of limit cycles, transverse contraction was first introduced in [14]. In this paper, we propose a transverse contraction frame- work for analysis of hybrid limit cycles, building on the work of transversal surface construction in [11], and continuous transverse contraction of [14]. For the purposes of robustness analysis, an important advantage is that Lyapunov functions must be generally constructed around a known equilibrium, whereas a contraction metric derived herein implies existence of a stable equilibrium indirectly. This is vital if the equi- librium point may change location depending on unknown dynamics. This paper proceeds as follows. Problem formulation and preliminary notations are outlined in Section II. In Section III, the transverse contraction conditions guaranteeing stabil- ity of a limit cycle in a nonlinear hybrid system is presented. We then formulate convex criteria enforcing these conditions on the nonlinear system in Section IV, thereby enabling the search for stability certificates via convex optimisation tech- niques such as sum-of-squares programming. An illustrative example is given with a dynamic walking model in Section V. Concluding remarks are given in Section VI. II. P RELIMINARIES AND   P ROBLEM   F ORMULATION  We consider the following class of autonomous hybrid dynamical problems.  Problem 1:   Consider a hybrid system with impulse.   ̇ x   =   f   ( x )   x / ∈   S −  i   (1)  x +   =   g ( x )   x   ∈   S −  i   (2)  arXiv:1403.5374v1 [math.OC] 21 Mar 2014
where   f   ,   g   are smooth, and   S i   for   i   = 1 ,   2 ...   is defined as a “switching surface.” We assume that   x   ∈   R n   and   f   :   K   →   R n , where a set   K  is a compact subset of   R n   and strictly forward invariant under  f   , such that any solutions of (1), (2) starting in   x (0)   ∈   K  is in the interior of   K   for all   t >   0 . We will refer to the Jacobian of   f   as   A ( x ) :=   ∂f ∂x   . We denote the solution curve of the system as   Φ( x 0 , t ) , such that   x ( t ) = Φ( x 0 , t )   is the solution at time   t >   0   of the dynamical system with initial state   x (0) =   x 0 . Suppose the system exhibits a non-trivial   T   -periodic orbit, i.e., for a periodic solution   x ∗ , there exists some   T >   0   such that   x ? ( t ) =   x ? ( t   +   T   )   for all   t . Such a solution cannot be asymptotically stable, as perturbations in phase are persistent. Instead,   orbital stability   is better posed [15]. The orbit of a periodic solution is the set   X   ?   :=   { x   ∈  R n   :   ∃ t   ∈   [0 , T   ) :   x   =   x ? ( t ) } . The solution is said to be   orbitally stable   if there exists a   b >   0   such that for any  x (0)   satisfying dist ( x (0) ,   X   ? )   < b , the unique solution exists and dist (Φ( x 0 , t ) ,   X   ? )   →   0   as   t   → ∞ . Further, the system is   exponentially orbitally stable   if the system is orbitally stable and there exists a   b >   0 , λ >   0 , k >   0   such that dist (Φ( x 0 , t ) ,   X   ? )   ≤   k dist (Φ( x 0 , t ) ,   X   ? ) e − λt . A switching surface defined by   S   :=   { x   |   c ( x ) = 0 }   is a   ( n   −   1) -dimensional hyperplane embedded in a manifold  M , where   c ( x )   is linear in   x . The tangent space of   M   at  x   ∈ M   is denoted as   T x M , and the tangent bundle of   M  is denoted as   T   M   =   ⋃  x ∈M { x } ×   T x M .  For simplicity, we assume that the region   K   is broken into a finite sequence of continuous “tubes” in the state space,   which   contains   no   equilibrium   point,   i.e.,   ∀ x   ∈  K, f   ( x ) T   f   ( x )   >   0 . Further, assume that all solutions start- ing in each continuous phase of   K   approach a particular switching surface and that   ∀ x   ∈   S i , f   ( x ) T   z i ( x )   6   = 0   where  z i ( x )   is the normal vector of   S i . Our goal is to verify that, in our hybrid system (1), (2), all solutions starting in the particular region   K   are orbitally stable and converge to a unique limit cycle. This is illustrated in Fig 1 for a system with two continuous phases and two switching surfaces. The verified region is shaded green, and the stable limit cycle in red. Continuous dynamics are shown in solid line with discrete impulse between switching surfaces shown in dotted line. For completeness, we now restate the transverse con- traction condition for continuous systems derived in [14]. Unless otherwise stated, we assume a Riemannian distance  V   ( x, δ x ) =   √ δ ′  x M   ( x ) δ x , with   M   ( x )   symmetric positive- definite for all   x .  Definition 1 (Transverse Contraction):   A continuous sys- tem    ̇ x   =   f   ( x )   is transverse contracting with rate   λ   if there exists a Riemannian metric   V   ( x, δ x )   satisfying  ∂V   ( x, δ x )  ∂x   f   ( x ) +   ∂V   ( x, δ x )  ∂δ x  ∂f   ( x )  ∂x   δ x   ≤ − λV   ( x, δ x )   (3) for all   δ x   6   = 0   such that   ∂V ∂δ x   f   ( x ) = 0 . S   1 S   2  Fig. 1.   Region of stability around a hybrid limit cycle  The latter condition requires   δ x   to be transverse to the flow of the system, i.e.,   δ x   and   f   ( x )   are orthogonal with respect to the metric   M   ( x ) . In the case that   V   ( x, δ x ) :=   √ δ ′  x M   ( x ) δ x , it is true if   δ ′  x M   ( x ) f   ( x ) = 0 . We now state a convex condition that is necessary and sufficient for transverse contraction derived in [14].  Definition 2 (Convex Criterion for Transverse Contraction):  A system    ̇ x   =   f   ( x )   is transverse contracting with rate   λ  and a metric   V   ( x, δ ) =   √ δ ′ M   ( x ) δ   if and only if there exists a function   W   ( x ) :=   M   ( x ) − 1   and   ρ ( x )   ≥   0   such that  W   ( x ) A ( x ) ′   + A ( x ) W   ( x ) −    ̇ W   ( x )+2 λW   ( x ) − ρ ( x ) Q ( x )   ≤   0  (4) where   Q ( x ) :=   f   ( x ) f   ( x ) ′  It   was   shown   in   [14]   that   if   a   system   is   transverse contracting, then all solutions starting with   x (0)   ∈   K   are stable under time reparameterization, or “Zhukovski stable” [16], and hence converge to a unique limit cycle. III.   C ONTRACTION   C ONDITIONS FOR   L IMIT   C YCLES IN  H YBRID   S YSTEMS  In this section, we derive the transverse contraction con- ditions for hybrid nonlinear limit cycles and show that these conditions guarantee the stability of a unique limit cycle within a particular set   K . For simplicity of expression we will consider the problem of single switching surface and a single set of continuous dynamics. However, extension to multiple switches and mul- tiple continuous phases is straightforward.  Condition 1 (Metric condition):   For a nonlinear system   ̇ x   =   f   ( x )   with a flat switching surface   S   :=   { x   |   c ( x ) = 0 }  if for all   x   ∈   S , the tangent space of the switching surface at   x ,   T x M   satisfies  f   ( x ) T   M   ( x ) δ x   = 0   (5) for all   δ x   ∈   T x M , then, according to the metric   M   ( x ) , all trajectories approach the switching surface orthogonally. This is true since the direction vector of for any trajectory on the switching surface   x   ∈   S   is given by   f   ( x ) .  Theorem 2:   If, firstly, the continuous part of the system in Eq. (1) is transverse contracting according to Definition 1;
secondly, Condition 1 is satisfied on the switching surface; and, thirdly, the discrete part of the system in Eq. (2) satisfies  ∂g ∂x  ′  M ∂g ∂x   −   M   ≤   0   (6) for all   δ ′  x M f   =   0 ; then, all solutions in   K   are locally Zhukovski stable. Hence, we can construct, for each   x   ∈   K  locally a ball,   B x , of constant radius centred around a given trajectory at   x , where trajectories starting in   B x   remain in  B x   under time reparameterization as   t   → ∞ .  Proof:   From the first condition in the Theorem and by the results of Definition 1, we know that the system is transverse contracting – i.e., all virtual displacements   δ x   in the subspace defined by the plane   δ ′  x M   ( x ) f   ( x ) = 0   are contracting. Suppose we define locally for each   x   ∈   K   a smooth change of coordinates   x   →   ( τ, x ⊥ ) , highlighting the dy- namics tangential and transverse to the flow of the system. The transformation,   Π( x ) , can be represented by a set of bases   { e 1 , ..., e n } , where   e 1   is in the direction of   f   ( x )   and  e 2 , ..., e n   are independent and lie on the contracting plane defined by   δ ′  x M   ( x ) f   ( x ) = 0 . By definition of orthogonality, this also implies   e 2 , ...e n   are orthogonal to   e 1 . This transfor- mation is given by:  Π( x ) x   =        e 1  e 2  . . .  e n         x   =  [   τ x ⊥  ]  (7) where   τ   is a 1-dimensional “phase” variable tangential to the flow of the system, and   x ⊥   is the   ( n   −   1) -dimensional transverse dynamics in the contracting subspace. The corresponding virtual displacements of the system in the new coordinate can be found via the Jacobian of the transformation.   [ δ τ  δ ⊥  ]  =  ̄ Θ δ x   (8) where    ̄ Θ   :=   ∂ Π  ∂x   .   δ ⊥   lies   in   the   subspace   defined   by  δ x M   ( x ) f   ( x ) = 0 . It is shown by [17], [11] that the differential system in the new coordinate becomes  d dt  [ δ τ  δ ⊥  ]  =  [ 0   ?  0   A ⊥ ( x )  ] [ δ τ  δ ⊥  ]  (9) where   A ⊥   is the transverse linearization. By the construction of   Π( x )   in (7), since   e 1   is orthogonal to all   e 2 , ..., e n , we have   e ′  1 M   ( x ) e k   = 0   for all   k   = 2 , ..., n . Therefore, we can separate the transverse and tangential components in the metric   M   ( x )   as below:  V   =   δ ′  x M   ( x ) δ x   (10)  V   =  [ δ τ  δ ⊥  ] T   [ M τ   ( x )   0 0   M ⊥ ( x )  ] [ δ τ  δ ⊥  ]  (11)  V   =   M τ   ( x ) | δ τ   | 2   +   δ ⊥ M ⊥ ( x ) δ ⊥   (12) Since   δ ⊥   lies   entirely   on   the   plane   defined   by  δ x M   ( x ) f   ( x ) = 0 , which is contracting by Definition 1, we yield:  δ ⊥  (  A ⊥ ( x ) ′ M ⊥ ( x ) +   M ⊥ ( x ) A ⊥ ( x ) +    ̇ M ⊥ ( x )  )  δ ⊥   <   0  (13) By the results of [17], we can construct a Lyapunov function   V ⊥   =   x ′ ⊥ M ⊥ x ⊥   and for sufficiently small   ‖ x ⊥ ‖  around the coordinate change,   d dt   V ⊥   <   0   can be guaranteed by (13). Now, by the second condition in the Theorem, during discrete instantaneous switching when   x   ∈   S , all trajectories approach the switching surface orthogonally by the results of Condition 1. Hence, the transversal plane   δ x M   ( x ) f   ( x ) = 0  for trajectories on the switching surface aligns with the switching   surface   itself.   Suppose   now   on   the   switching surface,   V   is non-increasing during the impulse:  δ ′ +  x   M   ( x ) δ +  x   ≤   δ ′−  x   M   ( x ) δ −  x   (14)  δ ′  x  (   ∂g ∂x  ′  M ∂g ∂x   −   M  )  δ x   ≤   0   (15) for all   δ x M   ( x ) f   ( x ) = 0 . Since the transversal plane for the trajectory coincide with the switching surface, by construction in (7),   δ ⊥   lies entirely on the switching surface. Hence (15) is equivalent to   V   +  ⊥   ≤  V   − ⊥   . We now prove local Zhukovski stability using construction similar to [17]. Let   σ   be a nearby solution to   x . Then, by the construction of   x ⊥ , there exists some compact region around   x   where   k 1 dist ( σ, x )   ≤ ‖ x ⊥ ‖ ≤   k 2 dist ( σ, x ) , for some   k 1 , k 2   >   0 . Hence,   k 3 dist ( σ, x )   ≤   V ⊥   ≤   k 4 dist ( σ, x )  for some   k 3 , k 4   >   0 . Since our conditions show   V ⊥   is uniformly decreasing for all   x   ∈   K   in both the continuous dynamics and across the switching surface, as   t   → ∞ ,   x ⊥   →   0   locally. Therefore, for all   x   ∈   K , there exists locally a ball,  B x , of constant radius centred around the trajectory at   x , where trajectories starting in   B x   remain in   B x   under time reparameterization as   t   → ∞ .  Theorem 3:   If the conditions of Theorem 2 is satisfied, for every pair of solutions   x 1   and   x 2   in   K , there exists time reparameterization   τ   ( t )   such that   x 1 ( t )   →   x 2 ( τ   ( t ))   as   t   → ∞ .  Proof:   Suppose there are two nearby trajectories   x 1 , x 2 . We define a smooth mapping   γ   : [0 ,   1]   →   K   where   γ (0) =  x 1   and   γ (1) =   x 2 , such that   ∂γ ∂s   6   = 0   for all   s . Using the Riemannian metric   M   ( x )   and associated dis- tance function   V   ( x, δ ) =   √ δ ′ M   ( x ) δ , the length of a smooth path between   x 1   and   x 2   is given by  L ( γ ) =  ∫   1 0  V  (  γ ( s ) , ∂ ∂s γ ( s )  )  ds   (16) Let   Γ( x 1 , x 2 )   be the set of all smooth paths between   x 1   and  x 2 . Then, the geodesic distance between   x 1   and   x 2   is given by  d ( x 1 , x 2 ) =   min  ι ∈ Γ( x 1 ,x 2 )   L ( ι )   (17)
and   ι ( s )   is the geodesic curve. We   prove,   by   contradiction,   that   Theorem   2   proves Zhukovski stability between any   x 1 , x 2   ∈   K . Suppose   x 1   and   x 2   diverge under some time reparameter- ization. Then there exists a supremum   k , where   k   ∈   [0 ,   1] , for which   ι ( k )   no longer converges to   x 1 . But by Theorem 2, since   ι ( k )   ∈   K , one can construct at   ι ( k )   a local ball of constant radius where immediately nearby trajectories would remain in that ball after the impulse, which contradicts the proposition. Therefore,   x 1 ( t )   →   x 2 ( τ   ( t ))   as   t   → ∞ .  Theorem 4:   If all conditions of Theorem 2 are satisfied, there exists a unique limit cycle that is orbitally stable.  Proof:   Since   K   is strictly forward invariant and com- pact, it follows that the omega-limit set,   Ω( x ) , exists and is a compact subset of   K . Further, an implication of Theorem 3 is that all points in   K   have the same   ω -limit set, which we denote   Ω( K ) . Pick a point   x ?   in   Ω( K ) , by strict forward invariance, this is an interior point of   K . Assume that   f   ( x ? )   6   = 0 , otherwise the results of [13] prove convergence to an equilibrium. Construct   a   hyperplane   orthogonal   to   f   ( x ? ) ,   which   we denote by   H . We prove convergence to a limit cycle by constructing a Poincar ́ e map on   H . Since   f   ( · )   is smooth, for   x   in some neighbourhood   B   of  x ?   we have that   f   ( x ) ′ f   ( x )   >   0 , so in   B H   :=   B   ∩ H   solution curves are transversal to   H   and pass through it in the same direction as at   x ? . Since   x ?   is in the   ω -limit set for all points in   K , and   B H  is transversal, the evolution of the system from any point  x ( t )   ∈   B H   eventually passes through   B H   again. That is,  x ( t   +   s )   ∈   B H   where   s >   0   depends on   x . This evolution can be represented by a Poincar ́ e map   T   :   B H   →   B H   . Take the distance between two points   d ( x 1 , x 2 )   in   B H  to be the Riemannian metric distance from Theorem 3. By Theorem 3, we have that   d ( T   ( x 1 ) , T   ( x 2 ))   <   d ( x 1 , x 2 ) . Hence,   T   is a contractive map from   B H   unto itself. By the Banach fixed point theorem it has a unique stable fixed point, which is its only limit point so must be   x ? . By standard results on Poincar ́ e maps this implies that   x ?   is a point on a limit cycle, to which all solutions converge, by Theorem 3. IV. C ONVEX   C RITERIA FOR   L IMIT   C YCLE   S TABILITY IN  H YBRID   S YSTEMS  In this section, we give convex conditions for transverse contraction of hybrid systems, enabling the search for the metric via sum-of-squares programming.  Condition 5 (Metric Condition linear in   W   ):   Suppose the normal vector of the switching surface is   z ( x ) . If the Riemannian   metric   M   and   the   continuous   dynamics   f  satisfies  α ( x ) f   ( x )   −   W   ( x ) z ( x ) =   β ( x ) c ( x )   (18) for some scalar function   α ( x )   ≥   0 , then Condition 1 is satisfied and all trajectories approach the switching surface orthogonally.  Proof:   For the orthogonality condition of (5) in Con- dition 1 to hold, we require   f   ′ M   ( x ) δ x   = 0   to hold for all  z ( x ) ′ δ x   = 0 . This is equivalent to requiring for some scalar  α ( x )   >   0 , the following holds  z ( x ) ′   =   α ( x ) f   ( x ) ′ M   ( x )   (19) for all   x   ∈   S . Reformulating this in terms of   W   :=   M   − 1  we yield the requirement  α ( x ) f   ( x ) =   W   ( x ) z ( x )   (20) for all   c ( x ) = 0 . Using an S-procedure formulation, we yield the equivalent equality constraint:  α ( x ) f   ( x )   −   W   ( x ) z ( x ) =   β ( x ) c ( x )   (21) which is the required condition for all   x .  Theorem 6 (Convex Conditions for Limit Cycle Stability in Hybrid Non  Suppose, firstly, there exists a Riemannian metric   M   ( x )  for nonlinear system (1) and (2) which satisfies Remark 5; secondly, the continuous dynamics of the system satisfies  W   ( x ) A ( x ) ′   + A ( x ) W   ( x ) −    ̇ W   ( x )+2 λW   ( x ) − ρ ( x ) Q ( x )   ≤   0  (22) where   W   ( x ) :=   M   − 1 , Q ( x ) :=   f   ( x ) f   ( x ) ′ , and thirdly, the discrete switching dynamics   g ( x )   of the system satisfies the following LMI  [  W   ( x ) +   ζ ( x ) Q ( x )   W   ( x )   ∂g ∂x T ∂g ∂x   W   ( x )   W   ( x )  ]  ≥   0   (23) For some   ζ ( x )   ≥   0 , then the overall hybrid system is contracting with respect to metric   M   .  Proof:   By   Remark   5,   just   prior   to   switching,   all trajectories approach the switching surface   S   orthogonally. During the discrete jump, we require the following condi- tion to be satisfied  δ ′  x  (   ∂g ∂x  ′  M ∂g ∂x   −   M  )  δ x   ≤   0   (24) for all   δ ′  x M f   = 0 . Reformulating in terms of the gradient of the metric, i.e.  η   :=   M   ( x ) δ x   such that   δ   =   M   − 1 η   :=   W η , we yield the equivalent condition:  η ′  (  W ∂g ∂x  T  W   − 1   ∂g ∂x W   −   W  )  η   ≤   0   (25) The transversality condition   δ ′ M f   = 0   becomes   η ′ f   ( x ) = 0 . Now, define matrix function   Q ( x ) :=   f   ( x ) f   ( x ) ′   which is rank-one and positive-semidefinite. Hence the sets   { η   :  η ′ f   ( x ) = 0 } ,   { η   :   η ′ Qη   = 0 }   and   { η   :   η ′ Q ( x ) η   ≤   0 }   are equivalent. Now, the transverse contraction of the discrete switching can be proved by the existence of   W   ( x )   such that:  η ′ Q ( x ) η   ≤   0   ⇒  η ′  (  W ∂g ∂x  T  W   − 1   ∂g ∂x W   −   W  )  η   ≤   0   (26)
By the S-procedure, the condition is only true if and only if there exists   ζ ( x )   ≥   0   such that  η ′  (  W   +   ζ ( x ) Q ( x )   −   W ∂g ∂x  T  W   − 1   ∂g ∂x  )  η   ≥   0   (27) By the Schur Complement, (27) is true if and only if (23) holds, which completes the proof. Note that these conditions are all linear in the unknown functions   W   ( x ) , α ( x ) , β ( x ) , ρ ( x ) ,   and   ζ ( x ) , i.e., it consists of a linear matrix inequality at each point   x . For polyno- mial systems, these conditions can be verified efficiently using   sum-of-squares   programming   and   positivstellensatz arguments [18]. V. A PPLICATION   E XAMPLE  The rimless wheel is a simple planar model of dynamic walking, exhibiting hybrid (switching) behaviour. It consists of a central mass, of mass   g,   with equally spaced spikes, of length   l   extending radially outwards. The system rolls down an incline of pitch   γ , as shown in Fig. 2.  Fig. 2.   The rimless wheel model  At any given moment, the rimless wheel rotates about the stance foot without slipping, behaving like an inverted pendulum. When the next foot contacts the ground, it is assumed that an elastic collision occurs such that the old stance foot lifts off and the system now rotates about the new stance foot. The Rimless Wheel state space   x   = [ θ,    ̇ θ ] ′   can be repre- sented with the following hybrid system dynamics:  d dt  [ θ   ̇ θ  ]  =   f   ( θ,    ̇ θ ) =  [    ̇ θ  g l   sin   θ  ]  for   θ   −   γ   −   α   6   = 0  (28)  [ θ +   ̇ θ +  ]  =   g (  ̇ θ − ) =  [   γ   −   α  cos(2 α )  ̇ θ −  ]  for   θ   −   γ   −   α   = 0  (29) On a sufficiently inclined slope, the system has a stable limit cycle, for which the energy lost in collision is perfectly compensated by the change in potential energy. The system has been studied extensively and its basin of attraction has ben computed exactly. [19] Figure 3 shows the phase portrait of the rimless wheel, with blue arrows indicating the direction of the continuous dynamics. The dotted line on the right of the graph indicates the collision surface that maps to the left edge of the graph −0.2   −0.1   0   0.1   0.2   0.3   0.4   0.5   0.6  −1  −0.5  0  0.5  1  1.5  2  2.5  3  3.5  4   Phase Portrait Trajectories of a Rimless Wheel  θ   [rad] d θ /dt [rad/s]  Fig. 3.   Phase diagram of the rimless wheel model  (or vice-versa, depending on the direction of dynamics). The grey line represents the homoclinic orbits of the system, and the red line represents the stable limit cycle. Using the convex conditions of Theorem 6, we formulate sum-of-squares (SOS) and Positivestellansatz conditions [20] which verifies transverse contraction for the hybrid system in a region around the limit cycle, defined by the switching surfaces and a B ́ ezier polynomial   b ( x ) . Let   H   =   A ( x ) W   ( x ) +   W   ( x ) A ( x ) ′   −    ̇ W   ( x ) + 2 λW   ( x ) , let   Σ n [ x ]   denote the set of   n   ×   n   matrices verified positive semidefinite. We approximate the continuous dynamics   f   ( x )  with a third order taylor series expansion. The conditions verified are given below.  W   ( x )   −   ( f   ( x ) T   f   ( x )   −    ) L 1 ( x )  − ( θ   −   ( γ   −   α )) L 2 ( x )  − (( α   +   γ )   −   θ ) L 3 ( x )  − (  ̇ θ   −   b ( x )) L 4 ( x )   ∈   Σ n [ x ]  (30)  αf   ( x )   −   W   ( x ) ∇ c   =   β ( x ) c ( x )   (31)  − H   −   ρ ( x ) f   ( x ) f   ( x ) ′  − ( f   ( x ) T   f   ( x )   −    ) L 5 ( x )  − ( θ   −   ( γ   −   α )) L 6 ( x )  − (( α   +   γ )   −   θ ) L 7 ( x )  − (  ̇ θ   −   b ( x )) L 8 ( x )   ∈   Σ n [ x ]  (32)  [  W   ( x ) +   ζ ( x ) Q ( x )   W   ( x )   ∂g ∂x T ∂g ∂x   W   ( x )   W   ( x )  ]  ∈   Σ 2 n [ x ]   (33)  L 1 , L 2 , L 3 , L 4 , L 5 , L 6 , L 7 , L 8 , β ( x )   ∈   Σ n [ x ]   (34)  α ( x ) , ρ ( x ) , ζ ( x )   ∈   Σ[ x ]   (35)
(30) verifies the positive-definiteness of   W   within the defined region; (31) verifies the condition of Remark 5; (32) and (33) verifies the conditions of Theorem 6; and, finally, (34) and (35) verifies positive semi-definiteness of the Lagrange multipliers and scalar functions. The above conditions were formulated in YALMIP [21], [22]   and   solved   by   commercial   SDP   solver   MOSEK v.7.0.0.103. The code has been made available online [23]. We found that these conditions could be verfified with  W   ( x )   and   β ( x )   a matrix of degree-four polynomials, and  L i ( x ) , α ( x ) , ζ ( x ) , ρ ( x )   degree-two. Figure 4 shows verified regions of stability coloured in green. −0.2   −0.1   0   0.1   0.2   0.3   0.4   0.5   0.6  −0.5  0  0.5  1  1.5  2  2.5  3  3.5  4  4.5  5   Phase Portrait of a Rimless Wheel  θ   (radians) d θ /dt [radians/sec]  System Trajectory  Stable Limit Cycle  Verified Region of Stability  Fig. 4.   Verified region of transverse contraction for Rimless Wheel  VI. C ONCLUSION  We have derived differential conditions guaranteeing the orbital stability of nonlinear hybrid limit cycles. These con- ditions are presented as pointwise linear matrix inequalities, enabling an efficient search for a stability certificate. The main advantages of this approach over traditional Lyapunov-based methods are two-fold. Firstly, the transverse contraction framework decouples the question of conver- gence from knowledge of a particular solution. This opens doors to robustness analysis when the exact location of the limit cycle is unknown due to uncertainty in the dynamics. Further, this method simplifies the search for stability certificate compared with previous Lyapunov-based method in   [4],   which   requires   a   separate   search   for   transversal surfaces and valid Lyapunov functions on those surfaces. By encapsulating the direction of transversal surfaces with the definition of orthogonality, this method allows the search for stability certificate by a single convex optimization problem – a search for a valid transverse contraction metric   M   ( x ) . The ability to efficiently compute stability certificates for hybrid systems in this work opens opportunities for control design with guaranteed stabilizability [24], and would enable the search for provably stable models in system identification [25], [26]. R EFERENCES [1]   S. Collins, A. Ruina, R. Tedrake, and M. Wisse, “Efficient bipedal robots based on passive-dynamic walkers,”   Science , vol. 27, no. 5712, p. 789, Oct. 2005. [2]   E. Westervelt, C. Chevallereau, B. Morris, J. Grizzle, and J. Ho Choi,   Feedback Control of Dynamic Bipedal Robot Locomotion , ser. Automation and Control Engineering.   CRC Press, Jun. 2007. [3]   A. Shiriaev, L. Freidovich, and I. Manchester, “Can we make a robot ballerina perform a pirouette? Orbital stabilization of periodic motions of underactuated mechanical systems,”   Annual Reviews in Control , vol. 32, no. 2, pp. 200–211, Dec. 2008. [4]   I. R. Manchester, M. M. Tobenkin, M. Levashov, and R. Tedrake, “Regions of Attraction for Hybrid Limit Cycles of Walking Robots,”  arXiv preprint arXiv: 1010.2247 , Oct. 2010. [5]   J. Guckenheimer and P. Holmes,   Nonlinear Oscillations, Dynamical Systems, and Bifurcations of Vector Fields .   Springer-Verlag, 1997. [6]   C. S. Hsu, “A Theory of Cell-to-Cell Mapping Dynamical Systems,”  Journal of Applied Mechanics , vol. 47, no. 4, pp. 931–939, 1980. [7]   A. Schwab and M. Wisse, “Basin of attraction of the simplest walking model,”   ASME Design Engineering Technical Conference , 2001. [8]   H. K. Khalil,   Nonlinear systems , 3rd ed.   Upper Saddle River: Prentice Hall, 2002. [9]   J. M. Gonc  ̧alves, A. Megretski, and M. A. Dahleh, “Global analysis of piecewise linear systems using impact maps and surface Lyapunov functions,”   Automatic Control, IEEE Transactions on , vol. 48, no. 12, pp. 2089–2106, 2003. [10]   J. Gonc  ̧alves, “Regions of stability for limit cycle oscillations in piecewise linear systems,”   Automatic Control, IEEE Transactions on , vol. 50, no. 11, pp. 1877–1882, 2005. [11]   I. Manchester, “Transverse dynamics and regions of stability for nonlinear hybrid limit cycles,”   arXiv preprint arXiv:1010.2241 , pp. 1–9, Oct. 2010. [12]   E. M. Aylward, P. A. Parrilo, and J.-J. E. Slotine, “Stability and robustness analysis of nonlinear systems via contraction metrics and SOS programming,”   Automatica , vol. 44, no. 8, pp. 2163–2170, Aug. 2008. [13]   W. Lohmiller and J.-J. E. Slotine, “On Contraction Analysis for Non- linear Systems,”   Automatica , vol. 34, no. 6, pp. 683–696, Jun. 1998. [14]   I. R. Manchester and J.-J. E. Slotine, “Transverse contraction criteria for existence, stability, and robustness of a limit cycle,”   Systems & Control Letters , vol. 63, pp. 32–38, Sep. 2014. [15]   J. K. Hale,   Ordinary Differential Equations .   New York: Robert E. Krieger Publishing Company, Jan. 1980. [16]   G. A. Leonov, “Generalization of the Andronov-Vitt Theorem,”   Reg- ular and Chaotic Dynamics , vol. 11, no. 2, pp. 281–289, 2006. [17]   J. Hauser and C. C. Chung, “Converse Lyapunov functions for expo- nentially stable periodic orbits,”   Systems & Control Letters , vol. 23, no. 1, pp. 27–34, Jul. 1994. [18]   W. Tan, “Nonlinear Control Analysis and Synthesis using Sum-of- Squares Programming,” Ph.D. dissertation, UC Berkeley, 2006. [19]   M. Coleman, “A Stability Study of a Three-Dimensional Passive- Dynamic Model of a Human Gait,” Ph.D. dissertation, Cornell Uni- versity, 1998. [20]   P. A. Parrilo, “Semidefinite programming relaxations for semialgebraic problems,”   Mathematical Programming , vol. 96, no. 2, pp. 293–320, May 2003. [21]   J. Lofberg, “YALMIP: A toolbox for modeling and optimization in MATLAB,”   IEEE International Symposium on Computer Aided Control Systems Design , 2004. [22]   ——, “Pre-and post-processing sum-of-squares programs in practice,”  Automatic Control, IEEE Transactions on , vol. 54, no. 5, pp. 1007– 1011, 2009. [23]   “Rimless   Wheel   transverse   contraction   example,”   2014. [Online].   Available:   http://www-personal.acfr.usyd.edu.au/ian/doku. php?id=wiki:software [24]   I. R. Manchester and J.-J. E. Slotine, “Control Contraction Metrics and Universal Stabilizability,”   arXiv , vol. 1311.4625, 2013. [25]   M. M. Tobenkin, I. R. Manchester, J. Wang, A. Megretski, and R. Tedrake, “Convex optimization in identification of stable non-linear state space models,”   49th IEEE Conference on Decision and Control (CDC) , vol. 1, no. 4, pp. 7232–7237, Dec. 2010. [26]   I. R. Manchester, M. M. Tobenkin, and J. Wang, “Identification of nonlinear systems with stable oscillations,”   IEEE Conference on Decision and Control and European Control Conference , pp. 5792– 5797, Dec. 2011.