A differential displacement, \( \delta{x}\) , is defined as an infinitesimal displacement at a fixed time as used in the calculus of variation [134, p. 107], and (1) yields the following differential dynamics:

where \( f({x}(t),t)\) is given in (1).

Suppose that a continuously differentiable function \( v\in\mathbb{R}_{\geq0}\mapsto\mathbb{R}\) satisfies the following differential inequality:

where \( \gamma \in \mathbb{R}_{>0}\) , \( c\in\mathbb{R}\) , and \( v_0\in\mathbb{R}\) . Then we have

Let \( \xi_0\) and \( \xi_1\) be the two arbitrary solutions of (1), and define the transformed squared length with \( M(x,t)\) of Theorem 1 as follows:

where \( x\) is a smooth path parameterized as \( x(\mu=0,t)=\xi_0(t)\) and \( x(\mu=1,t)=\xi_1(t)\) by \( \mu \in [0,1]\) . Also, define the path integral with the transformation \( \Theta(x,t)\) for \( M(x,t) = \Theta(x,t)^{\top}\Theta(x,t)\) as follows:

Then (25) and (26) are related as

where \( M(x,t)\succeq\underline{m}\mathrm{I}, \forall x,t\) for \( \exists\underline{m}\in\mathbb{R}_{>0}\) , and Theorem 1 can also be proven by using (25) and (26) as a Lyapunov-like function, resulting in incremental exponential stability of the system (1) (see Definition 2). Note that the shortest path integral \( V_\ell\) of (26) with a parameterized state \( x\) (i.e., \( \inf V_\ell=\sqrt{\inf V_{s\ell}}\) ) defines the Riemannian distance and the path integral of a minimizing geodesic [79].

If (29) is perturbed by \( {d}({x},t)\in \mathcal{L}_{pe}\) (i.e., \( \|(d)_{\tau}\|_{\mathcal{L}_p} < \infty\) for \( \tau\in\mathbb{R}_{\geq0}\) and \( p\in[1,\infty]\) , see Sec. 1.2.3) and Theorem 1 holds, then (29) is finite-gain \( \mathcal{L}_p\) stable with \( p\in [1,\infty]\) for an output \( {y}={h}({x},{d},t)\) with \( \int^{Y_1}_{Y_0}\|\delta {y}\|\leq \eta_0 \int^{\xi_1}_{\xi_0}\|\delta {x}\|+\eta_1\|{d}\|\) , \( \exists\eta_0,\eta_1\geq0\) , i.e., \( \forall \tau\in\mathbb{R}_{\geq 0}\)  [24]

where \( Y_0\) and \( Y_1\) denote the output trajectories of the original contracting system (1) and its perturbed system (29), respectively, \( \underline{m}\) is defined as \( {M}({x},t)\succeq \underline{m}{\mathrm{I}}\) , \( \forall x,t\) , as in (31), and \( V_\ell(0) = V_\ell(x(0),\delta x(0),0)\) for \( V_{\ell}\) in (26). Also, \( \zeta=1\) if \( p=\infty\) and \( \zeta=1/(\alpha p)^{1/p}\) if \( p\in[1,\infty)\) . The perturbed system (29) also exhibits ISS.

Consider the following hierarchically combined system of two contracting dynamics:

where \( F=F_{00}\) and \( F=F_{11}\) both satisfy the contraction condition (5), and suppose that it is subject to perturbation \( [{d}_0,{d}_1]^{\top}\) . Then the path length integral \( V_{\ell,i}(t)=\int^{\xi_1}_{\xi_0}\|\delta {z}_i\|\) with \( i=0,1\) between the original and perturbed dynamics trajectories, \( \xi_0\) and \( \xi_1\) , respectively, verifies [7]

where \( \alpha_i=\sup_{x,t}|\lambda_\mathrm{max}({F}_{ii})|, i=0,1\) . Hence, the error bounds of \( V_{\ell,0}(t)\) and \( V_{\ell,1}(t)\) can be obtained using Theorem 4 if \( \|{F}_{10}\|\) is bounded. In particular, if \( \|{\Theta}_i d_i\|\leq \sqrt{\overline{m}_i}\bar{d}_i\) , \( i=0,1\) , the relations (53) and (54) yield \( V_{\ell,0}(t) \leq V_{\ell,0}(0)e^{-\alpha_0t}+(\sqrt{\overline{m}_0}\bar{d}_0/{\alpha_0})(1-e^{-\alpha_0t})\) as in (32), and

where \( \bar{f}_{10}=\sup_{t\in\mathbb{R}_{\geq 0}}\|{F}_{10}\|\) and \( \bar{V}_{0}=V_{\ell,0}(0)+\sqrt{\overline{m}_0}\bar{d}_0/{\alpha_0}\) .

Also, similar to Theorem 6, we have \( \|(V_{\ell,0})_\tau\|_{\mathcal{L}_p} \leq V_{\ell,0}(0)\zeta_0+{\|({\Theta}_0{d}_0)_\tau \|_{\mathcal{L}_p}}/{\alpha_0}\) , and thus a hierarchical connection for finite-gain \( \mathcal{L}_p\) stability can be established as follows [24]:

where \( \zeta_i=1\) if \( p=\infty\) and \( \zeta_i=1/(\alpha_i p)^{1/p}\) if \( p\in[1,\infty)\) for \( i=0,1\) . By recursion, this result can be extended to an arbitrary number of hierarchically combined groups.

As demonstrated in Example 6, the Lagrangian virtual system (23) is contracting with respect to \( q\) , having \( q=\dot{\mathtt{q}}\) and \( q = \dot{\mathtt{q}}_r = \dot{\mathtt{q}}_d-\Lambda (\mathtt{q}-\mathtt{q}_d)\) as its particular solutions. Let \( q_0=q\) for such \( q\) . The virtual system of \( q_1\) which has \( [q_0^{\top},q_1^{\top}]^{\top}=[\dot{\mathtt{q}}^{\top},\mathtt{q}^{\top}]^{\top}\) and \( [\dot{\mathtt{q}}_r^{\top},\mathtt{q}_d^{\top}]^{\top}\) as its particular solutions is given as

resulting in \( \delta \dot{q}_1=\mathrm{I}\delta q_0-\Lambda\delta q_1\) . Since the virtual system of \( q_0\) is contracting as in Example 6 and the virtual system \( \delta \dot{q}_1=-\Lambda\delta q_1\) is contracting in \( q_1\) due to \( \Lambda \succ 0\) , (55) of Theorem 7 implies that the whole system (21) for \( [q_0^{\top},q_1^{\top}]^{\top}\) is hierarchically contracting and robust against perturbation in the sense of Theorem 7 (\( F_{10}=\mathrm{I}\) in this case). Also, see [40, 24] for the hierarchical multi-timescale separation of tracking and synchronization control for multiple Lagrangian systems.

Let \( {x}(k)=x_k\) and \( q(\mu,k)=q_k\) for any \( k\in\mathbb{N}\) for simplicity. If there exists a uniformly positive definite matrix \( {M}_k({x}_k,k)={{\Theta_k}}({x}_k,k)^{\top}{{\Theta_k}}({x}_k,k) \succ 0, \forall x_k,k\) , where \( {\Theta_k}\) defines a smooth coordinate transformation of \( \delta x_k\) , i.e., \( \delta{z}_k={\Theta_k}(x_k,k)\delta{x}_k\) , s.t. either of the following equivalent conditions holds for \( \exists \alpha \in (0,1)\) , \( \forall x_k,k\) :

then we have the following bound as long as we have \( \underline{m}\mathrm{I}\preceq M_x(x_k,k)\preceq \overline{m}\mathrm{I}, \forall x_k,k\) , as in (31):

where \( V_{\ell}(k)=\int_{\xi_{0}}^{\xi_{1}}\|\Theta_k(q_k,k)\delta q_k\|\) as in (26) for the unperturbed trajectory \( \xi_{0}\) , perturbed trajectory \( \xi_1\) , and virtual state \( q_k=q(k)\) given in (59).

Consider a general feedback controller \( u\) defined as \( u = k(x,x_d,u_d,t)\) with \( k(x_d,x_d,u_d,t)=u_d\) , where \( k:\mathbb{R}^n\times\mathbb{R}^n\times\mathbb{R}^m\times\mathbb{R}_{\geq 0}\mapsto\mathbb{R}^{m}\) . If \( k\) is piecewise continuously differentiable, then \( \exists K:\mathbb{R}^n\times\mathbb{R}^n\times\mathbb{R}^m\times\mathbb{R}_{\geq 0}\mapsto\mathbb{R}^{m\times n}\) s.t. \( u = k(x,x_d,u_d,t) = u_d-K(x,x_d,u_d,t)(x-x_d)\) .

Lemma 4 implies that designing optimal \( k\) of \( u = k(x,x_d,u_d,t)\) reduces to designing the optimal gain \( K(x,x_d,u_d,t)\) of \( u=u_d-K(x,x_d,u_d,t)(x-x_d)\) . We could also generalize this idea further using the CCM-based differential feedback controller \( \delta u = k(x,\delta x,u,t)\)  [79, 80, 81, 82, 83, 33] (see Theorem 18).

Let \( \beta\) be defined as \( \beta = 0\) for deterministic systems (63) and

for stochastic systems (64), respectively, where \( \bar{g}_c\) is given in (64), \( L_m\) is the Lipschitz constant of \( \partial M/\partial x_i\) for \( M\) of (77), and \( \alpha_G\in\mathbb{R}_{>0}\) is an arbitrary constant as in Theorem 5. Also, let \( W = M(x,x_d,u_d,t)^{-1}\) (or \( W=M(x,t)^{-1}\) , see Theorem 10), \( \bar{W}=\nu W\) , and \( \nu=\overline{m}\) . Then the following three matrix inequalities are equivalent:

where \( \zeta\) is as defined in (81). For stochastic systems with \( \beta=\alpha_s>0\) , these inequalities are also equivalent to

Note that \( \nu\) and \( \bar{W}\) are required for (86) and (87) and the arguments \( (x,x_d,u_d,t)\) for each matrix are suppressed for notational simplicity.

Furthermore, under these equivalent contraction conditions, Theorems 4 and 5 hold for the virtual systems (79) and (80), respectively. In particular, if \( \underline{m}\mathrm{I} \preceq M\preceq \overline{m}\mathrm{I}\) of (31) holds, or equivalently

holds for \( \chi=\overline{m}/\underline{m}\) , then we have the following bounds:

where \( V_{s\ell}=\int^x_{x_d}\delta q^{\top}M\delta q\) and \( V_{\ell}=\int^x_{x_d}\|\Theta\delta q\|\) are as given in Theorem 3 with \( M=\Theta^{\top}\Theta\) , the disturbance bounds \( \bar{d}_c\) and \( \bar{g}_c\) are given in (63) and (64), respectively, and \( C_C = \bar{g}_c^2({2}{\alpha_G}^{-1}+1)\) . Note that for stochastic systems, the probability that \( \|x-x_d\|\) is greater than or equal to \( \varepsilon\in\mathbb{R}_{> 0}\) is given as

Let \( (\bar{x},\bar{u})\) be a fixed point selected arbitrarily in \( \mathbb{R}^n\times\mathbb{R}^m\) , e.g.{}, \( (\bar{x},\bar{u})=(0,0)\) , and let \( A(x,t)\) be an SDC matrix constructed with \( ({s},\bar{s},\bar{u})=(x,\bar{x},\bar{u})\) in Lemma 3, i.e.,

Suppose that the contraction metric of Theorem 9 is designed by \( M(x,t)\) with \( A\) of (92), independently of the target trajectory \( (x_d,u_d)\) , and that the systems (63) and (64) are controlled by

with such \( M(x,t)\) , where \( R(x,t)\succ0\) is a weight matrix on \( u\) . If the function \( \phi(x,x_d,u_d,t)=A(\varrho,x,t)(x_d-\bar{x})+B(x,t)(u_d-\bar{u})\) is Lipschitz in \( x\) with its Lipschitz constant \( \bar{L}\) , then Theorem 9 still holds with \( \alpha\) of the conditions (84)–(87) replaced by \( \alpha+\bar{L}\sqrt{\overline{m}/\underline{m}}\) . The same argument holds for state estimation of Theorem 15 to be discussed in Sec. 2.3.2.

As demonstrated in [13], we could directly use the extra term \( 2\delta q^{\top}M({\partial\phi}/{\partial q})\delta q\) of (95) in (84)–(87) without upper-bounding it, although now the constraints of Theorem 10 depend on \( (x,q,t)\) instead of \( (x,t)\) . Also, the following two inequalities given in [13] with \( \bar{\gamma}=\nu \gamma\) , \( \gamma \in \mathbb{R}_{\geq 0}\) :

are combined as one LMI (87) in Theorems 9 and 10.

The inequalities in Theorem 9 can be interpreted as in the Riccati inequality in \( \mathcal{H}_{\infty}\) control. Consider the following system:

where \( A\in\mathbb{R}^{n\times n}\) , \( B_u\in\mathbb{R}^{n\times m}\) , \( B_w\in\mathbb{R}^{n\times w}\) , and \( C_z\in\mathbb{R}^{o\times n}\) are constant matrices, \( w\in\mathbb{R}^w\) is an exogenous input, and \( z\in\mathbb{R}^o\) is a system output. As shown in [25] and [21, p. 109], there exists a state feedback gain \( K=R^{-1}B_u^{\top}P\) such that the \( \mathcal{L}_2\) gain of the closed-loop system (97), \( \sup_{\|w\|\neq 0}{\|z\|}/{\|w\|}\) , is less than or equal to \( \gamma\) if

has a solution \( P \succ 0\) , where \( R\succ 0\) is a constant weight matrix on the input \( u\) . If we select \( B_w\) and \( C_z\) to have \( B_wB_w^{\top}\succeq (P^{-1})^2\) and \( C_z^{\top}C_z \succeq 2\alpha P\) for some \( \alpha >0\) , the contraction condition (85) in Theorem 9 can be satisfied with \( M=P\) , \( B = B_u\) , and \( \beta = 1/\gamma^2\) due to (98).

Let \( \mathcal{A}= A-BR^{-1}B^{\top}M\) , \( \mathcal{B}=M^{-1}\) , \( C=\sqrt{2\alpha}\Theta\) , and \( \mathcal{D}=0\) in (99), where \( M=\Theta^{\top}\Theta\) , and the other variables are as defined in Theorem 9. Then (100) with \( \mathcal{P}=M\) and \( \gamma=1/\sqrt{\beta}\) is equivalent to (87), and thus (87) implies \( \mathcal{L}_2\) gain stability of (99) with its \( \mathcal{L}_2\) gain less than or equal to \( \gamma\) .

If (104) holds, the LMI (100) for the bounded real lemma holds with \( \mathcal{P}=\gamma \mathcal{Q}\) , i.e., the system (99) is \( \mathcal{L}_2\) gain stable with its \( \mathcal{L}_2\) gain less than or equal to \( \gamma\) . Thus, for systems with \( \mathcal{A}\) , \( \mathcal{B}\) , \( \mathcal{C}\) , and \( \mathcal{D}\) defined in Theorem 11, the condition (104) guarantees the contraction condition (87) of Theorem 9.

Theorems 11 and 12 imply the following statements.

• If (100) holds, the system is finite-gain \( \mathcal{L}_2\) stable with \( \gamma\) as its \( \mathcal{L}_2\) gain.
• If (104) holds, the system is finite-gain \( \mathcal{L}_2\) stable with \( \gamma\) as its \( \mathcal{L}_2\) gain.
• (102) is equivalent to (104), and to the output strict passivity condition with dissipation \( 1/\gamma\)  [3, p. 231].
• If (104) holds (KYP lemma), (100) holds (bounded real lemma). The converse is not necessarily true.

If one of the matrix inequalities of Theorem 9 holds, we have the following bound for \( \chi={\overline{m}}/\underline{{m}}\) :

where \( c_0(\alpha,\alpha_G)=\bar{d}_c /\alpha\) for deterministic systems and \( c_0(\alpha,\alpha_G)=\bar{g}_c\sqrt{(2\alpha_G^{-1}+1)/(2\alpha)}\) for stochastic systems, with the variables given in Theorem 9.

Suppose that \( \alpha\) , \( \alpha_G\) , \( \bar{d}_c\) , and \( \bar{g}_c\) in (105) and the Lipschitz constant \( L_m\) of \( \partial M/\partial x_i\) in Theorem 9 are given. If the pair \( (A,B)\) is uniformly controllable, the non-convex optimization problem of minimizing the upper bound (105) is equivalent to the following convex optimization problem, with the convex contraction constraint (86) or (87) of Theorem 9 and \( \mathrm{I}\preceq \bar{W}\preceq\chi\mathrm{I}\) of (88):

where \( W = M(x,x_d,u_d,t)^{-1}\) for Theorem 9 or \( W=M(x,t)^{-1}\) for Theorem 10, \( \bar{W}=\nu W\) , \( \nu=\overline{m}\) , \( \chi=\overline{m}/\underline{m}\) , \( c_0\) is as defined in (105) of Theorem 13, \( c_1,c_2\in\mathbb{R}_{\geq 0}\) , and \( P\) is a performance-based convex cost function.

The weight \( c_1\) for \( \nu\) can be viewed as a penalty on the \( 2\) -norm of the feedback control gain \( K\) of \( u=u_d-K(x-x_d)\) in (77). Using non-zero \( c_1\) and \( c_2\) thus enables finding contraction metrics optimal in a different sense. Furthermore, the coefficients of the SDC parameterizations \( \varrho\) in Lemma 3, i.e., \( A = \sum\varrho_iA_i\) in (71), can also be treated as decision variables by convex relaxation, thereby adding design flexibility to mitigate the effects of external disturbances while verifying the system controllability.

The weights \( c_0\) and \( c_1\) of the CV-STEM control of Theorem 14 establish an analogous trade-off to the case of the LQR with the cost weight matrices of \( Q\) for state and \( R\) for control, since the upper bound of the steady-state tracking error (105) is proportional to \( \chi\) , and an upper bound of \( \|K\|\) for the control gain \( K\) in (77), \( u=u_d-K(x-x_d)\) , is proportional to \( \nu\) . In particular [14, 15],

• if \( c_1\) is much greater than \( c_0\) in (106) of Theorem 14, we get smaller control effort but with a large steady-state tracking error, and
• if \( c_1\) is much smaller than \( c_0\) in (106) of Theorem 14, we get a smaller steady-state state tracking error but with larger control effort.

Suppose the CV-STEM control of Theorem 14 is designed with its contraction condition replaced by the following set of convex constraints along with \( \mathrm{I}\preceq \bar{W}\preceq\chi\mathrm{I}\) of (88):

where \( B_{\bot}(x,t)\) is a matrix whose columns span the cokernel of \( B(x,t)\) defined as \( \mathrm{coker}(B)= \{a\in\mathbb{R}^n|B^{\top}a=0\}\) satisfying \( B^{\top}B_{\bot} = 0\) , \( b_i(x,t)\) is the \( i\) th column of \( B(x,t)\) , \( \bar{W}(x,t)=\nu W(x,t)\) , \( W(x,t) = M(x,t)^{-1}\) for \( M(x,t)\) that defines a contraction metric, \( \nu=\overline{m}\) , \( \chi=\overline{m}/\underline{m}\) , and \( \partial_{p} F = \sum_{k=1}^n(\partial F/\partial x_k)p_k\) for \( p(x,t)\in\mathbb{R}^n\) and \( F(x,t) \in \mathbb{R}^{n\times n}\) . Then the controlled system (63) with \( d_c=0\) , i.e., \( \dot{x}=f(x,t)+B(x,t)u\) , is

1. universally exponentially open-loop controllable
2. universally exponentially stabilizable via sampled-data feedback with arbitrary sample times
3. universally exponentially stabilizable using continuous feedback defined almost everywhere, and everywhere in a neighborhood of the target trajectory

all with rate \( \alpha\) and overshoot \( R = \sqrt{\overline{m}/\underline{m}}\) for \( \underline{m}\) and \( \overline{m}\) given in \( \underline{m}\mathrm{I}\preceq M\preceq\overline{m}\mathrm{I}\) of (31). Such a positive definite matrix \( M(x,t)\) defines a CCM.

Given a CCM, there exists a differential feedback controller \( \delta u = k(x,\delta x,u,t)\) that stabilizes the following differential dynamics of (63) with \( d_c=0\) along all solutions (i.e., the closed-loop dynamics is contracting as in Definition 3):

where \( A=\partial f/\partial x+\sum_{i=1}^m(\partial b_i/\partial x)u_i\) and \( \delta u\) is a tangent vector to a smooth path of controls at \( u\) . Furthermore, it can be computed as follows [79, 80, 33]:

where \( \partial_{\mu}\gamma=\partial \gamma/\partial \mu\) , \( \gamma\) is the minimizing geodesic with \( \gamma(0,t) = x_d\) and \( \gamma(1,t) = x\) for \( \mu\in[0,1]\) , \( (x_d,u_d)\) is a given target trajectory in \( \dot{x}_d=f(x_d,t)+B(x_d,t)u_d\) of (65), and the computation of \( k\) is given in [79, 80] (see Remark 11). Furthermore, Theorem 4 for deterministic disturbance rejection still holds and the CCM controller (129) thus possesses the same sense of oplimality as for the CV-STEM control in Theorem 14.

Let us consider a case where \( k\) of the differential feedback controller does not depend on \( u\) , and is defined explicitly by \( M(x,t) \succ 0\) as follows [119]:

where \( R(x,t) \succ 0\) is a given weight matrix and \( q\) is a smooth path that connects \( x\) to \( x_d\) . Since (131) yields \( \delta u=-R(x,t)^{-1}B(x,t)^{\top}M(x,t)\delta x\) , the contraction conditions (126) and (127) could be simplified as

yielding the convex optimization-based control synthesis algorithm of Theorem 18 independently of \( (x_d,u_d)\) similar to that of Theorem 14 (see [119] for details).

If \( B\) is of the form \( [0,\mathrm{I}]^{\top}\) for the zero matrix \( 0\in\mathbb{R}^{n_1\times m}\) and identity matrix \( \mathrm{I}\in\mathbb{R}^{n_2\times m}\) with \( n=n_1+n_2\) , the condition (133) says that \( M\) should not depend on the last \( n_2\) state variables [79].

We could consider stochastic perturbation in Theorem 18 using Theorem 5, even with the differential control law of the form (129) or (131) as demonstrated in [119]. Also, although the relation (127) or (133) is not included as a constraint in Theorem 14 for simplicity of discussion, the dependence of \( \dot{\bar{W}}\) on \( u\) in Theorem 14 can be removed by using it in a similar way to [119].

As stated in (129), the computation of the differential feedback gain \( k(x,\delta x,u,t)\) and minimizing geodesics \( \gamma\) is elaborated in [79, 80]. For example, if \( M\) is state-independent, then geodesics are just straight lines.

Let us consider the problem of learning a computationally-expensive (or unknown) feedback control law \( u^*(x,x_d,u_d,t)\) that tracks a target trajectory \( (x_d,u_d)\) given by \( \dot{x}_d=f(x_d,t)+B(x_d,t)u_d\) as in (65), assuming that \( f\) and \( B\) are known. If the dynamics is perturbed by \( d_c(x,t)\) as in (63), we have that

where \( u_L(x,x_d,u_d,t)\) denotes a learning-based control law that models \( u^*(x,x_d,u_d,t)\) . As long as \( u^*\) renders the closed-loop dynamics \( \dot{x}=f(x,t)+B(x,t)u^*(x,x_d,u_d,t)\) contracting with respect to \( x\)  [13, 79, 14, 80], we can define the functions of (134) as follows:

with \( d(\mu,\varpi,t)=\mu d_c(x,t)\) , where \( \varpi=[x^{\top},x_d^{\top},u_d^{\top}]^{\top}\) in this case. It can be easily verified that (134) indeed has \( q(\mu=0,t)=x_d\) and \( q(\mu=1,t)=x\) as particular solutions if \( u^*(x_d,x_d,u_d,t)=u_d\) . Similarly, we can use (135) with \( G=\mu G_c(x,t)\) if the dynamics is stochastically perturbed by \( G_c(x,t)d\mathscr{W}\) as in (64).

The learning objective here is to make \( \|u_L(x,x_d,u_d,t)-u^*(x,x_d,u_d,t)\|\) as small as possible, which aligns with the aforementioned objective in (136). Note that if \( \|B\|\) is bounded in the compact set \( \mathcal{S}\) of (136), we can bound \( \|\Delta_L\|\) of (139) with \( \epsilon_{\ell0}\neq 0\) and \( \epsilon_{\ell1}=0\) as to be explained in Remark 12. See Theorems 22, 24, and 27 for details.

Next, let us consider the problem of learning a computationally-expensive (or unknown) state estimator, approximating its estimation gain \( L(\hat{x},t)\) by \( L_L(\hat{x},t)\) . If there is no disturbance in (107) and (108), we have that

where \( \dot{x}=f(x,t)\) is the true system, \( y=h(x,t)\) is the measurement, and we assume that \( f\) and \( h\) are known. If \( L\) is designed to render \( \dot{q}=f(q,t)+L(\hat{x},t)(h(x,t)-h(q,t))\) contracting with respect to \( q\)  [14, 13, 16, 45], we could define the functions of (134) and (135) as follows:

where \( \varpi=[x^{\top},\hat{x}^{\top}]^{\top}\) . It can be seen that (134) and (135) with the relations (141) and (142) indeed has \( q(\mu=0,t)=x\) and \( q(\mu=1,t)=\hat{x}\) as its particular solutions when perturbed by deterministic and stochastic disturbances as in (107) and (108), respectively. We can bound \( \|\Delta_L\|\) of (142) in the compact set of (136) with \( \epsilon_{\ell0}=0\) and \( \epsilon_{\ell1}\neq0\) if \( h\) is Lipschitz, and with \( \epsilon_{\ell0}\neq0\) and \( \epsilon_{\ell1}=0\) if \( h\) is bounded in \( \mathcal{S}\) , using the techniques to be discussed in Remark 12. See Theorems 23 and 24 for details.

We can use the problem (136) with the systems (134) and (135) also if \( f_{\mathrm{true}}\) of the underlying dynamics \( \dot{x}^* = f_{\mathrm{true}}(x^*,t)\) is unknown and learned by \( \dot{x} = f_L(x,t)\) . Since we have

we could define \( g \) and \( \Delta_L\) of (134) and (135) as follows, to have \( q(\mu=0,t)=x\) and \( q(\mu=1,t)=x^*\) as its particular solutions:

where \( \varpi=x^*\) , as long as \( \dot{x} = f_L(x,t)\) is contracting with respect to \( x\)  [81, 133]. Since \( \Delta_L\) of (144) is the learning error itself, \( \|\Delta_L\|\) can be bounded in \( \mathcal{S}\) using the techniques of Remark 12. See Theorems 30 and 31 for details.

As seen in Examples 17–18, \( \Delta_L\) of (136) is typically given by a learning error, \( \phi(z)-\varphi(z)\) , multiplied by a bounded or Lipschitz continuous function, where \( z=(\varpi,t)\in\mathcal{S}\) for a compact set \( \mathcal{S}\) , and \( \phi\) is the true computationally-expensive/unknown function to be learned by \( \varphi\) as in (139), (142), (144), etc.

Let \( \mathcal{T}\) be a set of training data \( \{(z_i,\phi(z_i))\}_{i=1}^N\) sampled in \( \mathcal{S}\) . For systems with the true function \( \phi\) and its approximation \( \varphi\) being Lipschitz (e.g.{}, by spectral normalization [132] to be discussed in Definition 6 or by using contraction theory as in [62, 63, 64]), we can analytically find a bound for \( \|\phi(z)-\varphi(z)\|\) as follows if target data samples \( z\) are in \( B(r) = \{z\in\mathcal{S}|\sup_{z'\in\mathcal{T}}\|z-z'\|\leq r\}\) for some \( r \in\mathbb{R}_{\geq0}\) :

where \( L_{\phi}\) and \( L_{\varphi}\) are the Lipschitz constants of \( \phi\) and \( \varphi\) , respectively. The term \( \sup_{z'\in \mathcal{T}}\|\phi(z')-\varphi(z')\|\) can then be bounded by a constant, e.g.{}, by using a deep robust regression model as proven in [154, 155] with spectral normalization, under standard training data distribution assumptions.

Deep Neural Networks (DNNs) have been shown to generalize well to the set of unseen events that are from almost the same distribution as their training set [156, 157, 61, 158, 159], and consequently, obtaining a tighter and more general upper bound for the learning error as in (145) has been an active field of research [160, 161, 155, 154]. Thus, the condition (136) has become a common assumption in analyzing the performance of learning-based and data-driven control techniques [61, 162, 163, 164, 105, 104, 119, 165, 166].

Suppose that \( \Delta_L\) of (134) is learned to satisfy (136), and that \( g (q,\varpi,t)\) of (134) is Lipschitz with respect to \( q\) with its 2-norm Lipschitz constant \( L_{ g }\in\mathbb{R}_{\geq 0}\) , i.e.,

where \( q,q' \in \mathbb{R}^n\) and \( \mathcal{S}\) is the compact set of (136). Then we have the following upper bound for all \( (\varpi,t) \in \mathcal{S}\) :

where \( \mathtt{e}(t) = \xi_1(t)-\xi_0(t)\) , \( \xi_0(t)=q(0,t)\) and \( \xi_1(t)=q(1,t)\) for \( q\) of (134), and \( \bar{d}=\sup_{\mu,\varpi,t}\|\partial d/\partial \mu\|\) for \( d\) of (136).

Let us consider the virtual system of a smooth path \( q(\mu,t)\) in (134) and suppose that \( \Delta_L\) of (134) is learned to satisfy (136). If the condition \( \underline{m}\mathrm{I}\preceq M\preceq\overline{m}\mathrm{I}\) of (31) holds and the system (134) with \( \Delta_L=0\) and \( d=0\) is contracting, i.e.,

of Theorem 1 or 2 holds for \( M\succ 0\) that defines a contraction metric with the contraction rate \( \alpha\) , and if the learning error \( \epsilon_{\ell1}\) of (136) is sufficiently small to satisfy

then we have the following bound for all \( (\varpi,t) \in \mathcal{S}\) :

where \( \mathtt{e}(t)=\xi_1(t)-\xi_0(t)\) with \( \xi_0(t)=q(0,t)\) and \( \xi_1(t)=q(1,t)\) for \( q\) of (134), \( \bar{d}=\sup_{\mu,\varpi,t}\|\partial d/\partial \mu\|\) as in (134), \( \mathcal{S}\) is the compact set of Theorem 19, \( \epsilon_{\ell0}\) and \( \epsilon_{\ell1}\) are the learning errors of (136), and \( V_{\ell}(t) = \int_{\xi_0}^{\xi_1} \|\Theta(q(t),t)\delta q(t)\|\) for \( M=\Theta^{\top}\Theta\) as in Theorem 3.

Consider the virtual system of a smooth path \( q(\mu,t)\) in (135) and suppose \( \Delta_L\) of (135) is learned to satisfy (136). If \( \underline{m}\mathrm{I}\preceq M\preceq\overline{m}\mathrm{I}\) of (31) holds and the system (135) is contracting, i.e.,

of Theorem 5 holds, and if the learning error \( \epsilon_{\ell1}\) of (136) and an arbitrary constant \( \alpha_d \in \mathbb{R}_{>0}\) (see (156)) are selected to satisfy

then we have the following bound for all \( (\varpi,t) \in \mathcal{S}\) :

where \( \mathtt{e}(t)=\xi_1(t)-\xi_0(t)\) with \( \xi_0(t)=q(0,t)\) and \( \xi_1(t)=q(1,t)\) for \( q\) of (135), \( \mathcal{S}\) is the compact set given in Theorem 19, \( \epsilon_{\ell0}\) and \( \epsilon_{\ell1}\) are the learning errors of (136), \( V_{s\ell}(t) = \int_{\xi_0}^{\xi_1} \|\Theta(q(t),t)\delta q(t)\|^2\) for \( M=\Theta^{\top}\Theta\) as given in Theorem 3, and \( C = \bar{g}^2({2}{\alpha_G}^{-1}+1)+\epsilon_{\ell0}^2\alpha_d^{-1}\) for \( \bar{g}=\sup_{\mu,\varpi,t}\|\partial G/\partial \mu\|_F\) of (135) with an arbitrary constant \( \alpha_G \in \mathbb{R}_{>0}\) as in Theorem 5. Furthermore, the probability that \( \|\mathtt{e}\|\) is greater than or equal to \( \varepsilon\in\mathbb{R}_{> 0}\) is given as

Let \( \mathcal{M}\) define the NCM in Definition 5, and let \( u^*=u_d-R^{-1}B^{\top}M(x-x_d)\) for \( M\) and \( R\) given in Theorem 14. Suppose that the systems (157) and (158) are controlled by \( u_L\) , which computes the CV-STEM controller \( u^*\) of Theorem 14 replacing \( M\) by \( \mathcal{M}\) , i.e.,

for \( \mathtt{e}=x-x_d\) . Note that we could use the differential feedback control \( u=u_d+\int^1_0kd\mu\) of Theorem 18 for \( u^*\) and \( u_L\) . Define \( \Delta_L\) of Theorems 20 and 21 as follows:

where we use \( q(0,t)=\xi_0(t)=x_d(t)\) , \( q(1,t)=\xi_1(t)=x(t)\) , and \( \varpi=[x^{\top},x_d^{\top},u_d^{\top}]^{\top}\) in the learning-based control formulation (134) and (135). If the NCM is learned to satisfy the learning error assumptions of Theorems 20 and 21 for \( \Delta_L\) given by (164), then (150), (154), and (155) of Theorems 20 and 21 hold as follows:

where \( \mathtt{e}=x-x_d\) , \( C_C = \bar{g}_c^2({2}{\alpha_G}^{-1}+1)+\epsilon_{\ell0}^2\alpha_d^{-1}\) , the disturbance upper bounds \( \bar{d}_c\) and \( \bar{g}_c\) are given in (157) and (158), respectively, and the other variables are defined in Theorems 14, 20 and 21.

Let \( \mathcal{M}\) define the NCM in Definition 5, and let \( L(\hat{x},t)=M\bar{C}^{\top}R^{-1}\) for \( M\) , \( \bar{C}\) and \( R\) given in Theorem 17. Suppose that the systems (160) and (161) are estimated with an estimator gain \( L_L\) , which computes \( L\) of the CV-STEM estimator \( \dot{\hat{x}} = f(\hat{x},t)+L(\hat{x},t)(y-h(\hat{x},t))\) in Theorem 17 replacing \( M\) by \( \mathcal{M}\) , i.e.,

Define \( \Delta_L\) of Theorems 20 and 21 as follows:

where we use \( q(0,t)=\xi_0(t)=x(t)\) , \( q(1,t)=\xi_1(t)=\hat{x}(t)\) , and \( \varpi=[x^{\top},\hat{x}^{\top}]^{\top}\) in the learning-based control formulation (134) and (135). If the NCM is learned to satisfy the learning error assumptions of Theorems 20 and 21 for \( \Delta_L\) given by (171), then (150), (154), and (155) of Theorems 20 and 21 hold as follows:

where \( \mathtt{e}=\hat{x}-x\) , \( \bar{d}_{a} = \epsilon_{\ell0}+\bar{d}_{e0}\) , \( \bar{d}_{b}=\bar{\rho}\bar{c}\bar{d}_{e1}\) , \( C_E = (\bar{g}_{e0}^2+\bar{\rho}^2\bar{c}^2\bar{g}_{e1}^2\overline{m}^2)({2}{\alpha_G}^{-1}+1)+\epsilon_{\ell0}^2\alpha_d^{-1}\) , the disturbance upper bounds \( \bar{d}_{e0}\) \( \bar{d}_{e1}\) , \( \bar{g}_{e0}\) , and \( \bar{g}_{e1}\) are given in (160) and (161), and the other variables are defined in Theorems 17, 20 and 21.

As discussed in Examples 17 and 18 and in Remark 12, we have \( \epsilon_{\ell1}=0\) for the learning error \( \epsilon_{\ell1}\) of (136) as long as \( B\) and \( h\) are bounded a compact set by the definition of \( \Delta_L\) in (164) and (171). If \( h\) and \( u^*\) are Lipschitz with respect to \( x\) in the compact set, we have non-zero \( \epsilon_{\ell1}\) with \( \epsilon_{\ell0}=0\) in (136) (see Theorem 24).

Let \( \mathcal{M}\) define the NCM in Definition 5, and let \( \mathcal{S}_s\subseteq\mathbb{R}^n\) , \( \mathcal{S}_u\subseteq\mathbb{R}^m\) , and \( t\in\mathcal{S}_t\subseteq\mathbb{R}_{\geq 0}\) be some compact sets. Suppose that the systems (157)–(161) are controlled by (163) and estimated by (169), respectively, as in Theorems 22 and 23. Suppose also that \( \exists \bar{b},\bar{c},\bar{\rho}\in\mathbb{R}_{\geq 0}\) s.t. \( \|B(x,t)\|\leq\bar{b}\) , \( \|C(\varrho_c,x,\hat{x},t)\|\leq\bar{c}\) , and \( \|R^{-1}\|\leq\bar{\rho}, \forall x,\hat{x}\in\mathcal{S}_s\) and \( t\in\mathcal{S}_t\) , for \( B\) in (157), \( C\) in Theorem 17, and \( R\) in (163) and (169). If we have for all \( x,x_d,\hat{x}\in\mathcal{S}_s\) , \( u_d\in\mathcal{S}_u\) , and \( t\in\mathcal{S}_t\) that

where \( M(x,x_d,u_d,t)\) and \( M(\hat{x},t)\) (or \( M(x,t)\) and \( M(\hat{x},t)\) , see Theorem 10) are the CV-STEM solutions of Theorems 14 and 17 with \( \exists\epsilon_{\ell}\in\mathbb{R}_{\geq 0}\) being the learning error, then the bounds (165)–(167) of Theorem 22 and (172)–(174) of Theorem 23 hold with \( \epsilon_{\ell0} = 0\) and \( \epsilon_{\ell1} = \bar{\rho}\bar{b}^2\epsilon_{\ell}\) for control, and \( \epsilon_{\ell0} = 0\) and \( \epsilon_{\ell1} = \bar{\rho}\bar{c}^2\epsilon_{\ell}\) for estimation, as long as \( \epsilon_{\ell}\) is sufficiently small to satisfy the conditions (149) and (153) of Theorems 20 and 21 for deterministic and stochastic systems, respectively.

The left-hand side of Fig. 4 shows the state estimation errors of the following Lorenz oscillator perturbed by process noise \( d_0\) and measurement noise \( d_1\)  [14]:

where \( x = [x_1,x_2,x_3]^{\top}\) , \( \sigma = 10\) , \( \beta = 8/3\) , and \( \rho = 28\) . It can be seen that the steady-state estimation errors of the NCM and CV-STEM are below the optimal steady-state upper bound of Theorem 17 (dotted line in the left-hand side of Fig. 4) while the EKF has a larger error compared to the other two. In this example, training data is sampled along \( 100\) trajectories with uniformly distributed initial conditions (\( -10 \leq x_i \leq 10, i = 1,2,3\) ) using Theorem 17, and then modeled by a DNN as in Definition 5. Additional implementation and network training details can be found in [14].

Consider a robust feedback control problem of the planar spacecraft perturbed by deterministic disturbances, the unperturbed dynamics of which is given as follows [167]:

where \( x = [p_x,p_y,\phi,\dot{p}_x,\dot{p}_y,\dot{\phi}]^T\) with \( p_x\) , \( p_y\) , and \( \phi\) being the horizontal coordinate, vertical coordinate, and yaw angle of the spacecraft, respectively, with the other variables defined as \( m = \text{mass of spacecraft}\) , \( I = \text{yaw moment of inertia}\) , \( u = \text{thruster force vector}\) , \( \ell=\text{half-depth of spacecraft}\) , and \( b=\text{half-width of spacecraft}\) (see Fig. 8 of [167] for details).

As shown in the right-hand side of Fig. 4, the NCM keeps their trajectories within the bounded error tube of Theorem 14 (shaded region) around the target trajectory (dotted line) avoiding collision with the circular obstacles, even under the presence of disturbances, whilst requiring much smaller computational cost than the CV-STEM. Data sampling and the NCM training are performed as in Example 20 [14].

We could also simultaneously synthesize a controller \( u\) and the NCM, and Theorem 22 can provide formal robustness and stability guarantees even for such cases [120, 119]. In Sec. 3.3, we generalize the NCM to learn contraction theory-based robust control using a DNN only taking \( x\) , \( t\) , and a vector containing local environment information as its inputs, to avoid the online computation of \( x_d\) for the sake of automatic guidance and control implementable in real-time.

Suppose (157) and (158) are controlled by \( u = u_d-K^*(x,x_d,u_d,t)\mathtt{e}\) , where \( \mathtt{e}=x-x_d\) , and that \( K^*\) is designed by the following convex optimization for given \( (x,x_d,u_d,t)\) :

where \( \mathcal{M}\) defines the NCM of Definition 5 that models \( M\) of the CV-STEM contraction metric in Theorem 14 as in (177), \( \mathcal{C}_K\) is a convex set containing admissible \( K\) , \( \alpha\) is the contraction rate, and \( \beta\) is as defined in theorem 9 (i.e., \( \beta=0\) for deterministic systems and \( \beta = \alpha_s\) for stochastic systems). Then Theorems 4 and 5 still hold for the Lyapunov function defined in (183), yielding the following bounds:

where \( \bar{p}=\sup_{x,x_d,u_d,t}p\) , \( \underline{m}\mathrm{I}\preceq\mathcal{M}\preceq\overline{m}\mathrm{I}\) as in (31), \( \alpha_d\in\mathbb{R}_{>0}\) is an arbitrary positive constant (see (189)) selected to satisfy \( \alpha_p=\alpha-\alpha_d/2>0\) , \( C_C = \bar{g}_c^2({2}{\alpha_G}^{-1}+1)\) , \( \alpha_G\in\mathbb{R}_{>0}\) is an arbitrary constant as in Theorem 5, \( \varepsilon\in\mathbb{R}_{\geq0}\) , and the disturbance terms \( \bar{d}_c\) and \( \bar{g}_c\) are given in (157) and (158), respectively. Note that the problem (184) is feasible due to the constraint relaxation variable \( p\) .

Furthermore, if \( \epsilon_{\ell} = 0\) in (177), (184) with \( p=0\) is always feasible, and the optimal feedback gain \( K^*\) minimizes its Frobenius norm under the contraction constraint (185).

Suppose \( u=u^*\) of (157) is designed by the following convex optimization for given \( (x,x_d,u_d,t)\) :

where \( \mathcal{M}\) defines the NCM of Definition 5 that models the CCM of Theorem 18, and \( \mathcal{C}_{u}\) is a convex set containing admissible \( u\) . Then the bound (186) of Theorem 25 holds.

Furthermore, if \( \epsilon_{\ell}=0\) for the NCM learning error \( \epsilon_{\ell}\) in (177), then the problem (192) with \( p=0\) is always feasible, and \( u^*\) minimizes the deviation of \( u\) from \( u_d\) under the contraction constraint (193).

As discussed in Remark 11, Theorem 26 can also be formulated in a stochastic setting as in Theorem 25 using the technique discussed in [119].

If \( \mathcal{C}_K=\mathbb{R}^{m\times n}\) and \( \mathcal{C}_u=\mathbb{R}^{m}\) in (184) and (192), respectively, they can both be expressed as the following quadratic optimization problem:

by defining \( {c}\in\mathbb{R}^{\mathtt{v}}\) , \( \varphi_0\in\mathbb{R}\) , \( \varphi_1\in\mathbb{R}^{\mathtt{v}}\) , and \( v\in\mathbb{R}^\mathtt{v}\) appropriately. Applying the KKT condition [153, pp. 243-244] to this problem, we can show that

This implies that the controller utilizes \( u_d\) with the feedback term \( (({\varphi_0-\varphi_1^{\top}c})/{\varphi_1^{\top}\varphi_1})\varphi_1\) only if \( \dot{x}=f(x,t)+B(x,t)u\) with \( u=u_d\) is not contracting, i.e., \( \varphi_0-\varphi_1^{\top}c > 0\) .

When \( q\) is parameterized linearly by \( \mu\) in Theorem 25, i.e., \( q = \mu x+(1-\mu)x_d\) , we have \( V_{\mathrm{NCM}}=(x-x_d)^{\top}\mathcal{M}(x-x_d)\) as in the standard Lyapunov function formulation [70, 71, 3]. Also, a linear input constraint, \( u_{\min} \leq u \leq u_{\max}\) can be implemented by using

in (184) and (192), respectively.

For the CV-STEM and NCM state estimation in Theorems 17, 23, and 24, we could design a similar Lyapunov function-based estimation scheme if its contraction constraints of Theorem 17 do not explicitly depend on the actual state \( x\) . This can be achieved by using \( A(\varrho_A,\hat{x},t)\) and \( C(\varrho_C,\hat{x},t)\) instead of \( A(\varrho_A,x,\hat{x},t)\) and \( C(\varrho_C,x,\hat{x},t)\) in Theorem 17 as in Theorem 10 [16], leading to exponential boundedness of system trajectories as derived in [13, 16, 17].

 Learning-based Autonomous Guidance with RObustness and Stability guarantees (LAG-ROS) is a DNN model for approximating a function which maps state \( x\) , local environment information \( o_{\ell}\) , and time \( t\) to an optimal robust feedback control input \( u^*\) given by (77) of Theorem 14, i.e., \( u^*=u_d-K(x-x_d)\) , or (129) of Theorem 18, i.e., \( u=u_d+\int_{0}^1kd\mu\) , where its training data is sampled as explained in Sec. 3.3.3 (see Figures 11, 12, and 14).

Let \( u_L=u_L(x,o_{\ell}(x,o_g),t)\) be the LAG-ROS in Definition 7, and let \( u^*\) be given by (77) of Theorem 14, i.e., \( u^*=u_d-K(x-x_d)\) , or (129) of Theorem 18, i.e., \( u^*=u_d+\int^1_{0}kd\mu\) . Define \( \Delta_L\) of Theorems 20 and 21 as follows:

where \( B=B(x,t)\) is the actuation matrix given in (157) and (158). Note that we use \( q(0,t)=\xi_0(t)=x_d(o_g,t)\) and \( q(1,t)=\xi_1(t)=x(t)\) in the learning-based control formulation of Theorems 20 and 21.

If the LAG-ROS \( u_L\) is learned to satisfy (136) with \( \epsilon_{\ell1}=0\) , i.e., \( \|\Delta_L\| \leq \epsilon_{\ell0}\) for all \( x\in\mathcal{S}_s\) , \( o_g\in\mathcal{S}_o\) , and \( t\in\mathcal{S}_t\) , where \( \mathcal{S}_s\subseteq\mathbb{R}^n\) , \( \mathcal{S}_o\subseteq\mathbb{R}^{g}\) , and \( \mathcal{S}_t\subseteq\mathbb{R}_{\geq0}\) are some compact sets, and if the control-affine nonlinear systems (157) and (158) are controlled by \( u=u_L\) , then (150) of Theorem 20 and (154) of Theorem 21 hold with \( \bar{d}\) and \( \bar{g}\) replaced as in (168), i.e., we have the following for \( \mathtt{e}(t)=x(t)-x_d(o_g,t)\) :

where \( x_d=x_d(o_g,t)\) and \( u_d=u_d(o_g,t)\) are as given in (159), \( \bar{d}_{\epsilon_{\ell}}=\epsilon_{\ell0}+\bar{d}_c\) , \( C_{\epsilon_{\ell}} = \bar{g}_c^2({2}{\alpha_G}^{-1}+1)+\epsilon_{\ell0}^2\alpha_d^{-1}\) , and the other variables are as defined in Theorems 20 and 21 with the disturbance terms in (168) of Theorem 22.

Let us consider the cart-pole balancing task in Fig. 8, perturbed externally by deterministic perturbation, to demonstrate the differences of (I)–(I) summarized in Table 5. Its dynamical system is given in [104, 172] with its values defined in [172].

The target trajectories \( (x_d,u_d)\) are sampled with the objective function \( \int_0^{T}\|u\|^2dt\) , and the LAG-ROS of Theorem 27 is modeled by a DNN independently of the target trajectory as detailed in [105]. The environment information \( o_g\) is selected as random initial and target terminal states, and we let \( V_{\ell}(0)=0\) and \( \epsilon_{\ell0}+\bar{d}_c = 0.75\) in (204) of Theorem 27, yielding the following exponential bound:

where \( \bar{d}_{\epsilon_{\ell}}=\epsilon_{\ell0}+\bar{d}_c\) .

Figure 8 shows the tracking errors of each motion planner averaged over \( 50\) simulations at each time instant \( t\) . The LAG-ROS (I) and robust tube-based planner (II) indeed satisfy the bound (207) (dotted black curve) for all \( t\) with a small standard deviation \( \sigma\) , unlike learning-based motion planner (I) with a diverging bound in Theorem 19 and increasing deviation \( \sigma\) . This example demonstrates one of the major advantages of contraction theory, which enables such quantitative analysis on robustness and stability of learning-based planners.

If the solution of (208) gives \( (x_d,u_d)\) , the LAG-ROS \( u_L\) of Theorem 27 ensures the perturbed solution \( x(t)\) of \( \dot{x}=f(x,t)+B(x,t)u_L+d_c(x,t)\) in (157) to stay in the admissible state space \( \mathcal{X}(o_g,t)\) , i.e., \( x(t) \in \mathcal{X}(o_g,t)\) , even with the learning error \( \epsilon_{\ell0}\) of \( u_L\) in Theorem 27.

For stochastic systems, we could sample \( x\) around \( x_d\) using a given probability distribution of the disturbance using (205) and (206) of Theorem 27.

Let us again consider the spacecraft dynamics of Example 21, but now for designing a real-time robust motion planner using the LAG-ROS in a multi-agent setting. In this example, each agent is assumed to have access only to the states of the other agents/obstacles located in a circle with a limited radius, \( 2.0\) m, centered around it. The performance of LAG-ROS is compared with (I), (II), and the following centralized motion planner (I) which is not computable in real-time:

1. Centralized robust motion planner:
   \( (x,x_d,u_d,t) \mapsto u^*\) , offline centralized solution of (II).

The objective function for sampling target trajectories of (159) is selected as \( \int_0^{T}\|u\|^2dt\) with their target terminal positions \( x_f=[p_{xf},p_{yf},0,0]^{\top}\) , where \( (p_{xf},p_{yf})\) is a random position in \( (0,0) \leq (p_{xf},p_{yf}) \leq (5,5)\) . For the state constraint satisfaction in Lemma 7, the following error tube (204):

is used with an input constraint \( u_i\geq 0, \forall i\) , to avoid collisions with a random number of multiple circular obstacles and of other agents, even under the learning error and disturbances. See [105] for the LAG-ROS training details [105].

Figure 13 shows one example of the trajectories of the motion planners (I)–(I) under deterministic external disturbances. It implies the following:

• For (I), the tracking error accumulates exponentially with time due to the lack of robustness as proven in Theorem 19.
• (II) only yields locally optimal \( (x_d,u_d)\) of (159), as its time horizon has to be small enough to make the problem solvable online with limited computational capacity, only with local environment information. This renders some agents stuck in local minima as depicted in Fig. 13.
• LAG-ROS (I) tackles these two problems by providing formal robustness and stability guarantees of Theorem 27, whilst implicitly knowing the global solution (only from the local information \( o_{\ell}\) as in [115]) without computing it online. It indeed satisfies the state constraints due to Lemma 7 as can be seen from Fig. 13.

See [105] for the in-depth discussion on this simulation result.

The matrix \( B\) in (212) does not depend on \( \theta\) , and \( \exists \Pi(x) \in \mathbb{R}^{{c}\times n}\) s.t. \( \Pi(x)^{\top}\vartheta = f(x,\theta_n)-f(x,\theta)\) , where \( \vartheta=\theta_n-\theta\) and \( \theta_n\) is a nominal guess of the uncertain parameter \( \theta\) .

Suppose that Assumption 1 holds and let \( \mathcal{M}\) defines the NCM of Definition 5, which models \( M\) of the CV-STEM contraction metric in Theorem 14 for the nominal system (212) with \( \theta = \theta_n\) , constructed with an additional convex constraint given as \( \partial_{b_i(x)}\bar{W}+\partial_{b_i(x_d)}\bar{W}=0\) , where \( \partial_{b_i(q)}\bar{W}=\sum_{i}(\partial \bar{W}/\partial q_i)b_i(q)\) for \( B=[b_1,…,b_m]\) (see [79, 85]). Suppose also that the matched uncertainty condition [85] holds, i.e., \( (\Pi(x)-\Pi(x_d))^{\top}\vartheta \in \textrm{span}(B(x))\) for \( \Pi(x)\) , and that (212) is controlled by the following adaptive control law:

where \( u_L\) is given by (163) of Theorem 24, i.e., \( u_L=u_d-R^{-1}B^{\top}\mathcal{M}(x-x_d)\) , \( \Gamma \in \mathbb{R}^{{c}\times {c}}\) is a diagonal matrix with positive elements that governs the rate of adaptation, \( \sigma \in \mathbb{R}_{\geq 0}\) , and \( (\Pi(x)-\Pi(x_d))^{\top}\vartheta=B(x)\varphi(x,x_d)^{\top}\vartheta\) .

If \( \exists \underline{\gamma},\overline{\gamma},\bar{b},\bar{\rho},\bar{\phi},\bar{\vartheta} \in \mathbb{R}_{>0}\) s.t. \( \underline{\gamma} \mathrm{I} \preceq \Gamma \preceq \overline{\gamma}\mathrm{I}\) , \( \|B(x)\| \leq \bar{b}\) , \( \|R^{-1}(x,x_d,u_d)\| \leq \bar{\rho}\) , \( \|\varphi(x,x_d)\| \leq \bar{\phi}\) , and \( \|\vartheta\| \leq \bar{\vartheta}\) , and if \( \Gamma\) and \( \sigma\) of (215) are selected to satisfy the following relation for the learning error \( \|\mathcal{M}-M\| \leq \epsilon_{\ell}\) in some compact set \( \mathcal{S}\) as in (177) of Theorem 24:

for \( \exists \alpha_{a} \in \mathbb{R}_{>0}\) , \( \alpha_{\ell}\) given in Theorem 24, and \( \underline{m}\) and \( \overline{m}\) given in \( \underline{m}\mathrm{I}\preceq M \preceq \overline{m}\mathrm{I}\) of (31), then the system (212) is robust against bounded deterministic and stochastic disturbances with \( \sigma\neq0\) , and we have the following bound in the compact set \( \mathcal{S}\) :

where \( \mathtt{e}=x-x_d\) , and \( V_{\ell}=\int^{\xi_1}_{\xi_0}\|\Theta\delta q\|\) is defined in Theorem 3 with \( M=\Theta^{\top}\Theta\) replaced by \( \mathrm{diag}(M,\Gamma^{-1})\) for \( \xi_0=[x_d^{\top},\vartheta^{\top}]^{\top}\) and \( \xi_1=[x^{\top},\hat{\vartheta}^{\top}]^{\top}\) . Furthermore, if the learning error \( \epsilon_{\ell}=0\) (CV-STEM control), (215) with \( \sigma = 0\) guarantees asymptotic stability of \( x\) to \( x_d\) in (212).

Although Theorem 28 is for the case where \( f(x)\) is affine in its parameter, it is also useful for the following types of systems with an uncertain parameter \( \theta\in\mathbb{R}^{c}\) and a control input \( u\) (see [104]):

where \( {p}\in \mathbb{R}^n\) , \( u\in \mathbb{R}^n\) , \( h:\mathbb{R}^{n}\mapsto \mathbb{R}^{n}\) , \( H:\mathbb{R}^{n}\mapsto \mathbb{R}^{n\times n}\) , \( \Pi:\mathbb{R}^{n}\mapsto \mathbb{R}^{n\times {c}}\) , \( x=[({p}^{(n-2)})^{\top},…,({p})^{\top}]^{\top}\) , and \( {p}^{(k)}\) denotes the \( k\) th time derivative of \( {p}\) . In particular, adaptive sliding control [176] designs \( u\) to render the system of the composite variable \( s\) given as \( s = {p}^{(n-1)}-{p}_r^{(n-1)}\) to be contracting, where \( {p}_r^{(n-1)} = p_d^{(n-1)}-\sum_{i=0}^{n-2}\lambda_{i}\mathtt{e}^{(i)}\) , \( \mathtt{e}=p-p_d\) , \( p_d\) is a target trajectory, and \( \kappa^{n-1}+\lambda_{n-2}\kappa^{n-2}+…+\lambda_0\) is a stable (Hurwitz) polynomial in the Laplace variable \( \kappa\) (see Example 6). Since we have \( \mathtt{e}^{(n-1)}=s-\sum_{i=0}^{n-2}\lambda_{i}\mathtt{e}^{(i)}\) and the system for \( [\mathtt{e}^{(0)},…,\mathtt{e}^{(n-2)}]\) is also contracting if \( s=0\) due to the Hurwitz property, the hierarchical combination property [10, 9] of contraction in Theorem 7 guarantees \( \lim_{t\to\infty}\|p-p_d\| = 0\)  [5, p. 352] (see Example 7).

Using the technique of Remark 19, we can construct adaptive control for the Lagrangian system in Example 6 as follows:

where \( \hat{\mathcal{H}}\) , \( \hat{\mathcal{C}}\) , and \( \hat{\mathcal{G}}\) are the estimates of \( \mathcal{H}\) , \( \mathcal{C}\) , and \( \mathcal{G}\) , respectively, and the other variables are as defined in Example 6. Suppose that the terms \( \mathcal{H}\) , \( \mathcal{C}\) , and \( \mathcal{G}\) depend linearly on the unknown parameter vector \( \theta\) as follows [5, p. 405]:

Updating the parameter estimate \( \hat{\theta}\) using the adaptation law, \( \dot{\hat{\theta}} = -\Gamma Y(\mathtt{q},\dot{\mathtt{q}},\dot{\mathtt{q}}_r,\ddot{\mathtt{q}}_r)^{\top}(\dot{\mathtt{q}}-\dot{\mathtt{q}}_r)\) , as in (215) where \( \Gamma \succ 0\) , we can define the following virtual system which has \( q = \xi_0 = [\dot{\mathtt{q}}_r^{\top},{\theta}^{\top}]^{\top}\) and \( q = \xi_1 = [\dot{\mathtt{q}}^{\top},\hat{\theta}^{\top}]^{\top}\) as its particular solutions:

where the relation \( u=-\mathcal{K}(t)(\dot{\mathtt{q}}-\dot{\mathtt{q}}_r)+Y(\mathtt{q},\dot{\mathtt{q}},\dot{\mathtt{q}}_r,\ddot{\mathtt{q}}_r)\hat{\theta}\) is used. Thus, for a Lyapunov function \( V=\delta q^{\top}\left[\begin{smallmatrix}\mathcal{H} & 0\\ 0 & \Gamma^{-1}\end{smallmatrix}\right]\delta q\) , we have that

which results in asymptotic stability of the differential state \( \delta q\) (i.e., semi-contraction [11], see also Barbalat’s lemma [5, p. 405-406]).

The dynamical system (212) is multiplicatively separable in terms of \( x\) and \( \theta\) , i.e., \( \exists\) \( Y_f:\mathbb{R}^n\mapsto\mathbb{R}^{n\times {c}_z}\) , \( Y_{b_i}:\mathbb{R}^n\mapsto\mathbb{R}^{n\times {c}_z}\) , and \( Z:\mathbb{R}^{c}\mapsto\mathbb{R}^{{c}_z}\) s.t.

where \( B(x,\theta)=[b_1(x,\theta),…,b_m(x,\theta)]\) . We could define \( \theta\) as \( [\theta^{\top},Z(\theta)^{\top}]^{\top}\) so we have \( Y_f(q)\theta = f(q,\theta)\) and \( Y_{b_i}(q)\theta = b_i(q;\theta)\) . Such an over-parameterized system could be regularized using the Bregman divergence as in [177] (see Example 29), and thus we denote \( [\theta^{\top},Z(\theta)^{\top}]^{\top}\) as \( \theta\) in the subsequent discussion.

The adaptive Neural Contraction Metric (aNCM) in Fig. 15 is a DNN model for the optimal parameter-dependent contraction metric, given by solving the adaptive CV-STEM, i.e., (106) of Theorem 14 (or Theorem 18 for differential feedback) with its contraction constraint replaced by the following convex constraint:

for deterministic systems, and

for stochastic systems, where \( W=M(x,x_d,u_d,\hat{\theta})^{-1}\succ0\) (or \( W=M(x,\hat{\theta})^{-1}\succ0\) , see Theorem 10), \( \bar{W} = \nu W\) , \( \nu = \overline{m}\) , \( R=R(x,x_d,u_d)\succ0\) is a weight matrix on \( u\) , \( \Xi=({d}/{dt})|_{\hat{\theta}}{\bar{W}}\) is the time derivative of \( \bar{W}\) computed along (212) with \( \theta = \hat{\theta}\) , \( A\) and \( B\) are given in (229), \( \alpha\) , \( \underline{m}\) , \( \overline{m}\) , and \( \alpha_s\) are as given in Theorem 14, and the arguments \( (x,x_d,u_d,\hat{\theta})\) are omitted for notational simplicity.

Suppose that Assumption 2 holds and let \( \mathcal{M}\) define the aNCM of Definition 8. Suppose also that the true dynamics (212) is controlled by the following adaptive control law:

where \( d\mathcal{M}_q = [({\partial \mathcal{M}}/{\partial q_1})\mathtt{e},…,({\partial \mathcal{M}}/{\partial q_n})\mathtt{e}]^{\top}/2\) , \( \mathtt{e}=x-x_d\) , \( \Gamma \succ 0\) , \( \sigma \in \mathbb{R}_{\geq 0}\) , and \( Y\) , \( Y_d\) , \( \tilde{Y}\) are given in (230). Suppose further that the learning error in \( \|\mathcal{M}-M\|\leq\epsilon_{\ell}\) of Theorem 24 additionally satisfies \( \|d\mathcal{M}_{x_d} -dM_{x_d} \| \leq \epsilon_{\ell}\) and \( \|d\mathcal{M}_x -dM_x \| \leq \epsilon_{\ell}\) in some compact set \( \mathcal{S}\) .

If \( \exists \bar{b},\bar{\rho},\bar{y} \in \mathbb{R}_{>0}\) s.t. \( \|B(x,\hat{\theta})\| \leq \bar{b}\) , \( \|R^{-1}(x,x_d,u_d)\| \leq \bar{\rho}\) , \( \|Y\|\leq \bar{y}\) , \( \|Y_d\|\leq \bar{y}\) , and \( \|\tilde{Y}\|\leq \bar{y}\) in (233) and (234), and if \( \Gamma\) and \( \sigma\) of (234) are selected to satisfy the following as in Theorem 28:

for \( \exists \alpha_a\in\mathbb{R}_{>0}\) , \( \alpha_{\ell}\) given in Theorem 24, and \( \underline{m}\) and \( \overline{m}\) given in \( \underline{m}\mathrm{I}\preceq M \preceq \overline{m}\mathrm{I}\) of (31), then the system (212) is robust against bounded deterministic and stochastic disturbances, and we have the exponential bound (217) in the compact set \( \mathcal{S}\) . Furthermore, if \( \epsilon_{\ell}=0\) (adaptive CV-STEM control), (234) with \( \sigma = 0\) guarantees asymptotic stability of \( x\) to \( x_d\) in (212).

If \( \hat{\theta}\) is given by (236) with \( \hat{\theta}(0) \in \Omega_{\theta} = \{\theta\in\mathbb{R}^c|p(\theta) \leq 1\}\) for a convex function \( p(\theta)\) , then \( \hat{\theta}(t) \in \Omega_{\theta}, \forall t\geq 0\) (e.g.{}, \( \|\hat{\theta}\| \leq \theta_{\max}\sqrt{1+\epsilon_{\theta}}\) if \( p(\hat{\theta})=(\hat{\theta}^{\top}\hat{\theta}-\theta_{\max}^2)/(\epsilon_{\theta}\theta_{\max}^2)\) ).

The dependence on \( u\) and \( \dot{\hat{\theta}}\) in \( (d/dt)|_{\hat{\theta}}M\) can be removed by using \( \partial_{b_i(x)}M+\partial_{b_i(x_d)}M=0\) as in Theorem 28, and by using adaptation rate scaling introduced in [86]. In essence, the latter multiplies the adaptation law (234) by any strictly-increasing and strictly-positive scalar function \( v(2\rho)\) , and update \( \rho\) as

where \( v_{\rho}=\partial v/\partial \rho\) , \( \eta \in\mathbb{R}_{>0}\) , and \( V_{\mathtt{e}} = \mathtt{e}^{\top}M(x,x_d,u_d,\hat{\theta})\mathtt{e}\) for \( \mathtt{e}=x-x_d\) and \( M\) given in Definition 8, so the additional term due to (239) cancels out the term involving \( \dot{\hat{\theta}}\) in \( (d/dt)|_{\hat{\theta}}\bar{W}\) of (231) (see [86] for details). Its robustness property follows from Theorem 29 also in this case.

Using the Bregman divergence-based adaptation in [177], we could implicitly regularize the parameter estimate \( \hat{\theta}\) as follows:

where \( d_{\psi}\) is the Bregman divergence defined as \( d_{\psi}(x\|y)=\psi(x)-\psi(y)-(x-y)^{\top}\nabla\psi(y)\) for a convex function \( \psi\) , and \( \mathcal{A}\) is a set containing only parameters that interpolate the dynamics along the entire trajectory. If \( \hat{\theta}(0)=\mathrm{arg}\min_{b\in \mathbb{R}^p}\psi(b)\) , we have \( \lim_{t\to\infty}\hat{\theta}=\mathrm{arg}\min_{\vartheta\in \mathcal{A}}\psi(\vartheta)\) , which regularizes \( \hat{\theta}\) to converge to a parameter that minimizes \( \psi\) . For example, using \( 1\) -norm for \( \psi\) would impose sparsity on the steady-state parameter estimate \( \hat{\theta}\)  [177].

Note that the results presented earlier in this section do not necessarily mean parameter convergence, \( \lim_{t\to\infty}\|\tilde{\theta}\|=0\) , as the adaptation objective is to drive the tracking error \( \|x-x_d\|\) to zero [5, p. 331], not to find out the true parameter \( \theta\) out of the many that achieve perfect tracking.

Asymptotic parameter convergence, \( \lim_{t\to\infty}\|\tilde{\theta}\|=0\) for \( \theta\) of (212) under Assumption 2, could be guaranteed if there is no disturbance and learning error with \( \sigma=0\) , and we have \( \exists T,\alpha_{PE}\in\mathbb{R}_{>0}\text{ s.t. }\int_{t}^{t+T}\tilde{Y}^{\top}\tilde{Y}d\tau \succeq \alpha_{PE}\mathrm{I}, \forall{t}\) for \( \tilde{Y}\) given in (230) (the persistent excitation condition [5, p. 366]). We could also utilize the Bregman divergence-based adaptation to regularize the behavior of \( \lim_{t\to\infty}\|\hat{\theta}\|\) as in Example 29.

Define a mapping \( \mathcal{F}\) as \( \mathcal{F}(u) = u^*(x,t)-r_L(x,u,t)\) , where \( u^*\) and \( r_L\) are given in (248). If \( r_L\) is Lipschitz in \( u\) with a 2-norm Lipschitz constant \( L_u < 1\) , i.e., \( \|r_L(x,u,t)-r_L(x,u',t)\| \leq L_u\|u-u'\|, \forall u,u'\) , then \( \mathcal{F}\) is a contraction mapping for fixed \( x,t\) . Therefore, if \( x,t\) are fixed, discrete-time nonlinear control \( u_k\) defined as

converges to a unique solution \( u\) given by \( u=\mathcal{F}(u)\) .

Let \( x\) be the trajectory of the following ideal system without the unknown part \( r\) of the dynamics (247):

where \( u^*\) is a stabilizing controller that renders (251) contracting as in Theorem 30 for \( M\) which satisfies \( \underline{m}\mathrm{I}\preceq M \preceq \overline{m}\mathrm{I}\) of (31). Note that such \( u^*\) can be designed by using, e.g.{}, Theorems 14 and 18. Suppose that the true dynamics (247) is controlled by (249) and

for \( x^*\) in (247) [61]. If \( \exists\bar{b}\in\mathbb{R}_{\geq 0}\) s.t. \( \|B(x,t)\|\leq\bar{b}\) , and if \( r_L\) is modeled by a spectrally-normalized DNN of Definition 6 to have

for all \( x\in\mathcal{S}_s\) , \( u\in\mathcal{S}_u\) , and \( t\in\mathcal{S}_t\) , where \( \mathcal{S}_s\subseteq\mathbb{R}^n\) , \( \mathcal{S}_u\subseteq\mathbb{R}^m\) , and \( \mathcal{S}_t\subseteq\mathbb{R}_{\geq 0}\) are some compact sets, then \( r_L\) is Lipschitz continuous, and the controller (249) applied to (247) gives the following bound in the compact set:

as long as the Lipschitz constant of \( r_L\) is selected to have

where \( V_{\ell}=\int^{x}_{x^*}\|\Theta\delta x\|\) as in Theorem 3 with \( M=\Theta^{\top}\Theta\) . Furthermore, the system (247) with the controller (249) is robust against deterministic and stochastic disturbances.

Let us consider the following Lagrangian system of Example 6 perturbed externally by unknown control non-affine residual disturbance:

where \( x=[\mathtt{q}^{\top},\dot{\mathtt{q}}^{\top}]^{\top}\) and the other variables are as given in Example 6. Using the result of Theorem 31, we can design a discrete-time nonlinear controller by augmenting the exponentially stabilizing controller of Example 6 with a learned residual part \( r_L(x,u)\) as follows:

where \( \dot{\mathtt{q}}_r=\dot{\mathtt{q}}_d(t)-\Lambda(t)(\mathtt{q}-\mathtt{q}_d(t))\) , \( \mathcal{K}: \mathbb{R}_{\geq 0}\mapsto\mathbb{R}^{n\times n}\) , \( \Lambda: \mathbb{R}_{\geq 0}\mapsto\mathbb{R}^{n\times n}\) , and \( (\mathtt{q}_d,\dot{\mathtt{q}}_d)\) is the target trajectory of the state \( (\mathtt{q},\dot{\mathtt{q}})\) , and \( \mathcal{K},\Lambda \succ 0\) are control gain matrices (design parameters). Again, note that \( \dot{\mathcal{H}}-2\mathcal{C}\) is skew-symmetric with \( \mathcal{H} \succ 0\) by construction. This gives us the following closed-loop virtual system of a smooth path \( q(\mu,t)\) parameterized by \( \mu \in [0,1]\) , which has \( q(\mu=0,t)=\dot{\mathtt{q}}_r\) and \( q(\mu=1,t)=\dot{\mathtt{q}}\) as its particular solutions as in Example 6, but now with non-zero perturbation due to \( r(x,u)\) :

After some algebra as in the proof of Theorem 31, we can show that

where \( \mathcal{H}=\Theta^{\top}\Theta\) , \( h_{\ell}\mathrm{I}\preceq\mathcal{H}\preceq h_u\mathrm{I}\) , \( k_{\ell}\mathrm{I}\preceq\mathcal{K}\) , \( \|u_k-u_{k-1}\| \leq \rho\|\dot{\mathtt{q}}-\dot{\mathtt{q}}_r\|\) , \( \epsilon_{\ell}\) is the learning error of \( r_L\) , and \( L_u\) is the Lipschitz constant of \( r_L\) (assuming \( r_L\) is modeled by a spectrally-normalized DNN of Definition 6). This indeed indicates that the tracking error of the Lagrangian system (258) is exponentially bounded as proven in Theorem 31.

Suppose that \( \mathcal{S}_s \subseteq \mathbb{R}^n\) is full-dimensional, \( n \geq 2\) , and \( \dot{x} = f(x)\) is contracting in a metric defined by a uniformly bounded \( M(x)\succ 0\) , i.e.,

of (6) and \( \underline{m}\mathrm{I}\preceq M \preceq \overline{m}\mathrm{I}\) of (31) hold for \( M\) . Suppose also that a matrix-valued function \( \mathcal{M}(x)\) with \( \mathcal{M}(x)\succeq \underline{m}_L \mathrm{I}\) is learned to satisfy \( \nu(Z_b)\leq\epsilon_{\ell}\) , where the contraction violation set \( Z_b\) for \( p(x,\delta x) = (f(x),(\partial f/\partial x)\delta x)^{\top}\) is defined as

with \( \mathcal{V}: \mathcal{T}S \mapsto \mathbb{R}\) being a learned differential Lyapunov function \( \mathcal{V}(x,\delta x) = \delta x^{\top}\mathcal{M}(x)\delta x\) for \( S = \cup_{t\in \mathcal{S}_t}\varphi_t(\mathcal{S}_s)\) . Define \( \bar{B}\) , \( B_H\) , \( B_{\nabla q}\) , and \( B_{\nabla \mathcal{V}}\) as

where \( q = \nabla \mathcal{V}(x)^{\top} f(x)\) . Also, we define \( r_{\epsilon_{\ell}}\) and \( r_b\) as

where \( \eta \in (0,1)\) . Finally, define \( \tilde{X}_t(r_b)\) as

for \( t \in \mathcal{S}_t\) , where \( \tilde{X} = \{\xi \in \mathcal{S}_s: \mathbb{B}^n_2(\xi,r_{\epsilon_{\ell}}) \subset \mathcal{S}_s\}\) . Then the system is contracting in the learned metric defined by \( \mathcal{M}(x)\) at the rate \( (1-\eta)\alpha\) with \( \eta \in (0,1)\) , for every \( x \in \tilde{S}(r_b) = \cup_{t\in \mathcal{S}_t}\varphi_t(\tilde{X}_t(r_b))\) .