Part IV · 第 13 章反馈控制 · 13.2

状态空间与稳定性

PID 一次只看一个误差信号；状态空间把机器人的全部"记忆"写成一个向量 $\mathbf x$，把动力学写成 $\dot{\mathbf x}=A\mathbf x+B\mathbf u$。这一节在一个摆上走完整条链路：在平衡点处线性化，用特征值判断稳定性，用能量当作李雅普诺夫证书，用可控性判断"推得动推不动"，用极点配置把倒立摆立起来，再看把控制器放进计算机（离散化）和只测到一部分状态（观测器）时会发生什么。

约 45 分钟前置：1.11 特征值 · 2.5 常微分方程 · 3.3 雅可比 · 3.6 向量场与流 · 13.1 PID 6 个交互

1为什么机器人需要它

13.1 节的 PID 把一个关节当成"误差进、力矩出"的黑盒，调三个旋钮。它在单个关节上很好用，可是一旦遇到下面这些问题，就需要一个更完整的描述：

多个量互相牵连：倒立摆、平衡车、四旋翼的姿态和位置互相耦合，推一下小车，摆杆也跟着动。单独给每个量调 PID，常常调不出来。
平衡点不止一个：同一个摆，挂着是稳定的，倒立是不稳定的。同一个模型在不同工作点表现完全不同，控制器要知道自己在哪个点附近工作。
想要证明而不只是"试一试"：控制器会不会在某些初始状态下发散？能离平衡点多远而不失控？这需要一个不用把方程解出来也能判断稳定性的工具。
控制器跑在计算机上，传感器只测到一部分：编码器测角度不测角速度，IMU 测角速度不测位置；控制周期是 1 ms 还是 50 ms 也会改变稳定性。

状态空间给这些问题一个统一的语言：把系统此刻的全部"记忆"写成状态向量 $\mathbf x$，把动力学写成一阶方程组 $\dot{\mathbf x}=\mathbf f(\mathbf x,\mathbf u)$（2.5 节的标准形式）。在平衡点附近把它线性化，就得到全书控制部分的出发点：

$$ \dot{\mathbf x}=A\,\mathbf x+B\,\mathbf u,\qquad \mathbf y=C\,\mathbf x,\qquad \mathbf x\in\R^n,\ \mathbf u\in\R^m,\ \mathbf y\in\R^p . $$

四个问题依次是：它自己稳不稳（$A$ 的特征值，李雅普诺夫函数），输入能不能推动每个方向（可控性，$A$ 和 $B$），怎样用反馈 $\mathbf u=-K\mathbf x$ 把它变成我们想要的样子（极点配置），以及测量 $\mathbf y$ 够不够还原整个 $\mathbf x$（可观性与观测器，$A$ 和 $C$）。13.3 节的 LQR 和 13.4 节的 MPC 都建立在这四个问题的答案上。全篇用同一个摆：质量 $m=1$ kg 的质点挂在长 $l=1$ m 的轻杆上（转动惯量 $ml^2=1$ kg·m²），关节粘性阻尼 $b$，电机力矩 $u$，角度 $\theta$ 从竖直向下量起。

2直觉：状态是一个点，动力学是一片箭头

要预测摆下一刻往哪走，只知道角度不够：同一个角度，可能正往上甩，也可能正往下落。角度 $\theta$ 加上角速度 $\omega=\dot\theta$ 才够，这两个数就是状态 $\mathbf x=(\theta,\omega)$。牛顿定律 $ml^2\ddot\theta=-mgl\sin\theta-b\dot\theta+u$ 写成状态形式：

$$ \dot{\mathbf x}=\mathbf f(\mathbf x,u)=\begin{bmatrix}\omega\\[2pt] -\dfrac{g}{l}\sin\theta-\dfrac{b}{ml^2}\,\omega+\dfrac{u}{ml^2}\end{bmatrix}. $$

$\mathbf f$ 在 $(\theta,\omega)$ 平面（相平面）的每一点放一支箭头，这正是 3.6 节的向量场：状态顺着箭头流动。箭头为零的点叫平衡点：$\mathbf f(\mathbf x^*,u^*)=\mathbf 0$，即 $\omega^*=0$ 且 $mgl\sin\theta^*=u^*$。不加力矩（$u^*=0$）时有两个：挂着的 $\theta^*=0$ 和倒立的 $\theta^*=\pi$。

平衡点附近，向量场长什么样？3.3 节说过，放大到足够近，光滑的映射都像它的雅可比。令 $\delta\mathbf x=\mathbf x-\mathbf x^*$，$\delta u=u-u^*$，一阶泰勒展开给出

$$ \delta\dot{\mathbf x}\approx \underbrace{\frac{\partial\mathbf f}{\partial\mathbf x}\Big|_{*}}_{\hl{A}{A}}\delta\mathbf x+\underbrace{\frac{\partial\mathbf f}{\partial u}\Big|_{*}}_{B}\delta u,\qquad A=\begin{bmatrix}0&1\\ -\dfrac gl\cos\theta^* & -\dfrac{b}{ml^2}\end{bmatrix},\quad B=\begin{bmatrix}0\\ \dfrac1{ml^2}\end{bmatrix}. $$

只有 $\cos\theta^*$ 一项随平衡点变化，但它决定了一切。挂着时 $\cos\theta^*=1$，重力像一根把摆拉回去的弹簧；倒立时 $\cos\theta^*=-1$，同一个重力变成了把摆推开的"负弹簧"。代入本课的数值（$g/l=9.81\ \text{s}^{-2}$，$b=0.5$）：

$$ \begin{gathered} A_{\downarrow}=\begin{bmatrix}0&1\\-9.81&-0.5\end{bmatrix},\ \ \lambda=-0.25\pm3.12\,j;\\[6pt] A_{\uparrow}=\begin{bmatrix}0&1\\9.81&-0.5\end{bmatrix},\ \ \lambda=2.89,\ -3.39 . \end{gathered} $$

1.11 节和 2.5 节已经告诉我们怎么读特征值：线性系统的解是 $\sum_ic_ie^{\lambda_it}\mathbf v_i$，实部全为负就衰减。挂着的点是一个稳定螺旋（实部 $-0.25$，每秒振 $3.12/2\pi\approx0.5$ 圈）；倒立点有一个正特征值 $2.89$，沿它的特征向量，偏差每秒放大 $e^{2.89}\approx18$ 倍，这就是鞍点。3.6 节用迹–行列式图把二维线性系统分了类；这里关心的是另一件事：线性化在多大范围内可信，以及平衡点怎样随输入移动。

3动手玩：摆的相平面与它的线性化

在左边的相平面上点一下，就从那个状态放开摆（横轴是角度，纵轴是角速度）。琥珀色实线是真实的非线性运动，紫色虚线是在最近的平衡点处线性化后的预测 $\delta\mathbf x(t)=e^{At}\delta\mathbf x_0$。下方的曲线把两者对时间画出来，右上角的小图是摆本身。拖动"恒定力矩 $u_0$"：平衡点会移动，$A$ 跟着变。

FIG 1摆的相平面：平衡点、线性化与它的适用范围

相平面 (θ, ω)：点击放开摆

偏离平衡点 θ − θ* (rad) 对时间 (s)

小箭头是向量场 $\mathbf f$。绿色实心点 = 稳定平衡点，红色圆环 = 不稳定平衡点；在实特征值的平衡点上画出特征向量方向：青色 = 衰减方向（$\lambda<0$），红色 = 放大方向（$\lambda>0$）。蓝色虚线是分界线（倒立点的稳定流形）：它们之间的状态最终回到同一个挂着的平衡点。角度按 $2\pi$ 折回到 $[-\pi,\pi)$，所以 $\theta=\pm\pi$ 是同一个倒立点。仿真用 RK4，步长 2 ms，每条轨迹 8 s；"线性化误差"是 8 s 内 $\max|\delta\theta-\delta\theta_L|$（$\delta\theta_L$ 是线性预测）除以 $\max|\delta\theta_L|$。

试一试

在挂着的平衡点附近放开两次：一次幅值小于 20°（0.35 rad），线性化误差低于 2%；一次幅值超过 60°（1.05 rad），误差超过 15%。线性化在多大范围内"够用"？
加大阻尼 $b$，直到稳定平衡点的两个特征值变成实数（螺旋变成结点，出现两条青色方向）。在 $u_0=0$ 时这个阈值是 $b=2\sqrt{g/l}\approx6.26$；加上力矩 $u_0$ 后阈值变低还是变高？

再把 $u_0$ 推过 $mgl=9.81$ N·m：两个平衡点相互靠近、合并、然后消失（"鞍结分岔"），摆从此一圈圈地转下去。平衡点不是系统固有的，它是"状态 + 输入"的一对 $(\mathbf x^*,u^*)$，线性化也要在这一对上做。

核心洞见

线性化是一张局部地图：在平衡点附近它几乎完美（20° 以内误差不到 2%），离远了就会出错，而且它完全看不见全局结构：分界线、第二个平衡点、翻过顶点的旋转。"$A$ 的特征值实部都为负"只说明平衡点局部稳定（李雅普诺夫间接法），要回答"多远都能回来"，需要下一节的工具。

4第二视角：不解方程也能判断稳定（李雅普诺夫）

有阻尼的摆最后一定会停在下面，为什么这么肯定？因为它的能量只减不增：

$$ \begin{gathered} E(\mathbf x)=\tfrac12ml^2\omega^2+mgl\,(1-\cos\theta),\\[6pt] \dot E=ml^2\omega\dot\omega+mgl\sin\theta\,\omega=\omega\,(u-b\,\omega)\;\overset{u=0}{=}\;-b\,\omega^2\le0 . \end{gathered} $$

李雅普诺夫（1892）把这个论证抽象出来：找一个函数 $V(\mathbf x)$，它在平衡点为零、别处为正（像一只碗），而且沿着系统的运动不增加，$\dot V=\nabla V\cdot\mathbf f\le0$。那么状态就只能留在碗里，平衡点是稳定的；如果 $\dot V<0$ 处处成立（平衡点除外），状态会一路滑到碗底，平衡点是渐近稳定的。关键在于 $\dot V=\nabla V\cdot\mathbf f$ 只用到向量场本身，不需要解微分方程。

但 $V$ 要自己去找，而且找错了不说明任何问题。下图给出三个候选：能量 $E$；一个"看起来很自然"的 $\tfrac12(\theta^2+\omega^2)$；以及解李雅普诺夫方程 $A^\top P+PA=-I$ 得到的二次型 $\mathbf x^\top P\mathbf x$（第 6 节推导）。红色区域是 $\dot V>0$ 的地方，在那里状态会往碗外爬；绿色区域是这个 $V$ 能担保的吸引域：最大的、内部没有红色的碗。

FIG 2李雅普诺夫函数：哪个碗能证明稳定？

相平面：点击放开摆

沿轨迹的 V(t)

左：青色曲线是 $V$ 的等值线（碗的"等高线"），红色阴影是 $\dot V>0$ 的区域，绿色阴影是被证明的吸引域 $\{V<c^*\}$，绿色曲线是它的边界，$c^*$ 是红色区域里 $V$ 的最小值（能量的情形另行处理，见第 6 节的 LaSalle 原理）。能量候选下的琥珀色虚线是 $\dot V=0$ 的集合 $\omega=0$。蓝色曲线是你放开的轨迹（10 s，RK4，步长 2 ms），小圆点是平衡点。右：沿最新一条轨迹的 $V(t)$，绿色虚线是 $c^*$。窗口 $\theta\in[-4,4]$，$\omega\in[-7,7]$；$c^*$ 在这个窗口内的 $301\times301$ 网格上搜索。

试一试

选"朴素二次型" $\tfrac12(\theta^2+\omega^2)$，放开一条轨迹，让 $V(t)$ 在途中比起点高出 20% 以上。摆明明在稳定地停下来，这个 $V$ 却在上升：它什么也证明不了。
选"李雅普诺夫方程"（$b>0$），从绿色区域外面（$V(\mathbf x_0)>1.05\,c^*$）放开一条不翻过顶点、最终停在下面的轨迹。证书是充分条件：它没担保的地方，系统照样可以稳定。

核心洞见

能量是这个摆最好的证书：它担保的是整个 $E<2mgl$ 的"眼睛"（无阻尼摆的分界线所围的区域，有阻尼时真实的吸引域还要更大），因为它来自物理。由线性化解出的二次型只在 $|\theta|\lesssim44^\circ$ 的小椭圆里有效（$b=0.5$ 时），它对非线性一无所知。机器人里的稳定性证明大多是这种思路：8.3 节证明"重力补偿 + PD"稳定，用的就是"动能 + 弹簧势能"这个 $V$。

5第三视角：输入推得动哪些方向？然后把极点放到想要的位置

可控性

只有一个输入，却要控制两个（或更多）状态，这可能吗？对小车来说可以：力改变速度，速度再改变位置，输入通过动力学传到了位置上。$B$ 是输入直接推的方向，$AB$ 是"推了之后，动力学把它带向"的方向，$A^2B$ 是再下一步……只要 $B, AB, \dots, A^{n-1}B$ 张成整个空间，输入就能到达每个方向。这就是卡尔曼（1960）的秩判据：$(A,B)$ 可控，当且仅当下面这个可控性矩阵的秩等于 $n$。

$$ \hl{Ctrb}{\mathcal C=\begin{bmatrix}B & AB & \cdots & A^{n-1}B\end{bmatrix}},\qquad \operatorname{rank}\mathcal C=n . $$

秩判据只回答"能不能"，不回答"费多大劲"。可控性格拉姆矩阵 $W_T$（第 6 节）回答后者：从原点出发、在 $T$ 秒内、用输入能量 $\int_0^Tu^2dt\le1$ 能到达的状态，恰好组成椭圆 $\{\mathbf x:\mathbf x^\top W_T^{-1}\mathbf x\le1\}$。椭圆越扁，那个方向越"难推"；扁成一条线，就不可控了。下图拖动 $B$（琥珀色箭头）和目标点（红色十字）。

FIG 3可控性：一个输入能推到哪里，要花多少能量

状态平面：拖动 B 和目标

最小能量输入 u(t)

琥珀色箭头 = $B$，紫色箭头 = $AB$，青色虚线 = $A$ 的实特征向量方向（$B$ 靠近它们 5° 以内时会吸附上去）。蓝色椭圆 = $T$ 秒内用单位输入能量能到达的状态集合。绿色曲线 = 从原点到目标的最小能量轨迹，右图是对应的输入 $u(t)$。"所需能量"是 $\int_0^Tu^2\,dt=\mathbf x_T^\top W_T^{-1}\mathbf x_T$；$W_T$ 用 RK4 对 $\dot W=AW+WA^\top+BB^\top$ 积分 400 步得到。

试一试

选"两个一阶模态"（$A=\operatorname{diag}(-1,-3)$），把 $B$ 拖到一条青色特征方向上：椭圆塌成一条线，$\det\mathcal C=0$。这时哪个模态完全不受输入影响？
选"小车"，目标保持在"静止于 1 m 外"附近（位置 $\ge0.8$，速度绝对值 $\le0.2$），把 $T$ 从 2 s 减到 1 s 或更短：所需能量从 1.5 涨到 12 以上。为什么时间减半、能量涨 8 倍？

再试试"两个相同模态"（$A=-2I$）：无论 $B$ 指向哪里，椭圆都是一条线。一个输入没法分别控制两个一模一样、互不相干的系统，因为它对两者的作用永远成同一比例。

极点配置：把倒立摆立起来

倒立点不稳定（$\lambda=+3.03$，这里取 $b=0.2$），但系统可控。用状态反馈 $u=-K\mathbf x=-k_1\delta\theta-k_2\omega$，闭环变成 $\dot{\mathbf x}=(A-BK)\mathbf x$。对本例，闭环特征多项式是

$$ \begin{gathered} \det\big(s I-(A-B\hl{K}{K})\big)=s^2+(b+k_2)\,s+\big(k_1-\tfrac gl\big)\;\overset{!}{=}\;(s-p_1)(s-p_2)\\[6pt] \Longrightarrow\quad k_1=\tfrac gl+p_1p_2,\qquad k_2=-(p_1+p_2)-b . \end{gathered} $$

两个增益正好对应两个系数，所以任意想要的极点（复数要成对）都能配出来：这就是极点配置定理，可控是它的充要条件。在下图左侧的 $s$ 平面上拖动闭环极点，右侧的摆会用对应的 $K$ 从初始偏角站起来。仿真是非线性的摆，电机力矩有上限，紫色虚线是线性模型的预测。

FIG 4极点配置：倒立摆，带力矩上限

s 平面：拖动闭环极点

倒立摆（δ = θ − π）

偏角 δ (°) 对时间 (s)

电机力矩 u (N·m)

$s$ 平面：紫色叉 = 开环极点 $3.03$ 与 $-3.23$，琥珀色实心点 = 你拖动的闭环极点（空心圆 = 它的共轭）；把极点拖到实轴上（离实轴 0.4 以内会吸附），会出现第二个可以沿实轴拖动的极点。红色阴影是右半平面。摆和曲线：琥珀色 = 非线性仿真（RK4，1 ms），紫色虚线 / 空心圆 = 线性预测（无饱和）；力矩箭头和力矩曲线是绿色，饱和时变红，红色虚线是 $\pm u_{\max}$。"调节时间"是 $|\delta|$ 最后一次超过初始偏角 2% 的时刻；偏角超过 90° 判为摔倒。

试一试

初始偏角 $\ge20^\circ$、力矩上限 $\le20$ N·m，配出一组极点：调节时间 < 0.8 s，超调 < 5%，而且全程不饱和。提示：一对重极点 $-6$ 需要 0.97 s；试试略带虚部的一对复极点。
保持两个闭环极点都在左半平面，把初始偏角调到 40°，再降低力矩上限，直到摆摔倒。线性分析说"稳定"，为什么它还是倒了？在 40° 处，最少要多大力矩才可能托住它？

核心洞见

可控时，极点"想放哪就放哪"，但代价写在增益里：重极点 $-p$ 需要 $k_1=g/l+p^2$，从 20° 起步的初始力矩随 $p^2$ 增长（$p=6$ 时 16.0 N·m，$p=7$ 时 20.5 N·m）。极点放得越快，力矩越大、越先撞上饱和，对噪声和采样延迟也越敏感。饱和以后系统不再是线性的，极点承诺的东西就不再成立。13.3 节的 LQR 用代价函数把"快"和"省力"放在同一把秤上，自动选出极点。

6再看数学

提示：把鼠标悬停在公式中带底色的符号上，对应的图元素会被高亮：$A$ 加粗 FIG 1 的特征方向，$V$ 与 $\dot V$ 加粗 FIG 2 的等值线与红色区域，$\mathcal C$ 与 $W_T$ 加粗 FIG 3 的 $B$、$AB$ 箭头与椭圆，$K$ 放大 FIG 4 的闭环极点，$A_d$ 与 $I+AT$ 放大 FIG 5 的琥珀色与紫色极点，$L$ 加粗 FIG 6 的估计曲线。

状态空间模型与它的解

线性时不变（LTI）系统 $\dot{\mathbf x}=A\mathbf x+B\mathbf u$，$\mathbf y=C\mathbf x+D\mathbf u$（机器人里 $D$ 几乎总是 0）。对 $e^{-At}\mathbf x$ 求导，得 $\tfrac{d}{dt}(e^{-At}\mathbf x)=e^{-At}B\mathbf u$，两边积分即得"常数变易"公式：第一项是自由响应，第二项是输入引起的强迫响应。

$$ \mathbf x(t)=e^{At}\mathbf x_0+\int_0^t e^{A(t-\tau)}B\,\mathbf u(\tau)\,d\tau,\qquad e^{At}=\sum_{k\ge0}\frac{(At)^k}{k!} . $$

第一项是系统"自己"怎么走（3.6 节的线性流），第二项是输入的贡献：每一时刻的输入 $B\mathbf u(\tau)$ 被系统的自由运动 $e^{A(t-\tau)}$ 一路带到时刻 $t$。状态的选法不唯一：换坐标 $\mathbf x=T\mathbf z$ 得到 $(T^{-1}AT,\,T^{-1}B,\,CT)$，特征值不变，可控性、可观性也不变，所以这些都是系统本身的性质，而不是坐标的。

线性化与特征值稳定性

对非线性系统 $\dot{\mathbf x}=\mathbf f(\mathbf x,\mathbf u)$，在平衡对 $(\mathbf x^*,\mathbf u^*)$（$\mathbf f(\mathbf x^*,\mathbf u^*)=0$）处取 $A=\partial\mathbf f/\partial\mathbf x$、$B=\partial\mathbf f/\partial\mathbf u$，误差项是 $O(\|\delta\mathbf x\|^2+\|\delta\mathbf u\|^2)$（2.4 节）。LTI 系统 $\dot{\mathbf x}=A\mathbf x$ 的稳定性完全由特征值决定：

所有 $\operatorname{Re}\lambda_i<0$（$A$ 是 Hurwitz 的）$\iff$ 渐近稳定，而且是指数稳定：$\|\mathbf x(t)\|\le c\,e^{-\alpha t}\|\mathbf x_0\|$。
某个 $\operatorname{Re}\lambda_i>0$ $\Rightarrow$ 不稳定。
实部为零的特征值是临界情况：它们若是单重的（或每个都有足够多的特征向量），解有界但不衰减（无阻尼的摆，"中心"）；若是重根而特征向量不够（约当块），解会按 $t^k$ 增长。双积分器 $\ddot q=0$ 的两个特征值都是 0，可 $q(t)=q_0+\dot q_0t$ 一直跑远：不稳定。

李雅普诺夫间接法把这些结论带回非线性系统：若线性化 $A$ 是 Hurwitz 的，平衡点局部指数稳定；若 $A$ 有正实部特征值，平衡点不稳定；若有实部为零的特征值，线性化无法判断，要看高阶项（3.6 节提到的 Hartman–Grobman 定理只覆盖双曲平衡点）。无阻尼摆（$b=0$）的下平衡点正是这种情况：线性化是中心，判断要靠下面的能量论证。

李雅普诺夫直接法与 LaSalle 不变原理

平衡点 $\mathbf x^*=0$ 稳定：对任意 $\varepsilon>0$ 存在 $\delta>0$，使 $\|\mathbf x_0\|<\delta$ 的解永远满足 $\|\mathbf x(t)\|<\varepsilon$；渐近稳定：稳定且附近的解趋于 0。定理（李雅普诺夫，1892）：设在包含原点的区域 $D$ 上有连续可微的 $V$，满足

$$ \hl{V}{V(\mathbf 0)=0,\quad V(\mathbf x)>0\ \ (\mathbf x\ne\mathbf 0)},\qquad \hl{Vdot}{\dot V(\mathbf x)=\nabla V(\mathbf x)^{\top}\mathbf f(\mathbf x)\le0}\qquad(\mathbf x\in D) . $$

则原点稳定；若 $\dot V<0$（$\mathbf x\ne0$），则渐近稳定；若再有 $D=\R^n$ 且 $\|\mathbf x\|\to\infty$ 时 $V\to\infty$（径向无界），则全局渐近稳定。直觉：$V$ 的下水平集 $\Omega_c=\{V\le c\}$ 是一只只嵌套的碗，$\dot V\le0$ 说明状态不会爬出任何一只碗。

摆的能量只给出 $\dot E=-b\omega^2\le0$，在整条 $\omega=0$ 线上都等于零，不满足严格不等式。LaSalle 不变原理（1960）补上了这一步：若 $\Omega$ 是有界的正不变集且其中 $\dot V\le0$，则从 $\Omega$ 出发的解趋于 $\{\dot V=0\}$ 中最大的不变集。在 $\omega=0$ 线上能"待住"的解必须一直 $\omega\equiv0$，于是 $\dot\omega=-(g/l)\sin\theta=0$，只剩平衡点。取 $\Omega=\{E<2mgl,\ |\theta|<\pi\}$（能量低于翻过顶点所需的值），其中只有 $\theta=0$ 这一个平衡点，所以里面的每条轨迹都回到挂着的位置：这正是 FIG 2 中能量候选的绿色"眼睛"。

线性系统的二次型证书。取 $V=\mathbf x^\top P\mathbf x$，则 $\dot V=\mathbf x^\top(A^\top P+PA)\mathbf x$。定理：$A$ 是 Hurwitz 的，当且仅当对任意 $Q\succ0$，李雅普诺夫方程

$$ A^{\top}P+PA=-Q $$

有唯一解且 $P\succ0$（此时 $P=\int_0^\infty e^{A^\top t}Qe^{At}dt$）。它是 $n(n+1)/2$ 个未知数的线性方程组。对挂着的摆，$A=\begin{bmatrix}0&1\\-a&-b\end{bmatrix}$（$a=g/l$），$Q=I$，逐项解出 $p_{12}=\tfrac1{2a}$，$p_{22}=\tfrac{1+2p_{12}}{2b}$，$p_{11}=a\,p_{22}+b\,p_{12}$。$b=0.5$ 时 $P=\begin{bmatrix}10.83&0.051\\0.051&1.10\end{bmatrix}$。把它用到非线性摆上，$\dot V$ 等于 $-\mathbf x^\top\mathbf x$ 加上高阶项，只在原点附近为负；FIG 2 搜出的 $c^*=6.49$ 对应的椭圆在 $\theta$ 方向只伸到 0.77 rad（约 44°）。这就是"用线性化证明"的典型代价：可靠，但保守。

可控性：秩判据、格拉姆矩阵与 PBH 判据

从 $\mathbf x_0=0$ 出发，$\mathbf x(T)=\int_0^Te^{A(T-\tau)}B\mathbf u(\tau)d\tau$。由凯莱–哈密顿定理，$A^n$ 是 $I,A,\dots,A^{n-1}$ 的线性组合，所以 $e^{At}B$ 永远落在 $\operatorname{span}\{B,AB,\dots,A^{n-1}B\}$ 里：能到达的状态不会超出 $\mathcal C$ 的列空间。反过来，若 $\mathcal C$ 满秩，下面的格拉姆矩阵可逆，就能显式构造出到达任意目标的输入：

$$ \begin{gathered} \hl{W}{W_T=\int_0^T e^{A\tau}BB^{\top}e^{A^{\top}\tau}\,d\tau},\qquad \mathbf u^*(t)=B^{\top}e^{A^{\top}(T-t)}\,W_T^{-1}\mathbf x_T,\\[6pt] \int_0^T\|\mathbf u^*\|^2dt=\mathbf x_T^{\top}W_T^{-1}\mathbf x_T . \end{gathered} $$

$\mathbf u^*$ 是所有到达 $\mathbf x_T$ 的输入中能量最小的一个（它是一个最小范数最小二乘解，参见 1.13 节的思路）。能量公式说明单位能量能到达的集合就是椭圆 $\mathbf x^\top W_T^{-1}\mathbf x\le1$，半轴是 $W_T$ 特征值的平方根。例：小车（双积分器）$A=\begin{bmatrix}0&1\\0&0\end{bmatrix}$，$B=\begin{bmatrix}0\\1\end{bmatrix}$，$e^{A\tau}B=(\tau,1)$，得 $W_T=\begin{bmatrix}T^3/3&T^2/2\\T^2/2&T\end{bmatrix}$。从静止移动距离 $d$ 再停下（$\mathbf x_T=(d,0)$），所需能量 $12d^2/T^3$：$T=2$ s 时 1.5，$T=1$ s 时 12（FIG 3 的任务）。时间减半，最高速度要翻倍，而且要在一半的时间里达到，所以力要变成 4 倍；能量 $\int u^2dt$ 于是变成 $16\times\tfrac12=8$ 倍，即 $1/T^3$。

PBH 判据（Popov–Belevitch–Hautus）给出几何意义：$(A,B)$ 不可控，当且仅当存在 $A$ 的左特征向量 $\mathbf w^\top A=\lambda\mathbf w^\top$ 满足 $\mathbf w^\top B=0$。也就是说，有一个模态 $\mathbf w^\top\mathbf x$ 按 $\dot{(\mathbf w^\top\mathbf x)}=\lambda\,\mathbf w^\top\mathbf x$ 自己演化，输入完全碰不到它。FIG 3 的 $A=\operatorname{diag}(-1,-3)$：$\det\mathcal C=\det[B,\ AB]=-2b_1b_2$，$B$ 落在某根坐标轴上时，另一个坐标就是碰不到的模态。若不可控的模态本身是稳定的（如这里），系统仍可镇定（stabilizable），这对 LQR 已经够用。

注意线性化不可控不等于非线性系统不可控。独轮车在静止处的线性化 $\dot{\mathbf q}=B(\theta)\mathbf u$ 只有两列，秩为 2 < 3，可 9.1 节用李括号证明了它可控：侧移要靠"前进–转向–后退–转回"这种二阶效应，线性化看不见。

极点配置与 Ackermann 公式

定理：单输入系统 $(A,B)$ 可控，则对任意 $n$ 次实系数首一多项式 $\varphi(s)$，存在唯一的 $K$ 使 $\det(sI-(A-BK))=\varphi(s)$；不可控时，不可控模态的特征值无论 $K$ 取什么都不会动。构造性的写法是 Ackermann 公式：

$$ \hl{K}{K}=\begin{bmatrix}0&\cdots&0&1\end{bmatrix}\mathcal C^{-1}\varphi(A),\qquad \varphi(s)=\prod_{i=1}^{n}(s-p_i) . $$

对 $B=(0,\dots,0,1)^\top$ 形式的"可控标准形"，$A-BK$ 的最后一行恰好就是特征多项式的系数，于是每个增益直接对应一个系数，这正是第 5 节两行推导的来源。一个熟悉的特例：8.3 节的计算力矩把每个关节变成双积分器 $\ddot e=a$，PD 律 $a=-K_pe-K_d\dot e$ 就是 $K=[K_p,K_d]$，把两个极点放在 $s^2+K_ds+K_p=s^2+2\zeta\omega_ns+\omega_n^2$ 的根上。多输入时满足要求的 $K$ 有无穷多个，多出来的自由度可以用来改善鲁棒性，或交给 LQR 去选。选极点的经验：欠阻尼极点对 $\sigma\pm j\omega_d$ 的 2% 调节时间大约是 $4/|\sigma|$；极点越靠左，增益越大（FIG 4 中 $k_1=g/l+p_1p_2$）。

离散化：零阶保持与欧拉法

数字控制器每隔 $T$ 秒读一次传感器、算一次 $u_k$，然后在整个周期里保持这个值（零阶保持，ZOH）。在一个周期内 $\mathbf u$ 是常数，把常数变易公式从 $kT$ 积到 $(k+1)T$，得到精确的离散模型：

$$ \begin{gathered} \mathbf x_{k+1}=\hl{Ad}{A_d}\,\mathbf x_k+B_d\,\mathbf u_k,\qquad A_d=e^{AT},\quad B_d=\int_0^Te^{A\sigma}d\sigma\,B,\\[6pt] \exp\!\Big(\begin{bmatrix}A&B\\0&0\end{bmatrix}T\Big)=\begin{bmatrix}A_d&B_d\\0&I\end{bmatrix}. \end{gathered} $$

右边的分块指数是实际计算时的技巧（FIG 5 就用它）。$A_d$ 的特征值是 $z=e^{\lambda T}$：左半平面整个映到单位圆内，与 1.11 节"离散时间 $|\lambda|<1$"的判据完全对上。前向欧拉法把导数换成差分 $\dot{\mathbf x}\approx(\mathbf x_{k+1}-\mathbf x_k)/T$：

$$ \hl{Ae}{A_d^{\mathrm E}=I+AT},\qquad B_d^{\mathrm E}=BT,\qquad z=1+\lambda T,\qquad |z|<1\iff\big|\lambda+\tfrac1T\big|<\tfrac1T . $$

欧拉法只是 $e^{\lambda T}$ 的一阶近似，它的稳定区只是左半平面里一个圆心在 $-1/T$、半径 $1/T$ 的圆盘。弱阻尼极点 $-\zeta\omega_n\pm j\omega_n\sqrt{1-\zeta^2}$ 给出 $|z|^2=1-2\zeta\omega_nT+\omega_n^2T^2$，所以只要 $T>2\zeta/\omega_n$，欧拉模型就会把一个稳定的振子预测成发散的；快速实极点 $\lambda$ 则要求 $T<2/|\lambda|$，这就是 2.5 节讲过的显式欧拉步长限制。其它常用的映射还有双线性（Tustin）变换 $z=\tfrac{1+\lambda T/2}{1-\lambda T/2}$，它把整个左半平面映进单位圆，常用于把连续设计的控制器转成离散的。

采样控制回路：被控对象是连续的，控制器 $u_k=-K\mathbf x_k$ 每个周期更新一次，真实的闭环是 $\mathbf x_{k+1}=(A_d-B_dK)\mathbf x_k$。保持器让控制总是"慢半拍"，大致相当于 $T/2$ 的延迟，所以按连续时间设计的 $K$ 在 $T$ 变大时会失稳，而且比欧拉模型预测的早。

可观性与龙伯格观测器

状态反馈要用到整个 $\mathbf x$，可传感器只给出 $\mathbf y=C\mathbf x$。能不能从 $\mathbf y(t)$ 和已知的 $\mathbf u(t)$ 还原出 $\mathbf x$？把 $\mathbf y,\dot{\mathbf y},\ddot{\mathbf y},\dots$ 写出来（$\mathbf u=0$ 时 $\mathbf y^{(k)}=CA^k\mathbf x$），答案是

$$ \begin{gathered} \mathcal O=\begin{bmatrix}C\\CA\\\vdots\\CA^{n-1}\end{bmatrix},\qquad \operatorname{rank}\mathcal O=n,\\[6pt] \dot{\hat{\mathbf x}}=A\hat{\mathbf x}+B\mathbf u+\hl{L}{L}\big(\mathbf y-C\hat{\mathbf x}\big),\qquad \dot{\mathbf e}=(A-LC)\,\mathbf e . \end{gathered} $$

$(A,C)$ 可观，当且仅当 $\mathcal O$ 的秩为 $n$；这与"$(A^\top,C^\top)$ 可控"完全等价（对偶性）。龙伯格观测器（1964）是一个并行运行的模型，用"测量减预测"的差把它往正确的方向拉。估计误差 $\mathbf e=\mathbf x-\hat{\mathbf x}$ 的动态与输入无关，只由 $A-LC$ 决定；按对偶性，可观时 $A-LC$ 的特征值可以任意配置，配置方法就是对 $(A^\top,C^\top)$ 做极点配置。摆的两种传感器：

编码器 $C=[1\ \ 0]$：$\mathcal O=I$，总是可观。把两个观测器极点都放在 $-p$：$L_1=2p-b$，$L_2=p^2-a-bL_1$（$a=\tfrac gl\cos\theta^*$）。
陀螺仪（只测角速度）$C=[0\ \ 1]$：$\mathcal O=\begin{bmatrix}0&1\\-a&-b\end{bmatrix}$，$\det\mathcal O=a$。有重力作用在角度上（$a\ne0$）时可观：角度改变重力力矩，重力力矩又改变角速度，所以角速度里"藏着"角度。关节在水平面里转动时 $a=0$，不可观：给角度加任何常数偏移，陀螺仪读数都不变。这正是 IMU 只靠陀螺仪无法知道航向（偏航角）的原因（15.1 节）。

观测器极点越快，估计收敛越快，但 $L$ 越大，测量噪声也被放大得越多。4.8 节的卡尔曼滤波正是按噪声的统计量来选 $L$ 的观测器。把观测器和状态反馈 $\mathbf u=-K\hat{\mathbf x}$ 接在一起，闭环特征值恰好是 $A-BK$ 与 $A-LC$ 两组的并集（分离原理），可以分别设计；13.3 节的 LQG 会再回到这一点。

7机器人中的它：控制器跑在计算机上，传感器只看到一部分

离散化：同一个系统，两种离散模型

机器人的控制器以固定周期运行：力控机械臂常见 1 kHz，移动底盘的速度环可能只有几十 Hz，视觉伺服受相机帧率限制。模型预测控制（13.4 节）和 LQR 的离散版本都需要一个离散模型 $\mathbf x_{k+1}=A_d\mathbf x_k+B_d\mathbf u_k$。下图比较两种得到它的方法：精确的零阶保持（琥珀色）和"把导数换成差分"的欧拉法（紫色）。线是真实系统，点是两个离散模型每隔 $T$ 秒的预测。

FIG 5零阶保持 vs 欧拉：z 平面与采样响应

z 平面（单位圆 = 稳定边界）

左：琥珀色实心点 = 零阶保持模型的特征值 $e^{\lambda T}$，紫色空心点 = 欧拉模型的特征值 $1+\lambda T$；细的琥珀色曲线和紫色直线是 $T$ 从 0 增大到 0.3 s 时它们走过的轨迹；落到窗口外的点画在边缘并标出 $|z|$。右：青色实线 = 真实连续系统（RK4，0.5 ms），琥珀色实心点 = 零阶保持模型的递推，紫色空心点 = 欧拉模型的递推（最多画 60 个点，$T$ 小时每隔几个样本画一个），纵轴截断在 $\pm2$。"数字控制"预设里，青色线是真实的采样控制回路（连续的摆 + 每 $T$ 秒更新一次、其间保持不变的 $u=-K\mathbf x$，$K$ 把连续极点放在 $-5$ 重根），虚线是理想的连续时间闭环。

试一试

"柔性关节"（$\omega_n=6$ rad/s，$\zeta=0.1$）：加大 $T$，直到欧拉模型的极点跑出单位圆（$T>2\zeta/\omega_n\approx0.033$ s），而零阶保持模型和真实系统都依然稳定。
"数字控制"：加大控制周期 $T$，直到真实的采样回路失稳（琥珀色极点出圆，约 $T\ge0.211$ s），而欧拉模型还说"稳定"（它要到 $T=0.4$ s 才失稳）。

"刚性系统"预设里，快模态 $\lambda=-40$（比如电机电流）让欧拉模型在 $T>0.025$ s 时开始正负交替地振荡，$T>0.05$ s 时发散，而它对应的真实运动在几十毫秒内就衰减完了。所以：用于预测和优化的模型，用精确离散化（MPC 求解器通常用 RK4 或隐式积分器，不用一步欧拉）；控制周期要比闭环的时间常数短得多。FIG 5 的设计是 $-5$ 重根（时间常数 0.2 s），采样周期 $T$ 到 0.211 s，也就是大约一个闭环时间常数时，回路就失稳了。

观测器：从编码器还原角速度

编码器只给角度。直接对角度做差分得到的速度噪声很大（噪声除以一个很小的 $T$）。观测器换个思路：让一个模型和摆一起跑，用测量去修正它。下图中摆受已知的力矩 $u(t)=2\sin(1.3t)$ N·m 驱动，真实初始角 0.4 rad，观测器从 0 开始猜。选传感器和安装方式，调观测器极点和测量噪声。

FIG 6龙伯格观测器：能不能从测量还原状态

角度 θ (rad)

角速度 ω (rad/s) 对时间 (s)

琥珀色实线 = 真实状态，紫色虚线 = 观测器估计 $\hat{\mathbf x}$，青色小点 = 带噪声的测量（画在被测的那个量上）。右上：琥珀色摆 = 真实，紫色空心圆 = 估计的摆锤位置。小角度线性模型 $\ddot\theta=-a\theta-b\dot\theta+u$，$b=0.3$；"挂着"时 $a=g/l=9.81$，"水平面内"时 $a=0$。两个观测器极点都放在 $-p$（不可观时只能放可观的那一个，另一个停在 0）。系统和观测器一起用 RK4 积分，步长 1 ms，每步取一个新的噪声样本；"ω 估计 RMS"是 $t>2$ s 后 $\omega-\hat\omega$ 的均方根。

试一试

选"陀螺仪 ω"和"水平面内"，让观测器跑满 5 s：无论 $p$ 多大，角度误差都还停在初始值的 90% 以上。换成"挂着"，同样的陀螺仪为什么又能估出角度了？
用编码器，把测量噪声调到 $\sigma\ge0.02$ rad，再把观测器极点推到 $p\ge25$：ω 估计的 RMS 超过 0.03 rad/s。和 $p=5$ 比一比，快的观测器付出了什么？

真实系统里

在平衡点附近线性化是几乎所有经典机器人控制器的第一步：四旋翼在悬停点线性化后，小倾角下水平加速度约等于 $g$ 乘倾角，位置与姿态成为一串积分器，于是有了 13.1 节里的级联结构；平衡车和小车倒立摆在竖直点线性化后用状态反馈或 LQR（8.1 节的自动平衡）；9.2 节的线性倒立摆 $\ddot x=\tfrac{g}{z_c}(x-p)$ 是一个鞍点，特征值 $\pm\omega$（$\omega=\sqrt{g/z_c}$），"捕获点" $x+\dot x/\omega$ 正是它不稳定模态的坐标。机械臂经计算力矩变成双积分器后，$K_p=\omega_n^2,K_d=2\zeta\omega_n$ 就是极点配置（8.3 节）。工具：python-control（ss、ctrb、obsv、place、lyap、sample_system），SciPy（signal.place_poles、signal.cont2discrete 的 'zoh' 和 'euler' 方法、linalg.solve_continuous_lyapunov、linalg.expm），MATLAB 的 ss、ctrb、place、lyap、c2d；Drake 的 Linearize 给出系统在平衡点的 $(A,B,C,D)$，MuJoCo 的 mjd_transitionFD 用有限差分直接给出离散时间的 $A,B$。李雅普诺夫证书的自动化：平方和（SOS）规划可以搜索多项式 $V$，LQR-Trees（Tedrake 等，2010）用它验证每个 LQR 控制器的吸引域并拼成覆盖整个状态空间的"漏斗"；也有用神经网络学 $V$ 并验证的工作（Chang 等，2019）。观测器：许多伺服驱动器用观测器而不是差分来估计速度；卡尔曼滤波（4.8 节）是按噪声统计最优的观测器，视觉惯性里程计里全局位置和航向不可观（15.1 节）。

8常见误区

误区 1："线性化稳定"就是"稳定"

特征值只担保局部稳定。FIG 1 里离开平衡点 60° 线性化误差已超过 15%；FIG 2 里线性化导出的证书只覆盖约 44°；FIG 4 里闭环极点都在左半平面，力矩一受限，摆照样摔倒。设计完线性控制器，一定要在非线性模型上（带饱和、带延迟）仿真，并想清楚它要在多大的范围内工作。

误区 2：特征值决定一切

对固定的 $A$，特征值决定最终是否衰减，却不管中途。$A=\begin{bmatrix}-1&20\\0&-2\end{bmatrix}$ 的特征值是 $-1,-2$，可 $\|e^{At}\|$ 在 $t\approx0.69$ s 时达到约 5.0：状态先放大 5 倍再衰减（非正规矩阵的"瞬态增长"），足以撞上饱和或关节限位。对时变系统更糟：Markus–Yamabe（1960）的例子里 $A(t)$ 在每一时刻的特征值都是 $-0.25\pm0.66j$，解却按 $e^{t/2}$ 增长（10 s 后放大约 148 倍）。"冻结时间看特征值"对时变或非线性系统不是证明。

误区 3：李雅普诺夫候选失败 = 不稳定

李雅普诺夫定理是充分条件。FIG 2 的朴素二次型在稳定的摆上也会上升，那只说明这个 $V$ 选错了，什么也没证明。反过来，$\dot V\le0$ 只给出"稳定"，要得到"渐近稳定"还需要严格不等式或 LaSalle 原理；摆的能量在 $b=0$ 时 $\dot V\equiv0$，摆就永远摆下去，这时它稳定但不渐近稳定。要证明不稳定，需要另外的定理（如 Chetaev 定理）或一个正实部的线性化特征值。

误区 4："可控"就是"好控"，极点越快越好

秩判据是非黑即白的，可控程度却是连续的。FIG 3 的"两个一阶模态"里，默认的 $B$（45°）到达默认目标需要能量 0.81；$B$ 离特征方向只差 10° 时涨到约 14.4，差 5.1°（刚好在吸附范围外）时约 70。极点也一样：能放不代表该放，重极点 $-p$ 的增益按 $p^2$ 增长，很快撞上力矩上限、放大传感器噪声，并让采样延迟变得致命（FIG 4、FIG 5）。看 $W_T$ 的最小特征值或 $\mathcal C$ 的条件数，而不只看秩。

误区 5：用欧拉法离散化，或者忘了保持器

$A_d=I+AT$ 写起来最省事，但它会把稳定的弱阻尼振子预测成发散（FIG 5 中 $T>0.033$ s），把快模态预测成正负跳动；放进 MPC 里，优化器会认真地"利用"这些假动态。另一方面，即使模型精确，真实的数字回路也比连续设计更不稳定，因为保持器带来约 $T/2$ 的延迟：FIG 5 的控制器在 $T=0.211$ s 失稳，而欧拉模型说能撑到 0.4 s。离散化用 expm 或 cont2discrete(..., method='zoh')，并在真实的采样周期下检查 $A_d-B_dK$ 的特征值。

9小练习

1. 无阻尼的摆（$b=0$，$g/l=9.81\ \text{s}^{-2}$）在倒立点 $\theta^*=\pi$ 处线性化，特征值是什么？这个平衡点是什么类型？

$A=\begin{bmatrix}0&1\\ g/l&0\end{bmatrix}$（因为 $\cos\pi=-1$），特征方程 $\lambda^2-g/l=0$，$\lambda=\pm\sqrt{9.81}=\pm3.13$：一正一负，鞍点。选项 A 是挂着的点（$\cos0=1$）的答案。正特征值意味着倒立点不稳定，由李雅普诺夫间接法，这个结论对非线性摆也成立。

2. 有阻尼的摆（$b>0$）的能量满足 $\dot E=-b\omega^2\le0$，在 $\omega=0$ 的整条线上等于零。从 $E<2mgl$、$|\theta|<\pi$ 出发，能得出什么结论？

集合 $\{E<2mgl,|\theta|<\pi\}$ 有界且正不变（能量不增，又不够翻过顶点）。在 $\{\dot E=0\}=\{\omega=0\}$ 中能一直待住的解要求 $\omega\equiv0$，从而 $\dot\omega=-(g/l)\sin\theta=0$，只剩 $\theta=0$。LaSalle 原理于是给出收敛到 $\theta=0$。选项 B 是只用李雅普诺夫定理本身（不用 LaSalle）时能得到的较弱结论。

3. $A=\begin{bmatrix}-1&0\\0&-3\end{bmatrix}$，$B=\begin{bmatrix}1\\0\end{bmatrix}$。下列哪个说法正确？

$AB=(-1,0)^\top$，$\mathcal C=\begin{bmatrix}1&-1\\0&0\end{bmatrix}$，秩为 1：不可控。PBH：左特征向量 $\mathbf w=(0,1)$ 满足 $\mathbf w^\top B=0$，对应的模态 $\dot x_2=-3x_2$ 输入碰不到。但它本来就稳定，所以任何 $K$ 都只能动 $x_1$ 的极点，系统仍然可以被镇定。稳定性和可控性是两回事（选项 A），$B\ne0$ 也不够（选项 D）。

4. 双积分器 $\ddot q=u$，状态 $\mathbf x=(q,\dot q)$，用 $u=-K\mathbf x$ 把两个闭环极点都放在 $-2$。$K$ 是多少？

闭环 $\ddot q+k_2\dot q+k_1q=0$，特征多项式 $s^2+k_2s+k_1$，要等于 $(s+2)^2=s^2+4s+4$，所以 $k_1=4$，$k_2=4$。这就是 $\omega_n=2$、$\zeta=1$ 的 PD：$K_p=\omega_n^2$，$K_d=2\zeta\omega_n$。

5. 一个稳定模态 $\lambda=-50\ \text{s}^{-1}$，采样周期 $T=0.05$ s。零阶保持模型和欧拉模型的离散特征值分别是多少？

零阶保持是精确的：$z=e^{\lambda T}=e^{-2.5}\approx0.082$，在单位圆内。欧拉 $z=1+\lambda T=-1.5$，$|z|>1$：一个真实中很快衰减的模态，在欧拉模型里正负交替地发散。欧拉稳定要求 $T<2/|\lambda|=0.04$ s。

10延伸阅读

K. J. Åström & R. M. Murray, Feedback Systems: An Introduction for Scientists and Engineers, Princeton University Press, 2008 (2nd ed. 2021) — 动态行为、线性系统、状态反馈与输出反馈的几章是本节最好的入门读物（在线免费）
H. K. Khalil, Nonlinear Systems, 3rd ed., Prentice Hall, 2002, Ch. 4 — 李雅普诺夫稳定性、LaSalle 原理、线性化与吸引域估计
J.-J. E. Slotine & W. Li, Applied Nonlinear Control, Prentice Hall, 1991 — 相平面分析与李雅普诺夫方法，面向工程师
G. F. Franklin, J. D. Powell & M. L. Workman, Digital Control of Dynamic Systems, 3rd ed., Addison-Wesley, 1998 — 零阶保持离散化、采样率选择与数字控制器设计
R. Tedrake, Underactuated Robotics, MIT course notes — 摆、小车倒立摆、李雅普诺夫分析与平方和验证的机器人视角（在线免费）
A. M. Lyapunov, The General Problem of the Stability of Motion, 1892; English translation by A. T. Fuller, Int. J. Control 55(3), 1992 — 直接法与间接法的出处
R. E. Kalman, On the General Theory of Control Systems, Proc. 1st IFAC Congress, Moscow, 1960 — 可控性与可观性的概念及对偶性
D. G. Luenberger, Observing the State of a Linear System, IEEE Trans. Military Electronics 8(2), 1964 — 龙伯格观测器
R. Tedrake, I. R. Manchester, M. Tobenkin & J. W. Roberts, LQR-Trees: Feedback Motion Planning via Sums-of-Squares Verification, Int. J. Robotics Research 29(8), 2010 — 用 SOS 验证的吸引域把局部线性控制器拼成全局控制器