arXiv:1311.4625v2 [math.OC] 20 Nov 2013  Control Contraction Metrics and Universal Stabilizability  Ian R. Manchester   ∗   Jean-Jacques E. Slotine   ∗∗ ∗   ACFR, School of Aerospace, Mechanical and Mechatronic Engineering, University of Sydney, Australia. ian.manchester@sydney.edu.au  ∗∗   Nonlinear Systems Laboratory, Dept. of Mechanical Engineering, Massachusetts Institute of Technology  Abstract:   In this paper we introduce the concept of universal stabilizability: the condition that every solution of a nonlinear system can be globally stabilized. We give sufficient conditions in terms of the existence of a control contraction metric, which can be found by solving a pointwise linear matrix inequality. Extensions to approximate optimal control are straightforward. The conditions we give are necessary and sufficient for linear systems and certain classes of nonlinear systems, and have interesting connections to the theory of control Lyapunov functions.  1. INTRODUCTION Constructive control design for nonlinear systems remains a   challenging   problem,   even   in   the   case   of   full   state feedback [Slotine and Li 1991], [Isidori 1995], [Kokotovi ́ c and Arcak 2001]. The classical Lyapunov stability theory leads to necessary and sufficient conditions in terms of existence of control Lyapunov functions, however these may be difficult to find. A full state feedback linearization allows linear control tools to be applied; however finding such a transform can be non-trivial even when it can be shown that one exists [Isidori 1995]. Constructive design tools tend to be more limited in their scope of application. Backstepping and related methods provide a systematic and constructive procedure [Krstic et al. 1995], but are generally limited to systems of a particular “triangular” structure. For certain mechanical and electrical systems, energy and passivity methods can be used [van der Schaft 1999]. Recently, there has been increased interest in methods for systems analysis and control design based on convex op- timization, for example linear matrix inequalities, integral quadratic constraints, and sum-of-squares programming [Boyd et al. 1994], [Megretski and Rantzer 1997], [Parrilo 2000]. For nonlinear control design, the density functions of [Rantzer 2001], [Prajna et al. 2004] and related techniques of occupation measures [Lasserre et al. 2008] and control Lyapunov measures [Vaidya et al. 2010] explicitly address convexity of criteria. Another approach is to piece together stabilized trajectories, with regions of stability verified via sum-of-squares programming [Tedrake et al. 2010]. Constructing a control-Lyapunov or density function re- quires prior knowledge of the solution to be stabilized. However, in many applications with complex system archi- tectures, the role of feedback control is to track a setpoint  ⋆   This work was   supported in part by the   Australian Research Council.  or trajectory that changes in real-time. For nonlinear sys- tems the problem of stabilizing changing trajectories can be quite different to stabilizing a single   a priori   known trajectory. An alternative to explicit design is model predictive con- trol, based on real-time numerical optimization. For many years this method was limited to relatively slow linear processes, but rapid advances in computer technology and optimization algorithms mean it is emerging as a feasible tool for nonlinear problems (see, e.g. Diehl et al. [2009]). Despite some clear benefits, it generally remains difficult to predict or analyse performance of nonlinear MPC schemes by any method other than exhaustive simulations. In this paper, we introduce the concept of   universal sta- bilizability , i.e. the property that   every   solution of a sys- tem is globally stabilizable. Using the tools of contraction analysis [Lohmiller and Slotine 1998b], we give convex sufficient conditions for universal stabilizability, and a con- structive control design procedure. Useful properties such as exponential stability and approximate linear quadratic optimality are natural extensions, and for the case of linear systems, our conditions reduce to well-known necessary and sufficient conditions. Historically, basic convergence results on contracting sys- tems can be traced back to the results of Lewis [1949] in terms of Finsler metrics, and results of Hartman [1961] and Demidovich [1962]. Extensions to analysis of limit cycles have recently been developed by the authors [Manchester and Slotine 2013] building upon the early work of Borg [1960] and Hartman and Olech [1962]. Connections with Finsler structures and Lyapunov theory have recently been explored in detail by Forni and Sepulchre [2012]. The conditions we give are state-dependent linear matrix inequalities. The conditions are formally similar to those studied for a class of nonlinear   H ∞   problems in Lu and Doyle [1995] and to state-dependent Riccati equations (see, e.g., Cloutier [1997]), however the class of systems and   the   explicit   control   construction   we   consider   are
quite different. An interesting feature of our method is that it breaks the control-design problem into two stages: computation   of   a   metric,   which   can   be   performed   in advance of knowledge of the solution to be stabilized, and an on-line path integration to compute the control signal. Our methods build significantly upon the control design suggested in Lohmiller and Slotine [1998a]. In that paper, it was argued that requiring a nonlinear system to be feedback transformable to a stable linear system – as in feedback linearization – was too strong, since it requires an involutivity condition to be added to the rather nat- ural controllability condition [Isidori 1995]. If instead one required only that the system be feedback transformable to a “nice” stable nonlinear system (a contracting system), then this can be achieved using the controllability condi- tion alone. 2. UNIVERSAL STABILIZABILITY For most of this paper, we will consider a nonlinear time- dependant control-affine system  ̇ x ( t ) =   f   ( x ( t ) , t ) +   B ( t ) u ( t )   (1) where   x ( t )   ∈   R n , u ( t )   ∈   R m   are   state   and   control, respectively, at time   t   ∈   R +   :=   [0 ,   ∞ ). The function  f   :   R n   ×   R +   →   R n   is assumed to   be   smooth, and  B   :   R +   →   R n × m   is a time-dependent matrix. Contraction analysis is the study of (1) by way of the associated system of differential dynamics:  ̇ δ x ( t ) =   A ( x, t ) δ x ( t ) +   B ( t ) δ u ( t )   (2) where   A ( x, t ) =   ∂  ∂x   f   ( x, t ) is the Jacobian matrix. In Section 8 we will consider more general systems  ̇ x   =  f   ( x, u, t ) but for the moment we note that many systems not naturally appearing in the form (1) can be put in that form, either exactly or approximately, by change of variables or introducing new states. A solution of (1) is a pair of vector signals ( x ( t ) , u ( t )) satisfying (1) over the interval   R + . For simplicity, in this paper we will assume that (1) is such that solutions exist and are unique. As such, we will also use the notation  φ ( t, x 0 , u ) to denote the solution of   x   at time   t >   0 of (1) starting from initial condition   x (0) =   x 0   and under the application of the control input   u ( τ   ) , τ   ∈   [0 , t ]. A static state-feedback controller is a function   k   :   R n   →  R m . The main objective is to design such a function so that the behaviour of the closed-loop system  ̇ x ( t ) =   f   ( x ( t ) , t ) +   B ( t ) k ( x ( t ) , t )   (3) is in some sense desirable, e.g. globally stable or optimal in some sense. Solutions of the closed-loop system with a given state-feedback controller   k   will be denoted   φ k   ( t, x 0 ). As is standard, a solution ( x ⋆ , u ⋆ ) defined on [0 ,   ∞ ) is said to be globally asymptotically stabilized by a feedback controller   u   =   k ( x, t ) if both of the following hold: (1) For any   α   there exists an   ǫ   such that   | x 0   −   x ⋆ (0) |   < ǫ  implies   | φ k   ( t, x 0 )   −   x ⋆ ( t ) |   < α , (2) For any initial condition   x 0   ∈   R n , the closed loop solution satisfies   | φ k   ( t, x 0 )   −   x ⋆ ( t ) | →   0. Global   exponential   stabilization   refers   to   the   stronger condition that there exists a   K   and   λ   such that  | φ k ( t, x 0 )   −   x ⋆ ( t ) | ≤   Ke − λt | x 0   −   x ⋆ (0) |  for all   x (0). In this paper, we will study the following property:  Definition 1.   A   system   of   the   form   (1)   is   said   to   be  universally stabilizable by state feedback   if for any solution ( x ⋆ , u ⋆ ) defined on   t   ∈   [0 ,   ∞ ) there exists a state feedback controller   k   :   R n   →   R m   such that ( x ⋆ , u ⋆ ) is globally stabilized by   u   =   k ( x, t ). Analogously, we also consider the notion of   universally exponentially stabilizable with rate   λ . Note that universal stabilizability is a stronger condition than global stabilizability of a particular solution. 3. CONTROL CONTRACTION METRICS The following definition is central to this paper:  Definition 2.   A function   V   ( x, δ x , t ) =   δ ′  x M   ( x, t ) δ x , with  c 1 I   ≤   M   ( x, t )   ≤   c 2 I   for some   c 2   ≥   c 1   >   0, is said to be a  control contraction metric   for the system (1) if   ∂V  ∂x   B ( t ) = 0 and  ∂V  ∂δ x  B ( t ) = 0 = ⇒   ∂V  ∂t   +   ∂V  ∂x f   ( x, t )+   ∂V  ∂δ x  A ( x, t ) δ x   <   0 (4) for all   x, t . We will also make use of a Riemannian distance function between any two points at a given time   d ( x 1 , x 2 , t ) :   R n   ×  R n   ×   R +   →   R +   defined like so: let Γ( x 1 , x 2 ) denote the set of all smooth paths connecting   x 1   and   x 2 , where each   γ   ∈   Γ( x 1 , x 2 ) is parametrised by   s   ∈   [0 ,   1], i.e.  γ ( s ) : [0 ,   1]   →   R n . The path length of   γ   is then defined as  L ( γ, t ) :=  ∫   1 0  D  (  γ ( s ) , ∂  ∂s γ ( s ) , t  )  ds  where   D   is a Finsler function [Bao et al. 2000]. In this paper, we will primarily choose   D ( x, δ x , t ) =   √ V   ( x, δ x , t ) or   D ( x, δ x , t ) =   | δ x | . The distance between two points is then defined as  d ( x 1 , x 2 , t ) =   min  γ ∈ Γ( x 1 ,x 2 )   L ( γ, t ) The existence of a minimizing path, which we denote  γ x 2  x 1   ( t, s ), is implied by the Hopf-Rinow Theorem.  Remark 1.   In this conference paper we restrict ourselves to simple choices of   V   and   D . Other choices are possible, including non-quadratic metrics and matrix measures. See, e.g., Lohmiller and Slotine [1998a] and Forni and Sepulchre [2012] for detailed discussion of non-quadratic metrics in contraction analysis. The main utility of a control contraction metric is that it allows one to construct stabilizing controls by computing path integrals of a suitable   δ u .  Theorem 1.   If a control contraction metric exists for a system of the form (1), then the system is universally stabilizable by static state feedback.  Proof:   Consider an arbitrary initial condition   x (0)   ∈   R n  and a desired feasible trajectory ( x ⋆ ( t ) , u ⋆ ( t )) defined on  t   ∈   [0 ,   ∞ ). For each time   t , consider the distance   d ( x ⋆ ( t ) , x ( t ) , t ). By the construction of the distance given above, if   d   →   0 then
x ( t )   →   x ⋆ ( t ). So the objective of this proof is to show the existence of a control law such that   d   →   0 as   t   → ∞ . Consider also the minimal path   γ ( t, s ) =   γ x 2  x 1   ( t, s ). For a given time   t , at each point   s   ∈   [0 ,   1] along the path, consider two possibilities: firstly, suppose  ∂  ∂δ x  V   ( x, δ x , t ) B ( t ) = 0 when the left hand side is evaluated at   x   =   γ ( t, s ) , δ x   =  ∂γ ( s,t )  ∂s   . Then by the definition of a control contraction metric, one has  d  dt V   =   ∂V  ∂t   +   ∂V  ∂x f   ( x, t ) +   ∂V  ∂δ x  A ( x, t ) δ x   <   0 .  Secondly, suppose that  ∂  ∂δ x  V   ( x, δ x , t ) B ( t )   6   = 0 at   x   =   γ ( t, s ) , δ x   =   ∂γ ( s,t )  ∂s   . Then  d  dt V   =   ∂V  ∂t   +   ∂V  ∂x f   ( x, t ) +   ∂V  ∂δ x  ( A ( x, t ) δ x   +   B ( t ) δ u ) is affine in   δ u , and therefore there exists   δ u   making   d  dt   V <  0. Choose any such   δ u   and denote it   δ ⋆ u ( t, s ). Let us construct a control signal   u ( t, s ) for each point on the path   γ ( t, s ) like so:  ̄ u ( t, s ) =   u ⋆ ( t ) +  ∫   s  0  δ ⋆ u ( t,   s ) d s  then we have   u ( t,   0) =   u ⋆ ( t ) and   ∂    ̄ u ( t,s )  ∂s   =   δ ⋆ u ( t, s ) This implies that for all ( t, s ) we have  d  dt V  (  γ ( s ) , ∂  ∂s γ ( s ) , t  )  <   0 and therefore as   t   → ∞ , we have   V   ( γ ( s ) ,   ∂  ∂s   γ ( s ) , t )   →   0 for all   s , implying   D   ( γ ( s ) ,   ∂  ∂s   γ ( s ) , t )   →   0 for all   s , further implying   d ( x ⋆ ( t ) , x ( t ))   →   0. The specific control signal that is applied to the system is  ̄ u ( t,   1). Note that this procedure defines a static state- feedback function   u   =   k ( x, t ).  ✷  4. CONVEX CONDITION AND CONSTRUCTION OF EXPLICIT CONTROL The second main result of this paper is the following, where we use the notation   S n  +   to denote the convex cone of symmetric positive semidefinite matrices:  Theorem 2.   Consider the system (1) with differential dy- namics (2). If there exists a matrix function   W   ( x, t ) :   R n   ×  R   →   S n  + , a function   ρ ( x, t )   ≥   0, and constants   α 2   ≥   α 1   >   0 such that  α 1 I   ≤   W   ≤   α 2 I,   (5)  −    ̇ W   +   W A ′   +   AW   −   ρBB ′   <   0 ,   (6) for all   x, t , then the system (1) is universally stabilizable by static state feedback, and   δ x W   ( x, t ) − 1 δ x   is a control contraction metric. In the above condition,    ̇ W   ( x, t ) is a matrix with the   i, j  element given by   ∂W i,j  ∂t   +   ∂W i,j  ∂x   ( f   ( x, t ) +   B ( t ) u ).  Proof:   We will show that by taking   M   ( x, t ) =   W   ( x, t ) − 1 , the conditions of the theorem prove that   V   ( x, δ x , t ) =  δ ′  x M   ( x, t ) δ x   is a control contraction metric. Indeed, make the substitution of   M   for   W   and multiply (6) on either side by   M   , the conditions of the theorem are equivalent to  ̇ M   +   A ′ M   +   M A   −   ρM BB ′ M <   0 ,   (7)  ∀ x, u   with   1  α 2   I   ≤   M   ( x, t )   ≤   1  α 1   I . Note that    ̇ M   ( x, t ) =  − M   ( x, t )  ̇ W   ( x, t ) M   ( x, t ) follows from the total derivative of the identity   M   ( x, t ) W   ( x, t ) =   I   with respect to time. Let   γ t ( s ) be the minimizing path connecting   x ⋆ ( t ) to   x ( t ), parameterized by   s   ∈   [0 ,   1] and consider the control signal generated by  u ( t ) =   u ⋆ ( t )   −   1  2  ∫   1 0  ρ ( γ t ( s ) , t ) B ( t ) ′   M   ( γ t ( s ) , t )   ∂γ t  ∂s ds  is stabilizing to   x ⋆ ( t ). Now, along the path   γ , we have   δ u   =   −   1  2   B ′ M   ( x ) δ x   so  d  dt δ ′  x M δ x   = δ ′  x    ̇ M δ x   + 2 δ ′ M   [ Aδ   +   Bδ u ] = δ ′  x [  ̇ M   +   A ′ M   +   M A   −   ρM BB ′ M   ] δ x   <   0 (8) where the last inequality is implied by (7).  ✷  Remark 2.   We note that the conditions are pointwise lin- ear matrix inequalities in   W   ( x, t ) and   ρ ( x, t ), and are thus amenable to solution via convex optimization algorithms. We will address this further in Section 6. For systems of the form (1), the matrix    ̇ W   has elements  ̇ W i,j   =   ∂W i,j  ∂t   +   ∂W i,j  ∂x   ( f   ( x, t ) +   B ( t ) u ), which are affine in the control signal   u . Therefore for inequality (6) to hold for all   u , it is necessary that   ∂W i,j  ∂x   B ( t ) = 0, i.e.   W   ( x, t ) must be constant along the subspace   S u   of state space in which the control signal directly acts. This property deserves further investigation, but we re- mark that if   f   ( x, t ) is globally Lipschitz, then tracking errors in   S u   will generally be stabilized by sufficiently high- gain linear feedback as verified by a Lyapunov function that is quadratic with respect to   S u . This is closely con- nected to the method of sliding-mode control [Slotine and Li 1991]. Alternatively, dynamics in these directions can be directly cancelled, in a simple application of feedback linearization [Isidori 1995]. We also note that the restriction   ∂W i,j  ∂x   B ( t ) = 0 is removed in Section 8, at the expense of a more complex control computation.  4.1 Conditions for Exponential Stabilization  A stronger result is the following:  Theorem 3.   Consider the system (1) with differential dy- namics (2). If there exists a matrix function   W   ( x, t )   ∈   S n  + ,  ρ ( x, t )   ≥   0 and   α 2   ≥   α 1   >   0 such that  α 1 I   ≤   W   ≤   α 2 I,   (9)  −    ̇ W   +   W A ′   +   AW   −   ρBB ′   ≤ − 2 λW   ( x ) ,   (10) for all   x, u, t , then the system is universally exponentially stabilizable with rate   λ   by state feedback.
The proof uses the same controller construction and fol- lows the same outline as the first theorem, so we provide a sketch here: by taking   M   ( x, t ) =   W   ( x, t ) − 1   the conditions of the theorem are equivalent to 0   < M   ( x )   <   1  α   I   and  ̇ M   +   A ′ M   +   M A   −   M BB ′ M   ≤ − 2 λM.   (11) and, under the action of the control,  d  dt   ( δ ′  x M   ( x, t ) δ x )   ≤ − 2 λδ ′  x M   ( x, t ) δ x .  So there exists a   K 1   such that   V   ( s, t )   ≤   K 1 e − 2 λt V   ( s,   0) for all   s . Now, from the definition of   V   and   D   we can give  c 2   ≥   c 1   >   0 such that   c 1  √ V   ( t, s )   ≤   D ( t, s )   ≤   c 2  √ V   ( t, s ) it follows that   D ( t, s )   ≤   Ke − λt   with   K   =   c 2  c 1  √ K 1 . This further implies that  d ( x ⋆ ( t ) , x ( t ))   ≤   Ke − λt d ( x ⋆ (0) , x (0)) which completes the proof.  ✷  5. GUARANTEED QUADRATIC COST CONTROL Explicit solutions of nonlinear optimal control problems are generally extremely challenging to find. A general approach is to formulate solving a nonlinear partial dif- ferential equation, the Hamilton Jacobi Bellman equation, and attempt to find a so-called viscosity solution. Unfor- tunately, in most cases of practical interest this computa- tional problem is intractible. A less demanding approach, frequently used in robust control, is to find a   guaranteed cost control   [Petersen et al. 2000]. In this framework, one tries to find a controller that guarantees a conservative upper bound on the cost, which is usually obtained via dissipation inequalities. Then one can approach approximate optimal control by minimizing this upper bound. In this section we show that our method can be extended to guaranteed cost control of nonlinear systems, by making use of   differential   dissipation inequal- ities [Manchester and Slotine 2013].  Theorem 4.   Given   Q ( t )   >   0   and   R ( t )   >   0   and   the differential dynamics (2), suppose there exists a   W   ( x, t ) such that  [ (  ̇ W   −   W A ′   −   AW   +   BR − 1 B ′ )   W W   Q − 1  ]  ≥   0   (12) for all   x, u, t . Now, for any feasible trajectory ( x ⋆ ( t ) , u ⋆ ( t )), consider the LQ cost function  J   =  ∫   ∞  0  [ ( x ( t )   −   x ⋆ ( t )) ( u ( t )   −   u ⋆ ( t ))  ] ′   [ Q ( t )   0 0   R ( t )  ] [ ( x ( t )   −   x ⋆ ( t )) ( u ( t )   −   u ⋆ ( t ))  ]  dt.  then there exists a state-feedback controller such that  J ( x (0))   ≤  ∫   1 0  δ ′  x W   ( x ) − 1 δ x ds  with   x   =   γ ( t, s ) and   δ x   =   ∂γ ( t,s )  ∂s   over the path   γ x (0)  x ⋆ (0) ( t, s ).  Proof:   Follows from similar argument as above and the following differential dissipation inequality  d  dt   ( δ x M δ x )   ≤ −  [ δ x  δ u  ] ′   [ Q ( t )   0 0   R ( t )  ] [ δ x  δ u  ]  along any path   γ   joining   x ⋆ ( t ) and   x ( t ). This inequality is satisfied by design of controllers satisfy- ing the differential Riccati inequality  ̇ M   +   A ′ M   +   M A   −   M BR − 1 B ′ M   +   Q   ≤   0 which can be transformed to  −    ̇ W   +   W A ′   +   AW   −   BR − 1 B ′   +   W QW   ≤   0 which can be linearized via Schur complement to  [ (  ̇ W   −   W A ′   −   AW   +   BR − 1 B ′ )   W W   Q − 1  ]  ≥   0 .  To show that the differential dissipation inequality implies the cost bound, it is more convenient to consider the equivalent formulation  J   =  ∫   ∞  0  ∣ ∣ ∣ ∣ H ( t )  [ ( x ( t )   −   x ⋆ ( t )) ( u ( t )   −   u ⋆ ( t ))  ]∣ ∣ ∣ ∣  2  dt  where   H ( t ) satisfies  H ( t ) ′ H ( t ) =  [ Q ( t )   0 0   R ( t )  ]  e.g. by Cholesky factorization. The upper bound for the cost can be obtained by noting that the Cauchy-Schwarz inequality on   L 2 [0 ,   1] implies  ∫   1 0  ∣ ∣ ∣ ∣ H ( t )  [ δ x  δ u  ]∣ ∣ ∣ ∣  2  ds   ≥  ∣ ∣ ∣ ∣ ∫   1 0  H ( t )  [ δ x  δ u  ]  ds  ∣ ∣ ∣ ∣  2  and if the right-hand side is computed along any path, it is an upper bound for  ∣ ∣ ∣ ∣ H ( t )  [ ( x ( t )   −   x ⋆ ( t )) ( u ( t )   −   u ⋆ ( t ))  ]∣ ∣ ∣ ∣  2  .  with equality if a straight-line path is used, i.e. a geodesic with respect to   | δ x | .  ✷  As approximate optimal control, one could consider maxi- mizing, e.g. the smallest eigenvalue of   W   ( x ) or the trace of  W   − 1 ( x ), both of which are concave in   W   ( x ), to minimize “worst case” or “expected” value of  ∫   1 0  δ ′  x W   − 1 ( x, t ) δ x ds.  The above construction can be easily extended to cost functions with a cross term ( u ( t ) − u ⋆ ( t )) ′ S ( t )( x ( t ) − x ⋆   ( t )) as long as the cost matrix remains positive-definite. It is well known that for linear systems, stabilizability implies solvability of the LQR problem if   Q >   0 , R >   0 Hespanha [2009]. In fact, for nonlinear systems a similar claim can be made in terms of guaranteed cost controlla- bility.  Theorem 5.   Suppose a system of the form (1) is univer- sally stabilizable, as verified by the conditions of Theorem 2, then there exists a solution to the quadratic guaranteed cost control problem in Theorem 4.  Proof:   Since both   of the conditions are pointwise, it suffices to show at each state and time ( x, t ) that the existence of a   W   ( x, t ) and   ρ ( x, t ) satisfying (5)and  −    ̇ W   +   W A ′   +   AW   −   ρBB ′   <   0   (13) implies   the   existence   of   a   (possibly   different)    ̄ W   ( x, t ) satisfying  −    ̇  ̄ W   +  ̄ W A ′   +   A    ̄ W   −   BR − 1 B ′   +  ̄ W Q    ̄ W   ≤   0 ,   (14)
for a fixed   Q ( t )   >   0 , R ( t )   >   0. Note that we have taken the Schur complement of (12). Let us first consider the condition  δ ′  1  (  −    ̇  ̄ W   +  ̄ W A ′   +   A    ̄ W   −   BR − 1 B ′   +  ̄ W Q    ̄ W  )  δ 1   ≤   0 ,  when   δ 1   is in the nullspace of   B ′ , i.e.   B ′ δ 1   = 0. Clearly this is equivalent to the condition  δ ′  1  (  −    ̇  ̄ W   +  ̄ W A ′   +   A    ̄ W   +  ̄ W Q    ̄ W  )  δ 1   ≤   0 .   (15) Now, the conditions of Theorem 2 imply that for any   δ 1  in the null-space of   B ′ , we have  δ ′  1  (  −    ̇ W   +   W A ′   +   AW  )  δ 1   <   0 .   (16) Take    ̄ W   ( x, t )   =   μ ( x, t ) W   ( x, t ) for some scaling factor  μ ( x, t )   >   0, then (15) becomes  δ ′  1  (  μ ( −    ̇ W   +   W A ′   +   AW   ) +   μ 2 W QW  )  δ 1   ≤   0 .  From (16) it is clear that taking   μ   sufficiently small this condition is satisfied. On the other hand, suppose   δ 2   is not in the null-space of  B , then as   μ   →   0  δ ′  2  (  μ ( −    ̇ W   +   W A ′   +   AW   )   −   BR − 1 B ′   +   μ 2 W QW  )  δ 2  → − δ ′  2 BR − 1 B ′ δ 2   <   0 so again, sufficiently small   μ   will satisfy the conditions of Theorem 4. Since it is sufficient to verify the contraction condition for   δ   :   | δ |   = 1, a compact set, this implies the existence of a sufficiently small   μ ( x, t ) for all   δ .  ✷  Remark 3.   Taking    ̄ W   ( x, t ) =   μW   ( x, t ) very small in the proof above leads to a upper bound on the cost given by the path integral of   1  μ   δ ′ W   − 1 ( x, t ) δ . Therefore the smaller  μ   is taken, the more conservative the upper bound on cost. 6. COMPUTATIONAL APPROACHES The methods we have described naturally break the con- trol problem down into two stages: firstly, a pointwise LMI must be solved for   W   ( x ), giving the metric   δ ′  x M   ( x ) δ x , guaranteeing a particular form of stabilizability; secondly, the on-line computation involves computing a path in- tegral along a line or geodesic. A full discussion is not possible due to space restrictions, but in this section we briefly discuss some applicable techniques.  6.1 Finding a Control Contraction Metric  The main calculation is to find a matrix function   W   ( x, t ) and a scalar function   ρ ( x, t ) satisfying (5) and (6) for all  x, u, t , or the similar conditions in Sections 4.1 and 5. These conditions are convex, but the decision variables  W   and   ρ   are infinite dimensional and there are infinitely many LMI constraints (5) and (6), due to the dependence on   x   and   t . To be practically computable, we must reduce these to a finite-dimensional optimization problem. If   f   ( x, t )   =   f   ( x ) is a vector of polynomials in   x , one can search for polynomial   W   ( x ) and   ρ ( x ) satisfying the conditions of the theorem using sum-of-squares program- ming [Parrilo 2003]. In that framework, a matrix inequality  G ( x )   ≥   0 for all   x , where   G ( x ) is a   q   × q   symmetric matrix, is represented by introducing an auxiliary variable   y   ∈   R q  and search for a representation   y ′ G ( x ) y   =   ∑  i   g i ( x, y ) 2  for some polynomials   g i ( x, y ). When   G ( x ) is linearly pa- rameterized by unknowns, this can be represented as a finite-dimensional semidefinite program. This gives a con- servative test, since not all non-negative polynomials are sums of squares. Alternatively, over a compact set in state space one could specify a finite set of basis functions for   W   and   ρ , and check   the   inequalities   (5)   and   (6)   at   a   finite   grid   of points. This again leads to a finite-dimensional set of linear matrix inequalities. A similar method was investigated by Johansen [2000] for computation of Lyapunov functions.  6.2 Feedback Control via Geodesic Computation  In the case that   D ( x, δ x , t ) =   | δ x |   is used to compute the distance function, construction of the feedback control is relatively straightforward: shortest paths are straight lines, and the control is computed by a path integral along a line. In the simplest case that a   M   ( t ) and   ρ ( t ) can be found independent of   x   and satisfying the conditions, then the feedback control is simply the linear gain  u ( t ) =   u ⋆ ( t )   −   1  2   ρ ( t ) B ( t ) ′   M   ( t )( x ( t )   −   x ⋆ ( t )) .  It   is   to   be   expected   that   better   performance   can   be achieved if   D ( x, δ x , t ) =   √ V   ( x, δ x , t ) is used instead. In that case, the presented approach does not remove the need for on-line computation but it does potentially reduce its complexity, as compared to model predictive control. In particular, we have reduced the problem from an   optimal control   problem to a   shortest path   problem, with respect to a Riemannian distance function. The shortest path problem is easier to solve because of the special structure of the distance metric and the fact that the path need not satisfy any differential equations. Furthermore, it is computed only in the state space, rather than in the state/control space, as with MPC when using collocation or multiple shooting. Algorithms for computing geodesics have been developed in several fields, including computational physics, com- puter graphics, and robot motion planning M ́ emoli and Sapiro [2001], Boykov and Kolmogorov [2003], LaValle [2006]. The approach described in this subsection may also be applicable to the controller construction in [Lohmiller and Slotine 1998a].  6.3 Open-loop Control via Path Images  An alternative method giving an open-loop control signal that drives   x ( t )   →   x ⋆ ( t ) was suggested in [Lohmiller and Slotine 1998a]. This method also uses path integrals of a differential feedback gain   δ u   =   Kδ x , although   K  was constructed by a different method   related to the controllability conditions of Isidori [1995]. In particular, it was suggested to compute a path from   x ⋆ (0) to   x (0), and integrate   δ u   =   Kδ x   along this path for   u (0). For times  t >   0, the forward image of each point on the path under the flow of the system is simulated, generating a family of paths   γ ( s, t ), and controls   u ( s, t ).
This method can also be applied for the systems we study with   K   =   −   1  2   ρB ′ M   . One advantage is that it does not require recomputing a geodesic at each time to compute the control signal. A disadvantage is that it does not give a feedback law, and could be non-robust. E.g., if the true system behaviour does not match the simulated forward image of   x (0) then there are no guarantees of stability. One can also imagine a hybrid of this approach and the previous one, in which forward images of the initial path are computed, but with some corrective term based on gra- dient descent to converge towards geodesics. In this case, the resulting feedback system would be a dynamic state feedback with an infinite-dimensional controller state.  6.4 Two-stage Computation of Explicit Feedback  If   W   ( x, t ) has been found using the above-described con- vex optimization procedure, then one can fix   M   ( x, t ) =  W   ( x, t ) − 1   which implies the existence of matrix   K ( x ) such that  δ ′  x  (    ̇ M   ( x, t ) + 2 M   ( x, t ) [ A ( x, t ) +   B ( t ) K ( x, t )]  )  δ x   <   0 (17) for all   δ x , i.e., the existence of a differential state feedback  δ u   =   K ( x ) δ x . However, in order to compute an explicit control law   u   =   k ( x, t ), this differential relation must be integrable. Let   K i ( x ) refer to the   i th   row of the matrix   K ( x ), then the additional convex constraint  ∂  ∂x K i ( x ) =  (   ∂  ∂x K i ( x )  ) ′  , i   = 1 ,   2 , ...m   (18) ensures the existence of a function   k ( x ) satisfying   ∂k  ∂x   =  K ( x ) by the Poincar ́ e lemma. For example, consider a linearly parameterized class of controllers where the   i th  control element is given by  u i   =   k i ( x ) :=  p ∑  j =1  κ i,j   φ j   ( x ) where   φ j   ( x ) are smooth basis functions (e.g. a monomial basis for polynomial feedback). In this case,  K i ( x ) :=  p ∑  j =1  κ i,j  ∂  ∂x φ j   ( x ) and (18) imposes a linear constraint on the coefficients, as long as the set of basis functions is closed under differen- tiation with respect to   x . This can also be considered as a curl condition on a vector field, so recent work on curl-free wavelets by Deriaz and Perrier [2009] may be applicable when choosing basis functions for the control gain. For stabilizing particular equilibrium   x ⋆ , one can add the additional linear constraint 0 =   f   ( x ⋆ , t ) +   B ( t ) k ( x ⋆ , t ) ,  or, e.g. a periodic solution could be stabilized with family of constraints over [0 , T   ]  ̇ x ⋆   =   f   ( x ⋆ , t ) +   B ( x ⋆ , t ) k ( x ⋆ , t ) .  It is important to note that the conditions of the main theorems imply the existence of a   K ( x, t ) satisfying (17) but not necessarily satisfying the additional constraint (18). To our knowledge, whether or not this   stronger implication is true is an open question. 7. DISCUSSION  7.1 Necessity for the Case of Linear Systems  In general, the conditions we provide will be sufficient but not necessary. An interesting question is whether there is a class of nonlinear systems for which they are necessary, i.e. they completely describe stabilizability. As a basic requirement, one would expect this class to include all stabilizable linear time-invariant (LTI) systems:  ̇ x   =   Ax   +   Bu  and indeed this is the case. Our condition for LTI systems with a matrix   W >   0 independent of   x   is:  AW   +   W A ′   −   BB ′   <   0 where   ρ , also constant,has been absorbed into   W   . This is well-known Lyapunov condition for stabilizability of linear systems – see, e.g., [Hespanha 2009, Sec 14.5]. The condition we give on computing a guaranteed cost control similarly reduces to the Riccati inequality:  AW   +   W A ′   −   BR − 1 B ′   +   W QW   ≤   0 which is known to have a solution if and only if the algebraic Riccati equation:  AW   +   W A ′   −   BR − 1 B ′   +   W QW   = 0 has a solution, given by the solution of the Riccati in- equality with maximal trace. By multiplying the latter on either side by   M   =   W   − 1 , this is the exactly the Riccati inequality associated with the linear quadratic regulator problem, with the minor strengthening that   Q   must be invertible and hence positive-definite. In linear control design, changes of variables such as that from   M   to   W  are frequently used to construct linear matrix inequality conditions [Boyd et al. 1994].  7.2 Connection to Control Lyapunov Functions  The distance function   d ( x ⋆   ( t ) , x ( t ) , t ) is a control Lya- punov function for a solution   x ⋆ , thus the conditions we give imply the existence of a control Lyapunov function for   every   solution of the system. Control Lyapunov functions for systems of the form (1) can be characterized by the condition:  ∂  ∂x V   ( x, t ) B ( t ) = 0   ⇒   ∂  ∂x V   ( x, t ) f   ( x, t ) +   ∂  ∂t V   ( x, t )   <   0 for all   x   and   t . In general, it is computationally challenging to search for a   V   ( x, t ) satisfying such a condition impli- cation. If we instead consider a similar condition for the  differential   dynamics (2) and a function   δ ′  x M   ( x, t ) δ x :  M Bδ x   = 0   ⇒   δ ′  x (  ̇ M   +   A ′ M   +   M A ) δ <   0 Then by Finsler’s theorem [Uhlig 1979], or a version of the S-Procedure losslessness theorem [Yakubovich 1971], this is equivalent to the existence of a non-negative function  ρ ( x, t ) for which the following is true:  ̇ M   +   A ′ M   +   M A   −   ρM BB ′ M <   0 which, by pre and post-multiplication by   W   is our con- dition (6). Therefore our condition guarantees both the existence of a control Lyapunov function in the traditional sense, and a “differential” control Lyapunov function – i.e. a control contraction metric.
7.3 On Transverse Contraction  The results in this paper also suggest a new interpreta- tion of limit cycle stability. In [Manchester and Slotine 2013] the authors introduced the concept of transverse contraction, a generalisation of prior work by Borg [1960], Hartman and Olech [1962], Leonov et al. [1996], to study the existence of limit cycles in autonomous systems:  ̇ x   =   f   ( x ) .  Transverse contraction refers to the existence of a metric function   V   ( x, δ x ) such that  ∂V   ( x, δ )  ∂x   f   ( x ) +   ∂V   ( x, δ )  ∂δ ∂f   ( x )  ∂x   δ   ≤ − λV   ( x, δ ) ,   (19) for all   δ   6   = 0 such that   ∂V  ∂δ   f   ( x ) = 0. In [Manchester and Slotine 2013] it was shown that, for a Riemannian metric √ δ ′ M   ( x ) δ , this is equivalent to a pointwise LMI condition:  ̇ W   ( x )   ≥   W   ( x ) A ( x ) ′   +   A ( x ) W   ( x )   −   ρ ( x ) Q ( x ) ′   + 2 λW   ( x ) ,  (20) where   A ( x )   =   ∂f  ∂x   ,   Q ( x )   =   f   ( x ) f   ( x ) ′   and   W   ( x )   >  0. It was shown that this implies stability under time reparametrisation – a.k.a. Zhukovsky stability – which is known to imply existence of a stable limit cycle. In the context of this paper, that condition is equivalent to stabilizability of the system:  ̇ x   =   f   ( x ) +   f   ( x ) u.  The interpretation of this fact is that a control input acting “along” the flow of the system   f   ( x ) can speed up or slow down the trajectory, but not change its phase portrait. 8. EXTENSION TO MORE GENERAL NONLINEAR SYSTEMS In this section we briefly discuss extension to more general nonlinear systems  ̇ x ( t ) =   f   ( x ( t ) , u ( t ) , t ) with differential dynamics  ̇ δ x ( t ) =   A ( x, u, t ) δ x ( t ) +   B ( x, u, t ) δ u ( t ) .  The significant change from the previous sections is that the   A   and   B   matrices now depend on   u . In this case, we can also consider a matrix   W   ( x, u, t ) and constant   ρ ( x, u, t ). Then, suppose conditions (5) and (6) hold for all   x, u, t . The path integral equation for the control signal   u ( s ) now depends on   u ( s ), and can be constructed if a solution exists to the differential equation:  d  ds u   =   −   1  2   ρ ( γ, u ) B ( γ ) ′ M   ( γ, u )   ∂γ  ∂s   ( s ) with boundary condition   u (0) =   u ⋆ ( t ). If the dependency of these terms on   u   can be well controlled, then it may be possible to guarantee existence by, e.g., the Bellman- Gr ̈ onwall lemma. We defer detailed discussion of this issue for a later publication. 9. APPLICATION EXAMPLE In this section we show an example class of systems for which the above analysis can be applied: compositions of a stabilizable linear system with a contracting system. Consider the following hierarchical composition of sys- tems:  ̇ z   =  ̄ Az   +  ̄ Bu,    ̇ y   =   f 1 ( y, z ) ,  for which   f 1   is partially contracting with respect to   y , and the pair    ̄ A,    ̄ B   is stabilizable. That is, there exists  M 1 ( y )   ≥   αI   and   P   such that  ̇ M 1 ( y ) +   ∂f 1  ∂y  ′  M 1 ( y ) +   M 1 ( y )   ∂f 1  ∂y   < 0   (21)  ̄ AP   +   P    ̄ A ′   −    ̄ B    ̄ B ′   < 0 .   (22) Let   W 1 ( y ) =   M   − 1 1   ( y ), then we also have from (21) Π( y ) :=   −    ̇ W 1 ( y ) +   W 1 ( y )   ∂f 1  ∂y  ′  +   ∂f 1  ∂y W 1 ( y )   <   0   (23) with   W 1 ( y )   ≤   1  α   I . The combined system with state   x   := [ y ′ , z ′ ] ′   is:  ̇ x   =  [    ̇ y   ̇ z  ]  =  [ f 1 ( y, z )  ̄ Az  ]  +  [   0  ̄ B  ]  u   =:   f   ( x ) +   Bu  and the differential dynamics are (2) with  A ( y ) :=   ∂f  ∂x   =     ∂f  ∂y ∂f  ∂z  0    ̄ A     Let us consider a class of metrics parametrised by   β >   0:  W   ( y ) =  [ β ( y ) M 1 ( y ) − 1   0 0   P  ]  Noting that  ̇ W   =  [    ̇ W 1   0 0   0  ]  , BB ′   =  [ 0   0 0    ̄ B    ̄ B ′  ]  (24) we see that (6) reduces to      β ( y )Π( y )   β ( y ) W y  ∂f  ∂z β ( y )   ∂f  ∂z  ′  W y    ̄ A ′ P   +   P    ̄ A   −    ̄ B    ̄ B ′       <   0 .   (25) Now, since the top-left block is negative definite by as- sumption (23), negativity of the entire matrix is equiva- lent, by the Schur complement, to the following  ̄ A ′ P   +   P    ̄ A   −    ̄ B    ̄ B ′   −   β ( y ) W y  ∂f  ∂z   Π( y )   ∂f  ∂z  ′  W y   <   0 .   (26) Now, since  ̄ A ′ P   +   P    ̄ A   −    ̄ B    ̄ B ′   <   0 by assumption (22), this can clearly be satisfied for sufficiently small   β ( y )   >   0. Thus we have another class of systems for which the conditions of Theorem 2 are necessary and sufficient. This   situation is   clearly simpler than   Lyapunov-based control design: if a nonlinear system is Lyapunov stable at the origin, and driven by a stabilizable linear system, it is not so clear that the entire system is stabilizable. This is the fact that led to the development of backstepping [Krstic et al. 1995]. 10. CONCLUSIONS In this paper we have introduced the concept of   univer- sal stabilizability : the condition that   every   solution of a system is globally stabilizable. We have proven sufficient conditions in terms of the existence of a   control contraction metric . Unlike   control   Lyapunov   functions,   the   set   of   control contraction metrics for a given system can be parametrised
as   a   convex set   –   defined   by   pointwise   linear   matrix inequality constraints – and thus amenable to search via convex optimization methods. The conditions we give are necessary and sufficient for lin- ear systems and certain classes of interconnected nonlinear systems. Straightforward extensions allow one to construct convex upper bounds for a nonlinear quadratic regulator problem. REFERENCES Bao, D.D.W.,   Chern, S.S.,   and Shen,   Z.   (2000).   An introduction to Riemann-Finsler geometry , volume 200. Springer Verlag. Borg, G. (1960). A condition for the existence of orbitally stable solutions of dynamical systems.   Kungl. Tekn. H ̈ ogsk. Handl. , (153). Boyd, S., el Ghaoui, L., Feron, E., and Balakrishnan, V. (1994).   Linear matrix inequalities in system and control theory . Society for Industrial and Applied Mathematics (SIAM). Boykov,   Y.   and   Kolmogorov, V.   (2003).   Computing geodesics and minimal surfaces via graph cuts. In   Proc. IEEE International Conference on Computer Vision . Cloutier, J.R. (1997).   State-dependent riccati equation techniques: an overview.   In   Proc. American Control Conference . Demidovich, B. (1962). Dissipativity of nonlinear system of differential equations.   Ser. Mat. Mekh , 19–27. Deriaz, E. and Perrier, V. (2009). Orthogonal helmholtz decomposition in arbitrary dimension using divergence- free and curl-free wavelets.   Applied and Computational Harmonic Analysis , 26(2), 249–269. Diehl, M., Ferreau, H.J., and Haverbeke, N. (2009). Effi- cient numerical methods for nonlinear MPC and moving horizon estimation. In   Nonlinear Model Predictive Con- trol , 391–417. Springer. Forni,   F.   and   Sepulchre,   R.   (2012).   A   differential Lyapunov framework for contraction analysis.   arXiv preprint arXiv:1208.2943 . Hartman, P. (1961). On stability in the large for systems of ordinary differential equations.   Canadian Journal of Mathematics , 13(3), 480–492. Hartman, P. and Olech, C. (1962). On global asymptotic stability of solutions of differential equations.   Trans- actions of the American Mathematical Society , 104(1), 154–178. Hespanha, J.P. (2009).   Linear systems theory . Princeton university press. Isidori, A. (1995). Nonlinear control systems. Johansen, T.A. (2000).   Computation of Lyapunov func- tions for smooth nonlinear systems using convex opti- mization.   Automatica , 36(11), 1617–1626. Kokotovi ́ c, P. and Arcak, M. (2001). Constructive nonlin- ear control: a historical perspective.   Automatica , 37(5), 637–662. Krstic, M., Kanellakopoulos, I., and Kokotovic, P. (1995).  Nonlinear   and   adaptive   control   design ,   volume   222. Wiley. Lasserre, J.B., Henrion, D., Prieur, C., and Tr ́ elat, E. (2008). Nonlinear optimal control via occupation mea- sures and lmi-relaxations.   SIAM Journal on Control and Optimization , 47(4), 1643–1666. LaValle, S.M. (2006).   Planning algorithms .   Cambridge university press. Leonov, G., Burkin, I., and Shepeljavyi, A. (1996).   Fre- quency methods in oscillation theory . Kluwer Academic Publishers Dordrecht–Boston–London. Lewis, D. (1949).   Metric properties of differential equa- tions.   American Journal of Mathematics , 294–312. Lohmiller, W. and Slotine, J.J.E. (1998a).   Contraction analysis: A practical approach to nonlinear control ap- plications. In   Proc. IEEE International Conference on Control Applications . Lohmiller, W. and Slotine, J.J.E. (1998b). On contraction analysis for non-linear systems.   Automatica , 34(6), 683– 696. Lu, W.M. and Doyle, J.C. (1995).   H ∞   control of nonlinear systems: a convex characterization.   IEEE Transactions on Automatic Control , 40(9), 1668–1675. Manchester, I.R. and Slotine, J.J.E. (2013).   Transverse contraction criteria for existence, stability, and robust- ness of a limit cycle.   Systems & Control Letters , 62, 32–38. Megretski, A. and Rantzer, A. (1997).   System analysis via integral quadratic constraints.   Automatic Control, IEEE Transactions on , 42(6), 819–830. M ́ emoli, F. and Sapiro, G. (2001).   Fast computation of weighted distance functions and geodesics on im- plicit hyper-surfaces.   Journal of computational Physics , 173(2), 730–764. Parrilo, P.A. (2003).   Semidefinite programming relax- ations for semialgebraic problems.   Mathematical Pro- gramming , 96(2), 293–320. Parrilo, P.A. (2000).   Structured Semidefinite Programs and   Semialgebraic   Geometry   Methods   in   Robustness and Optimization . Ph.D. thesis, California Institute of Technology. Petersen, I., Ugrinovskii, V., and Savkin, A. (2000).   Robust control design using   H ∞   methods . Springer Verlag. Prajna, S., Parrilo, P.A., and Rantzer, A. (2004). Nonlin- ear control synthesis by convex optimization.   Automatic Control, IEEE Transactions on , 49(2), 310–314. Rantzer, A. (2001). A dual to lyapunov’s stability theorem.  Systems & Control Letters , 42(3), 161–168. Slotine,   J.J.E. and   Li,   W.   (1991).   Applied   nonlinear control . Prentice-Hall. Tedrake, R., Manchester, I.R., Tobenkin, M., and Roberts, J.W. (2010). LQR-Trees: Feedback motion planning via sums-of-squares verification.   The International Journal of Robotics Research , 29(8), 1038–1052. Uhlig, F. (1979).   A recurring theorem about pairs of quadratic   forms   and   extensions:   A   survey.   Linear Algebra and Its Applications , 25, 219–237. Vaidya, U., Mehta, P.G., and Shanbhag, U.V. (2010). Nonlinear stabilization via control lyapunov measure.  Automatic Control, IEEE Transactions on , 55(6), 1314– 1328. van der Schaft, A. (1999).   L2-Gain and Passivity in Nonlinear Control . Springer-Verlag New York, Inc. Yakubovich, V. (1971).   S-procedure in nonlinear control theory.   Vestnik Leningrad University , 1, 62–77.