130_鲁棒与随机MPC

专题 3.13：鲁棒与随机 MPC — 从名义模型走向真实世界¶

博士前数学路线图 · 第三批 · Phase E MPC 收官专题 前置：3.1–3.12（ISS/Lyapunov、CBF、MPC 四条件/Tube/ISS-MPC/Economic MPC、数值求解）面向：机器人算法工程师（RL motion control、具身智能、SLAM），C++/Python 时间预算：档位 3（8 节必学）约 40–55 小时；档位 4（A–F 选学）按需再加 25–40 小时

摘要与本专题定位¶

名义 MPC 的一切稳定性论证，建立在"模型就是真相"这一谎言之上。 专题 3.13 把 3.11 的四条件（terminal set, terminal cost, stage cost, horizon）推到"风会吹、摩擦会变、SLAM 会漂、质量没称准"的真实机器人世界。核心工具有三件：Tube MPC 家族（把不确定性塞进 Minkowski 管道）、Stochastic/Scenario MPC（机会约束 $\Pr(g\le 0)\ge 1-\varepsilon$）、Distributionally Robust MPC（Wasserstein 球内 worst-case）。再加上**学习型鲁棒 MPC**（GP tube、safety filter），构成机器人 MPC 在 2026 年前后的完整工程谱系。

作为第三批 13 个子专题（Phase E）的收官，本专题回答两个问题：(i) 给定一个 3.12 能跑实时的 NMPC，如何给它加一层安全外衣，使其在真实机器人上不出事？(ii) 这层外衣与第六批要学的 RL/MBRL、第四批的刚体动力学、第五批的 SLAM 估计之间，是什么关系？

BLUF：在机器人场景下，如果不确定性有确定边界（质量 ±10%、风 ≤5 m/s），首选 Rigid/Homothetic Tube MPC；如果扰动是概率性的并且你愿意接受 $\varepsilon=5\%$ 的违反率，用 Scenario MPC 或 analytic chance constraint；如果你只有数据没有模型，用 DeePC + L2 正则（隐式 DRO）；如果你在 sim-to-real 场景，tube MPC + domain randomization 是当前最稳妥的组合。Min-max MPC 虽数学优美，工程上几乎从不直接用，仅作为理论对照物。

§3.13.1 不确定性的分类与建模¶

定义（三类不确定性）¶

考虑离散时间系统 $$x_{k+1}=f(x_k,u_k,\theta,w_k),\qquad y_k=h(x_k,\theta)+v_k.$$

加法扰动：$w_k\in\mathcal{W}\subset\mathbb{R}^{n_w}$，$\mathcal{W}$ 紧凸包含原点。典型线性情形 $x_{k+1}=Ax_k+Bu_k+w_k$。

参数/乘法不确定性：$\theta\in\Theta\subset\mathbb{R}^{n_\theta}$，$\Theta$ 紧凸（通常多面体或椭球）。LPV 视图 $x_{k+1}=A(\theta)x_k+B(\theta)u_k$。

分布不确定性：$w_k\sim \mathbb{P}$，$\mathbb{P}\in\mathcal{P}_{\mathrm{amb}}$，$\mathcal{P}_{\mathrm{amb}}$ 为**模糊集 (ambiguity set)**。常见构造：矩约束 $\{\mathbb{Q}:\mathbb{E}[w]=\mu,\,\mathrm{Cov}[w]\preceq\Sigma\}$、Wasserstein 球 $\mathbb{B}_\varepsilon(\hat{\mathbb{P}}_N)=\{\mathbb{Q}:W_p(\mathbb{Q},\hat{\mathbb{P}}_N)\le\varepsilon\}$、KL 球 $\{\mathbb{Q}:\mathrm{KL}(\mathbb{Q}\|\hat{\mathbb{P}}_N)\le\eta\}$。

机器人物理含义映射表¶

不确定来源	数学建模	典型量级（机器人）
风扰（四旋翼）	加法 $w\in\mathcal{W}$，$\\|w\\|_\infty\le 5\,\mathrm{m/s}$ 对应加速度扰动	3–8 m/s²
地面摩擦 $\mu$	参数 $\theta\in[0.4,1.0]$	±50%
负载质量 $m$	参数 $\theta=m\in[m_0,m_0+\Delta m]$	±20%
SLAM 位姿误差	加法 $v_k$ 或分布 $\mathbb{P}\in\mathcal{P}_{\mathrm{amb}}$	位置 1–5 cm、yaw 0.5–2°
接触位置不确定	参数（接触点坐标）+ 加法（法向力扰动）	接触点 ±1 cm
执行器延迟	参数（延迟 $\tau$）	5–50 ms

核心观察：物理不确定性**很少是单一类型**——四足机器人"地形变化"同时包含摩擦（参数）、接触位置（参数）、外部推力（加法）。务实做法是**先分解再叠加**：把主导项归到一个可处理的类型，剩余项作为有界加法残差打包。

不确定性量化的工程方法论¶

问题：如何从一个真实机器人上获取 $\mathcal{W}$（扰动集）或 $\Theta$（参数不确定性集）的定量估计？这是整个鲁棒 MPC 链条中**最工程化、最不"数学"的一步**——却往往是成败的关键。

方法 1：物理分析法 - 风扰动：查询当地气象数据，阵风速度 $v_g \le 5$ m/s → 加速度扰动 $w_a = \frac{1}{2}\rho C_D A v_g^2 / m$ - 负载变化：已知工件质量范围 $m \in [m_0 - \Delta m, m_0 + \Delta m]$ - 执行器延迟：从数据手册或 step response 实验测量

方法 2：数据驱动法 - 收集闭环数据 $\{x_k, u_k, x_{k+1}\}$ - 计算残差 $e_k = x_{k+1} - f_{\text{nom}}(x_k, u_k)$ - $\mathcal{W}$ 取为 $e_k$ 的经验分布的 $(1-\delta)$ 置信凸包 - 更精细：区分状态依赖残差（参数不确定）和状态独立残差（加法扰动）

方法 3：蒙特卡洛仿真法 - 在仿真器中随机化所有不确定参数 - 跑大量轨迹，收集状态偏差统计 - 用 $3\sigma$ 规则或 Chebyshev 界构造 $\mathcal{W}$

方法 4：GP 后验法（§3.13.8 详述） - 用 GP 学习残差，$\sigma_g(x,u)$ 自动给出状态依赖的不确定性估计 - $\mathcal{W}(x,u) = \{w: \|w\| \le \beta_\delta \sigma_g(x,u)\}$——自适应扰动集

方法	精度	适用阶段	缺点
物理分析	粗略	初始设计	容易低估
数据驱动	中等	系统已运行	需要激励数据
蒙特卡洛	高	仿真阶段	仿真本身有误差
GP 后验	自适应	在线学习	计算复杂度高

工程经验法则：初始设计用物理分析（乘以安全系数 1.5-2），仿真验证用蒙特卡洛，实机部署用数据驱动方法动态更新。永远不要只用一种方法——交叉验证是唯一靠谱的做法。

从名义 MPC 到鲁棒 MPC 的决策流程¶

你的系统有不确定性吗？
├── 否 → 使用名义 MPC (3.11)
└── 是 → 不确定性有确定性边界吗？
    ├── 是 (|w| ≤ w_max) → 不确定性维度/量级？
    │   ├── 小扰动 → Rigid Tube MPC (§3.13.3)
    │   ├── 大扰动 → Homothetic/Elastic Tube (§3.13.3.3-4)
    │   └── 参数不确定 → LPV Tube MPC (§3.13.4)
    └── 否 (概率性扰动)
        ├── 知道分布（如 Gaussian） → 解析 Chance Constraint (§3.13.5)
        ├── 只能采样 → Scenario MPC (§3.13.6)
        └── 只有少量数据 → Wasserstein DRO (§3.13.7)

+ 学习残差？ → GP-MPC / Safety Filter (§3.13.8)
+ 无模型？ → DeePC (§3.13.D)
+ RL 策略安全部署？ → PSF + Tube (§3.13.E)

练习¶

练习 3.13.1.1 ⭐：对你的实验室机器人（或一个仿真模型），列出所有不确定性源，分类为加法/参数/分布型。估计每种不确定性的量级。

练习 3.13.1.2 ⭐⭐：设计一个实验方案，从实际数据中估计加法扰动集 $\mathcal{W}$。讨论：用多长的数据？如何处理异常值？如何验证估计的可靠性？

桥梁¶

→ §3.13.3/4/8：加法 → Tube；参数 → LPV-Tube / Adaptive；分布 → Stochastic/DRO。
→ 第四批刚体动力学：$M(q)\ddot q+C(q,\dot q)\dot q+g(q)=\tau+\tau_{\mathrm{ext}}$ 中 $\tau_{\mathrm{ext}}$ 是典型加法扰动。
→ 第五批 SLAM：SLAM 估计误差协方差 $\Sigma_{\hat x}$ 直接构造 $\mathcal{W}$ 或 Wasserstein 半径 $\varepsilon_N$。

§3.13.2 Min-max MPC：数学上的顶峰，工程上的陷阱¶

严格定义¶

Open-loop min-max： $$\min_{\mathbf{u}=(u_0,\dots,u_{N-1})}\,\max_{\mathbf{w}\in\mathcal{W}^N}\,J(x_0,\mathbf{u},\mathbf{w})\quad \text{s.t.}\ x_{k+1}=Ax_k+Bu_k+w_k,\,x_k\in\mathcal{X},\,u_k\in\mathcal{U},\ \forall\mathbf{w}\in\mathcal{W}^N.$$

Closed-loop (feedback) min-max（Scokaert–Mayne 1998）： $$\min_{\pi=\{\mu_0,\dots,\mu_{N-1}\}}\,\max_{\mathbf{w}\in\mathcal{W}^N}\,J(x_0,\pi,\mathbf{w}),\quad u_k=\mu_k(x_k).$$

关键区别：开环 min-max 要求**一条固定输入序列**对**所有扰动**可行；闭环 min-max 允许**每步根据已实现扰动做决策**（recourse），保守性大幅降低甚至可行性大幅扩张。

核心定理（Scokaert–Mayne 1998，IEEE TAC 43(8):1136–1142，DOI:10.1109/9.704989）¶

定理 3.13.1（反馈 min-max MPC 稳定性）。设 $(A,B)$ 可镇定，$\mathcal{W}$ 凸紧，终端约束集 $\mathcal{X}_f$ 为鲁棒正不变集，$V_f$ 为其上的鲁棒 Lyapunov 函数。则反馈 min-max MPC 律 $\kappa_N(x)=\mu_0^*(x)$ 使原点 ISS 稳定，且 $\forall k,\,x_k\in\mathcal{X},\,u_k\in\mathcal{U}$。

证明骨架：(i) 终端集上存在局部鲁棒控制器 $\kappa_f$ 使 $V_f(f(x,\kappa_f(x))+w)-V_f(x)\le -\ell(x,\kappa_f(x)),\forall w\in\mathcal{W}$；(ii) 取候选策略 $\pi^+=(\mu_1^*,\dots,\mu_{N-1}^*,\kappa_f)$，对任意 $w_0$ 仍可行；(iii) 最优值函数 $V_N^*$ 满足 Lyapunov 下降 $V_N^*(x^+)-V_N^*(x)\le -\ell(x,\kappa_N(x))+\alpha(\|w_0\|)$。

计算复杂度爆炸¶

若 $\mathcal{W}$ 为凸多面体含 $q$ 个顶点，最坏情形仅在顶点处出现（$J$ 对 $w$ 凸）。对闭环 min-max，每步扰动可取 $q$ 个分支 ⇒ 长度 $N$ 的 scenario tree 含 $q^N$ 条路径，每条路径独立决策变量 $\Rightarrow$ 变量数 $\mathcal{O}(q^N)$。

$N$	$q=3$（1D 扰动）	$q=8$（3D box 扰动）	$q=64$（6D box）
5	243	32,768	$10^9$
10	59,049	$10^9$	$10^{18}$
20	$3\times 10^9$	$10^{18}$	$10^{36}$

对机器人 $N\ge 20$、扰动维度 $n_w\ge 3$ 的场景，闭环 min-max 完全不可行。

Bemporad–Borrelli–Morari 2003（mp-LP 方法）¶

对 $\ell_1/\ell_\infty$ 代价 + 多面体不确定，**显式解**是 $x_0$ 的分段仿射函数 $u^*(x_0)=F_i x_0+g_i,\,x_0\in\mathcal{R}_i$（TAC 48(9):1600，DOI:10.1109/TAC.2003.816984）。但分区数仍随 $N$ 和顶点数指数增长，仅对低维（$n_x\le 4$）嵌入式系统可行。

工程含义¶

Min-max MPC 的工程价值：(i) 作为上界——任何 tube/scenario 方法都不会超过 min-max 的保守性；(ii) 离线分析——用 mp-LP 计算低维系统的精确鲁棒可行域；(iii) 理论地基——recourse 思想直接演化为 tube MPC 的 $K$ 反馈律和 multi-stage MPC 的 scenario tree。直接部署 ≈ 学术自杀。

为什么 min-max 在工程上不实用——深层原因分析¶

原因 1：保守性悖论。开环 min-max 要求一条固定输入序列对所有扰动可行。但物理直觉告诉我们：合理的控制器应该根据实际观测到的扰动做出调整。不允许 recourse（反馈）的控制器，就像一个蒙着眼睛走迷宫的人——必须对每个路口的每种可能性都预留余量。这种余量累积起来，导致可行域急剧缩小。

原因 2：维度爆炸的本质。闭环 min-max 允许反馈（recourse），但必须对**所有可能的扰动历史**提供一个决策。如果每步扰动有 $q$ 种可能值，$N$ 步后有 $q^N$ 种历史——对应 $q^N$ 个独立决策变量。这不是算法不够好的问题，而是**问题本身的信息论复杂度**：你需要为每种可能的信息状态准备一个响应。

原因 3：鞍点可能不存在。min-max 要求 $\min_u \max_w$ 有鞍点（即 $\min_u \max_w = \max_w \min_u$）。对非凸非线性系统，鞍点不一定存在——此时 min-max 问题本身定义不良。

不是 X 而是 Y：Min-max MPC 不是"最好的鲁棒 MPC"——它是"最保守且最贵的鲁棒 MPC"。Tube MPC 的核心洞察是：与其为每种扰动历史准备完整响应（$q^N$ 复杂度），不如用一个**简单但有效的线性反馈** $K$ 把误差压缩到小集合 $Z$ 内（$O(1)$ 复杂度），然后在收紧的约束内做名义优化。这是一种"用简单反馈换取计算效率"的工程哲学。

Bemporad–Borrelli–Morari 2003 的工程启示¶

虽然 explicit MPC / mp-LP 方法在高维不可用，但它揭示了一个重要结构：鲁棒 MPC 的闭环律是分段仿射 (piecewise affine) 函数。这意味着： - 在每个多面体区域 $\mathcal{R}_i$ 内，最优控制是 $u^*(x) = F_i x + g_i$（仿射） - 区域边界对应**活动约束的切换**（active set 变化） - 在线实现只需确定当前状态属于哪个区域（point location），然后查表

对低维嵌入式系统（$n_x \le 4$），这仍然是实用方案——例如直流电机的位置/速度控制、简单温控系统等。

桥梁¶

→ §3.13.3：Tube MPC 就是"用一个固定反馈 $K$ 作为 $\mu_k$ 的参数化"——把 $q^N$ 策略压成一个 $K$ 矩阵。
→ §3.13.C（多阶段鲁棒）：multi-stage 在 scenario tree 前 $N_r$ 步 branch，之后冻结——min-max 与 tube 之间的工程折中。

练习¶

练习 3.13.2.1 ⭐⭐：对一维系统 $x_{k+1} = x_k + u_k + w_k$，$\mathcal{X} = [-5,5]$，$\mathcal{U} = [-1,1]$，$\mathcal{W} = [-0.3,0.3]$，$N=3$，写出开环 min-max 的所有约束（enumerate $w$ 的顶点 $\pm 0.3$）。计算总约束数量，验证与 $q^N = 2^3 = 8$ 条路径的对应。

练习 3.13.2.2 ⭐⭐⭐：证明：对线性系统 + 凸代价 + 多面体约束，开环 min-max 的最坏情况代价 $\max_w J(u,w)$ 的最大值一定在 $\mathcal{W}$ 的顶点处取得（利用 $J$ 对 $w$ 的凸性）。

§3.13.3 Tube MPC 完整理论：机器人鲁棒控制的主力¶

核心思想：不确定性经线性反馈 $K$ 压缩为一个**误差管道** $Z\subset\mathbb{R}^{n_x}$（Minkowski 意义上的"误差不变集"），名义轨迹 $\bar x_k$ 在**收紧的约束** $\mathcal{X}\ominus Z$、$\mathcal{U}\ominus K Z$ 上优化，真实轨迹 $x_k=\bar x_k+e_k$ 自动满足 $x_k\in\mathcal{X}$。

3.13.3.1 预备：Pontryagin 差与 mRPI 集¶

不变集理论的来龙去脉¶

在正式进入 Tube MPC 之前，必须理解一个更基本的概念：不变集 (invariant set)。这是鲁棒控制中最核心的工具之一，其历史可以追溯到 1960 年代 Blanchini 和 Miani 的正不变集理论。

为什么需要不变集？ 考虑一个简单的问题：如果系统受到持续扰动 $w_k$，那么状态会不会无限漂移？如果不会，它最终会稳定在哪个区域内？不变集正是回答这个问题的工具——它定义了"无论扰动怎么作用，状态都不会离开的区域"。

正不变集 (Positively Invariant Set, PI)：集合 $\Omega \subseteq \mathbb{R}^{n_x}$ 对系统 $x_{k+1} = f(x_k)$ 是正不变的，当且仅当

\[x_k \in \Omega \Rightarrow x_{k+1} = f(x_k) \in \Omega\]

鲁棒正不变集 (Robust Positively Invariant Set, RPI)：集合 $\Omega$ 对不确定系统 $x_{k+1} = f(x_k, w_k)$, $w_k \in \mathcal{W}$ 是鲁棒正不变的，当且仅当

\[x_k \in \Omega \Rightarrow f(x_k, w_k) \in \Omega, \quad \forall w_k \in \mathcal{W}\]

与 Lyapunov 稳定性的深层联系：不变集可以看作 Lyapunov 函数水平集的离散版本。回顾 3.7 节 ISS 理论：如果存在 ISS-Lyapunov 函数 $V(x)$ 使得 $V(f(x,w)) - V(x) \le -\alpha(\|x\|) + \sigma(\|w\|)$，则任何水平集 $\{x: V(x) \le c\}$，当 $c$ 足够大时（$c \ge \alpha^{-1}(\sigma(\|w\|_{\max}))$），都是 RPI 集。这不是巧合——不变集理论和 Lyapunov 理论是同一枚硬币的两面。

本质洞察：最小 RPI 集 $Z_\infty$ 描述的是"在最坏情况下，误差最终会落入的最小区域"。Tube MPC 的全部工程智慧就在于：先计算这个最小区域，然后在设计名义控制时预留出这个区域的空间（约束收紧），从而保证真实状态永远不违反原始约束。

Pontryagin（Minkowski）差的完整推导¶

Pontryagin（Minkowski）差： $$A\ominus B:=\{x:x+B\subseteq A\}=\bigcap_{b\in B}(A-b).$$

为什么需要这个运算？ 直觉上，如果名义状态 $\bar x$ 在集合 $A \ominus B$ 内，而误差 $e$ 在集合 $B$ 内，那么真实状态 $x = \bar x + e$ 一定在集合 $A$ 内。这正是约束收紧的数学本质。

逐步推导 H-表示下的 Pontryagin 差：

设 $A = \{x : Fx \le g\}$（H-多面体，$F \in \mathbb{R}^{m \times n}$，$g \in \mathbb{R}^m$），$B$ 为紧凸集。

Step 1. 由定义，$x \in A \ominus B$ 等价于 $\forall b \in B, x + b \in A$。

Step 2. 代入 $A$ 的定义：$\forall b \in B, F(x + b) \le g$。

Step 3. 分离：$\forall b \in B, Fx + Fb \le g$，即 $Fx \le g - Fb, \forall b \in B$。

Step 4. 取最坏情况（约束最紧）：第 $i$ 行约束为 $F_i x \le g_i - \max_{b \in B} F_i b$。

Step 5. 识别 $\max_{b \in B} F_i b = h_B(F_i^\top)$——这正是 $B$ 在方向 $F_i^\top$ 上的**支持函数**。

因此： $$A \ominus B = \{x : Fx \le g - h_B(F^\top)\}$$

其中 $h_B(F^\top) \in \mathbb{R}^m$ 是逐行计算支持函数的向量。

几何含义：每个约束面向内平移了 $h_B(F_i^\top)$ 的距离——恰好是集合 $B$ 在该约束法向方向上的"厚度"。如果 $B$ 是以原点为中心的球 $\{b: \|b\| \le r\}$，则 $h_B(f) = r\|f\|$，所有约束面统一内缩 $r \cdot \|F_i\|$。

支持函数的性质（后续反复使用）：

性质	表达式	含义
正齐次	$h_B(\lambda f) = \lambda h_B(f), \lambda \ge 0$	方向缩放
次可加	$h_B(f_1 + f_2) \le h_B(f_1) + h_B(f_2)$	三角不等式
Minkowski 和	$h_{A \oplus B}(f) = h_A(f) + h_B(f)$	支持函数线性叠加
线性变换	$h_{MB}(f) = h_B(M^\top f)$	矩阵变换

如果不用 Pontryagin 差会怎样？ 你必须在 MPC 中对每一步、每一个扰动实现都检查约束满足——这就退化为 min-max MPC 的指数复杂度。Pontryagin 差把"对所有扰动可行"的无穷维约束**离线压缩**为一次性的约束收紧，使在线问题规模不变。

最小鲁棒正不变集 (mRPI) 的完整理论¶

最小鲁棒正不变集 (mRPI)：对误差动力学 $e_{k+1}=A_K e_k+w_k,\,A_K=A+BK$ Schur（所有特征值严格在单位圆内）， $$Z_\infty=\bigoplus_{i=0}^\infty A_K^i \mathcal{W}.$$

为什么 mRPI 是无穷和？ 这来自递归分析。如果 $e_0 = 0$（初始误差为零），则： - $k=1$: $e_1 = w_0 \in \mathcal{W}$ - $k=2$: $e_2 = A_K w_0 + w_1 \in A_K \mathcal{W} \oplus \mathcal{W}$ - $k=n$: $e_n = \sum_{i=0}^{n-1} A_K^i w_{n-1-i} \in \bigoplus_{i=0}^{n-1} A_K^i \mathcal{W}$ - $k \to \infty$: $e_\infty \in \bigoplus_{i=0}^{\infty} A_K^i \mathcal{W} = Z_\infty$

因为 $A_K$ Schur，$\|A_K^i\| \to 0$ 指数衰减，所以这个无穷 Minkowski 和是收敛的（即 $Z_\infty$ 有界）。

mRPI 的"最小性"：$Z_\infty$ 是所有 RPI 集的子集——任何 RPI 集 $\Omega$ 都满足 $Z_\infty \subseteq \Omega$。这是因为 $Z_\infty$ 恰好是"从零初始误差出发，在最坏扰动序列下所有可能到达的误差"的集合。

关键问题：$Z_\infty$ 是无穷和，无法精确计算。如何近似？

这就是 Raković et al. 2005 的贡献——提供了一个有限步可计算的**外近似**，且近似精度可控。

定理 3.13.2（Raković–Kerrigan–Kouramas–Mayne 2005, IEEE TAC 50(3):406–410, DOI:10.1109/TAC.2005.843854）¶

mRPI 的 $(s,\alpha)$ 外近似。设 $A_K$ Schur，$\mathcal{W}$ 含原点紧凸。若存在 $s\in\mathbb{N}_+$、$\alpha\in[0,1)$ 使 $A_K^s\mathcal{W}\subseteq\alpha\mathcal{W}$，则 $$Z_\infty\subseteq F(s,\alpha):=\frac{1}{1-\alpha}\bigoplus_{i=0}^{s-1}A_K^i\mathcal{W},\qquad A_K F(s,\alpha)\oplus\mathcal{W}\subseteq F(s,\alpha).$$

完整证明：

Step 1. 将 $Z_\infty$ 分成前 $s$ 项和尾部： $$Z_\infty = \bigoplus_{i=0}^{s-1} A_K^i \mathcal{W} \oplus \bigoplus_{i=s}^{\infty} A_K^i \mathcal{W}$$

Step 2. 尾部利用收缩性递推。由 $A_K^s \mathcal{W} \subseteq \alpha \mathcal{W}$，得 $A_K^{s+j} \mathcal{W} \subseteq \alpha A_K^j \mathcal{W}$。更一般地，$A_K^{ks} \mathcal{W} \subseteq \alpha^k \mathcal{W}$，因此： $$\bigoplus_{i=s}^{\infty} A_K^i \mathcal{W} \subseteq \frac{\alpha}{1-\alpha} \bigoplus_{j=0}^{s-1} A_K^j \mathcal{W}$$

Step 3. 合并前 $s$ 项与尾部： $$Z_\infty \subseteq (1 + \frac{\alpha}{1-\alpha}) \bigoplus_{i=0}^{s-1} A_K^i \mathcal{W} = \frac{1}{1-\alpha} \bigoplus_{i=0}^{s-1} A_K^i \mathcal{W} = F(s,\alpha)$$

Step 4. 验证 RPI 性质：利用支持函数的线性性和收缩条件，$h_{A_K F \oplus \mathcal{W}}(f) \le h_F(f)$ 逐行成立。

数值实例（二维双积分器）：$A = \begin{bmatrix}1 & 0.1\\ 0 & 1\end{bmatrix}$，$B = \begin{bmatrix}0.005\\ 0.1\end{bmatrix}$，$K$ 使 $\text{eig}(A_K) = \{0.6, 0.7\}$，$\mathcal{W} = [-0.1, 0.1]^2$。

$s$	$\alpha^\circ(s)$	近似误差	说明
5	0.51	12%	过粗
10	0.26	0.8%	可用
15	0.13	0.1%	精确
20	0.07	$<0.01\%$	过精

反事实推理：如果直接截断前 $s$ 项 $F_s^{\text{trunc}} = \bigoplus_{i=0}^{s-1} A_K^i \mathcal{W}$，则 $F_s^{\text{trunc}}$ 不是 RPI 集（$A_K F_s^{\text{trunc}} \oplus \mathcal{W} \not\subseteq F_s^{\text{trunc}}$），失去递归可行性保证。$1/(1-\alpha)$ 的放大是必要代价。

算法：给定精度 $\varepsilon>0$，从 $s=0$ 起递增；每步计算 $\alpha^\circ(s)=\max_j h_\mathcal{W}((A_K^s)^\top f_j)/g_j$；当 $\alpha^\circ(s)/(1-\alpha^\circ(s))\cdot\|F_s\|\le\varepsilon$ 即停。典型 $s\in[5,30]$。

$K$ 选取的 trade-off 分析：$A_K$ 特征值越小（$K$ 越激进），$Z$ 越小——约束收紧量减少，名义 MPC 可行域更大。但 $K$ 越大，$KZ$ 也越大——$\mathcal{U} \ominus KZ$ 缩小，可能导致控制输入约束收紧过度甚至不可行。最优 $K$ 是两者的平衡：使 $\text{vol}((\mathcal{X} \ominus Z) \cap \{x: Kx \in \mathcal{U} \ominus KZ\})$ 最大化。实践中以 LQR 增益为起点，调 $Q/R$ 的相对权重迭代。

3.13.3.2 Rigid Tube MPC（Mayne–Seron–Raković 2005）¶

完整出处：Mayne, Seron, Raković, "Robust model predictive control of constrained linear systems with bounded disturbances", Automatica 41(2):219–224, 2005, DOI:10.1016/j.automatica.2004.08.019。

构造：名义 $\bar x_{k+1}=A\bar x_k+B\bar u_k$，反馈 $u_k=\bar u_k+K(x_k-\bar x_k)$，误差 $e_{k+1}=A_Ke_k+w_k\in Z$（$Z$ 为 mRPI 外近似）。

在线 OCP： $$\min_{\bar x_0,\bar u_{0:N-1}}\sum_{i=0}^{N-1}\ell(\bar x_i,\bar u_i)+V_f(\bar x_N)$$ s.t. $\bar x_{i+1}=A\bar x_i+B\bar u_i$，$\bar x_i\in\mathcal{X}\ominus Z$，$\bar u_i\in\mathcal{U}\ominus KZ$，$\bar x_N\in\mathcal{X}_f\subseteq(\mathcal{X}\ominus Z)$，$x\in\{\bar x_0\}\oplus Z$（初始约束——该论文的关键创新：$\bar x_0$ 是决策变量而非 $=x$）。

定理 3.13.3（Rigid Tube MPC 鲁棒指数稳定）。若 $\mathcal{X}_f$ 是名义系统在局部律 $\bar u=K_f\bar x$ 下的正不变集、$V_f$ 是对应 Lyapunov 函数、$Z$ 是 mRPI 外近似，则闭环系统使 $\{0\}\oplus Z$ 鲁棒指数稳定，且 $\forall k\ge 0,\,x_k\in\mathcal{X},\,u_k\in\mathcal{U}$，对所有 $w\in\mathcal{W}^\infty$。

证明骨架：(i) 递归可行性：取候选 $(\bar u_1^*,\dots,\bar u_{N-1}^*,K_f\bar x_N^*)$ 与候选 $\bar x_0^+=\bar x_1^*$；验证 $x^+\in\{\bar x_1^*\}\oplus Z$（因 $e^+=A_Ke+w\in Z$）。(ii) 名义稳定性：$V_N^*(\bar x)$ 对名义系统是 Lyapunov 函数。(iii) 管道收缩 + 名义 ISS：联合给出 $x_k\to\{0\}\oplus Z$ 指数收敛。

工程含义：在线**只解一个 QP**（规模与名义 MPC 同阶）；$K,Z$ 离线一次算好。这是目前机器人部署最多的鲁棒 MPC 架构——ETH RSL 的 ANYmal、MIT Cheetah、多数四旋翼的 rate control 都用某种变体。

Rigid Tube MPC 完整设计流程（从零到部署）¶

以下是一个完整的 Rigid Tube MPC 设计实例，展示每一步的数学计算和工程决策。

Step 1：确定系统模型与不确定性集

考虑四旋翼高度通道的简化模型（二维双积分器）： $$x_{k+1} = \underbrace{\begin{bmatrix}1 & T_s\\ 0 & 1\end{bmatrix}}_{A} x_k + \underbrace{\begin{bmatrix}T_s^2/2\\ T_s\end{bmatrix}}_{B} u_k + w_k$$

其中 $x = [z, v_z]^\top$（高度和垂直速度），$u$ 是推力加速度偏差，$T_s = 0.05$ s。扰动 $w_k \in \mathcal{W} = \{w: \|w\|_\infty \le [0.005, 0.05]^\top\}$（风扰动对位置和速度的影响）。

Step 2：设计反馈增益 $K$ 并计算 $A_K$

选 LQR 权重 $Q = \text{diag}(100, 1)$，$R = 10$： $$K_{\text{lqr}} = -(R + B^\top P B)^{-1} B^\top P A = [-3.16, -1.77]$$ $$A_K = A + BK = \begin{bmatrix}0.996 & 0.046\\ -0.316 & 0.823\end{bmatrix}$$

验证：$\text{eig}(A_K) = \{0.76, 0.46\}$——均在单位圆内，$A_K$ Schur。

Step 3：计算 mRPI 外近似 $Z$

用 $(s, \alpha)$ 算法，取精度 $\varepsilon = 0.001$： - $s = 8$, $\alpha = 0.35$ - $Z \approx F(8, 0.35) = \frac{1}{0.65} \bigoplus_{i=0}^{7} A_K^i \mathcal{W}$

计算得 $Z$ 为一个 16 面多面体，在 $z$ 方向半径 $\approx 0.03$ m，$v_z$ 方向 $\approx 0.18$ m/s。

Step 4：收紧约束

原始约束：$\mathcal{X} = \{x: |z| \le 5, |v_z| \le 3\}$，$\mathcal{U} = \{u: |u| \le 5\}$。

收紧后： $$\mathcal{X}_{\text{tight}} = \mathcal{X} \ominus Z = \{x: |z| \le 4.97, |v_z| \le 2.82\}$$ $$\mathcal{U}_{\text{tight}} = \mathcal{U} \ominus KZ = \{u: |u| \le 4.68\}$$

Step 5：设计终端集和终端代价

终端集 $\mathcal{X}_f$：名义系统在 $K_f = K_{\text{lqr}}$ 下的最大正不变集，与 $\mathcal{X}_{\text{tight}}$ 取交集。终端代价 $V_f(\bar x) = \bar x^\top P_{\text{lqr}} \bar x$。

Step 6：在线 QP 的标准形式

决策变量：$\bar x_0, \bar u_0, \bar x_1, \bar u_1, \ldots, \bar x_N$（共 $(N+1) n_x + N n_u$ 个标量）。

\[\min \sum_{i=0}^{N-1} (\bar x_i^\top Q \bar x_i + \bar u_i^\top R \bar u_i) + \bar x_N^\top P \bar x_N\]

约束：动力学等式 + 收紧不等式 + 终端集 + 初始管道约束 $\|x_{\text{meas}} - \bar x_0\|_Z \le 1$。

Step 7：控制律执行

\[u_{\text{apply}} = \bar u_0^* + K(x_{\text{meas}} - \bar x_0^*)\]

本质洞察：Tube MPC 的控制律不是简单的 $u = \bar u_0^*$，而是**名义前馈 + 误差反馈**的叠加。名义前馈 $\bar u_0^*$ 负责"按计划行动"，误差反馈 $K(x - \bar x_0)$ 负责"把误差拉回管道内"。这与 RL 中 teacher-student 架构的思想异曲同工：teacher（名义 MPC）做规划，student（线性反馈 $K$）做跟踪。

Tube MPC 时变截面的 2024 前沿（Sieber et al. 2024, Wang et al. 2024 TAC）¶

2024 年的最新进展是**时变截面 Tube MPC**（Tube MPC with Time-Varying Cross-Sections）。核心思想：传统 Rigid Tube 用固定截面 $Z$，保守性来源于"用最坏情况 tube 约束所有时步"。时变截面允许 $Z_k$ 随预测时域变化：

\[x_k \in \{\bar x_k\} \oplus Z_k, \quad Z_{k+1} \supseteq A_K Z_k \oplus \mathcal{W}\]

Wang et al. (IEEE TAC, 2024) 证明了时变截面可以**严格缩小保守性**同时保持递归可行性。Sieber et al. (2024) 进一步将此与 data-driven 方法结合，用在线数据动态调整 $Z_k$。

3.13.3.3 Homothetic Tube MPC（Raković 2012）¶

出处：Raković, Kouvaritakis, Cannon, Panos, Findeisen, "Parameterized tube model predictive control", IEEE TAC 57(11):2746–2761, 2012, DOI:10.1109/TAC.2012.2191174。

构造：管道截面形状 $Z$ 固定，缩放因子 $\alpha_k\ge 0$ 随时间变化： $$x_k\in\{z_k\}\oplus\alpha_k Z,\quad u_k=v_k+K(x_k-z_k).$$

动力学递推利用 $A_K Z\oplus\mathcal{W}\subseteq\alpha^\star Z$（$\alpha^\star$ 为 mRPI 收缩率），给出标量不等式 $\alpha_{k+1}\ge\alpha^\star\alpha_k+\rho_\mathcal{W}$。收紧约束 $z_i\in\mathcal{X}\ominus\alpha_i Z,\,v_i\in\mathcal{U}\ominus\alpha_i KZ$。

价值：Rigid Tube 在 $k=0$ 就需要 $\bar x_0$ 留出 $Z$ 全尺寸的余量；Homothetic 允许 $\alpha_0$ 小，$\alpha_k$ 逐渐放大——可行域显著扩张。在线仍是 QP/LP，变量数仅多 $N$ 个标量。

3.13.3.4 Elastic Tube MPC（Raković–Levine–Açıkmeşe 2016, ACC）¶

出处：Raković, Levine, Açıkmeşe, "Elastic tube model predictive control", Proc. ACC 2016, pp. 3594–3599, DOI:10.1109/ACC.2016.7525471。

截面形状 $S(a)=\bigoplus_j a_j S_j$ 与局部反馈 $K(a)$ 均随 $a\in\mathbb{R}^p_{\ge 0}$ 在线可变。递推约束 $Az+Bv\oplus(A+BK(a))S(a)\oplus\mathcal{W}\subseteq z^+\oplus S(a^+)$ 退化为逐行线性不等式。

保守性最低，代价是在线变量从 $\mathcal{O}(N)$ 增至 $\mathcal{O}(N^2)$，仍为 LP/QP。

注：任务原文提到"Lorenzetti–Landry–Pavone 2020 elastic tube"实际对应的是 Lorenzetti & Pavone 的 tube-based 输出反馈 MPC（arXiv:1911.07360, ECC 2020），采用单一常数截面使在线复杂度与名义 MPC 同阶；而真正的 elastic tube MPC 原始文献是 Raković–Levine–Açıkmeşe 2016 ACC。两者学习时均应掌握。

3.13.3.5 Nonlinear Tube MPC（Mayne–Kerrigan–van Wyk–Falugi 2011）¶

出处：IJRNC 21(11):1341–1353, DOI:10.1002/rnc.1758。

思想：线性 $K$ 反馈在非线性系统上失效——改用 ancillary MPC（短视界辅助 NMPC）强迫 $x_k$ 跟 $\bar x_k$ 差在管道 $\mathcal{S}$ 内。在线需解 两个 NLP（名义 + 辅助），计算负担重。

现代替代：Köhler–Müller–Allgöwer 2020（见 §3.13.B）用 incremental Lyapunov / contraction 构造单层非线性 tube，在线只解一个 NMPC。

常见陷阱¶

陷阱	现象	正确做法
$K$ 选择过"软"	$A_K$ 特征值接近单位圆 → $Z$ 巨大 → 约束收紧后不可行	用 LQR 作初值，调权重矩阵 $Q,R$ 使 $Z$ 合理
低估 $\mathcal{W}$	仿真测试没出错，实车坏掉	$\mathcal{W}$ 应包含建模误差 + 延迟 + 量化误差的联合估计
忘记 $K_f$ 与 $K$ 的区别	终端集用错误反馈计算	终端局部律 $K_f$ 可独立于 tube 反馈 $K$，但 $\mathcal{X}_f\oplus Z\subseteq\mathcal{X}$ 必须满足
mRPI 近似精度不够	$(s,\alpha)$ 取得太粗 → 管道过大 → 保守	选 $\varepsilon\le 1\%\cdot\\|\mathcal{X}\\|$，$s$ 通常 10–20

Tube MPC 的 Python/CasADi 完整实现¶

以下给出一个可运行的 Rigid Tube MPC 完整实现，包含离线 mRPI 计算和在线 QP 求解。

import numpy as np
import casadi as ca
from scipy.linalg import solve_discrete_are

class RigidTubeMPC:
    """
    完整的 Rigid Tube MPC 实现

    离线计算:
        1. 设计反馈增益 K (LQR)
        2. 计算 mRPI 外近似 Z (Rakovic 算法)
        3. 收紧约束 X_tight = X ⊖ Z, U_tight = U ⊖ KZ
        4. 计算终端集 X_f 和终端代价 V_f

    在线求解:
        名义 QP (规模与普通 MPC 完全相同)
        + 反馈律 u = u_bar + K*(x - x_bar)
    """

    def __init__(self, A, B, Q, R, N, x_max, u_max, W_max):
        self.nx, self.nu = A.shape[0], B.shape[1]
        self.A, self.B = A, B
        self.Q, self.R = Q, R
        self.N = N

        # --- 离线 Step 1: LQR 增益 ---
        P_lqr = solve_discrete_are(A, B, Q, R)
        self.K = -np.linalg.solve(R + B.T @ P_lqr @ B, B.T @ P_lqr @ A)
        self.A_K = A + B @ self.K

        # 验证 A_K 稳定性
        eigs = np.linalg.eigvals(self.A_K)
        assert np.all(np.abs(eigs) < 1.0), f"A_K 不稳定! eig={eigs}"

        # --- 离线 Step 2: mRPI 外近似 (简化为 box 近似) ---
        self.Z = self._compute_mrpi_box(W_max, s_max=20, eps=1e-3)

        # --- 离线 Step 3: 约束收紧 ---
        self.x_tight = x_max - self.Z[:self.nx]
        self.u_tight = u_max - np.abs(self.K) @ self.Z[:self.nx]

        # --- 离线 Step 4: 终端代价 ---
        self.P_f = P_lqr

        # --- 在线 QP 构建 (CasADi) ---
        self._build_qp()

    def _compute_mrpi_box(self, W_max, s_max=20, eps=1e-3):
        """
        计算 mRPI 的 box 外近似 (适用于 box 扰动集)

        对 box W = [-W_max, W_max], 计算 Z = sum_{i=0}^{s-1} |A_K^i| * W_max / (1-alpha)
        """
        Z = np.zeros(self.nx)
        A_pow = np.eye(self.nx)

        for i in range(s_max):
            Z += np.abs(A_pow) @ W_max
            A_pow = A_pow @ self.A_K
            # 检查收敛
            alpha = np.max(np.abs(A_pow) @ W_max / (W_max + 1e-10))
            if alpha < 0.5:
                Z = Z / (1 - alpha)
                print(f"  mRPI 收敛: s={i+1}, alpha={alpha:.4f}")
                return Z

        print(f"  mRPI 警告: 未在 {s_max} 步内收敛")
        return Z * 2  # 保守放大

    def _build_qp(self):
        """构建 CasADi 的名义 MPC QP"""
        opti = ca.Opti('conic')

        # 决策变量: 名义状态和输入
        x_bar = opti.variable(self.nx, self.N + 1)
        u_bar = opti.variable(self.nu, self.N)

        # 参数: 当前测量状态
        x_meas = opti.parameter(self.nx)

        # 代价函数
        cost = 0
        for k in range(self.N):
            cost += ca.mtimes([x_bar[:, k].T, self.Q, x_bar[:, k]])
            cost += ca.mtimes([u_bar[:, k].T, self.R, u_bar[:, k]])
        cost += ca.mtimes([x_bar[:, self.N].T, self.P_f, x_bar[:, self.N]])
        opti.minimize(cost)

        # 约束
        for k in range(self.N):
            # 名义动力学
            opti.subject_to(x_bar[:, k+1] == 
                           self.A @ x_bar[:, k] + self.B @ u_bar[:, k])
            # 收紧的状态约束
            for i in range(self.nx):
                opti.subject_to(x_bar[i, k] <= self.x_tight[i])
                opti.subject_to(x_bar[i, k] >= -self.x_tight[i])
            # 收紧的输入约束
            for i in range(self.nu):
                opti.subject_to(u_bar[i, k] <= self.u_tight[i])
                opti.subject_to(u_bar[i, k] >= -self.u_tight[i])

        # 初始 tube 约束: x_meas ∈ {x_bar_0} ⊕ Z
        for i in range(self.nx):
            opti.subject_to(x_meas[i] - x_bar[i, 0] <= self.Z[i])
            opti.subject_to(x_meas[i] - x_bar[i, 0] >= -self.Z[i])

        # 存储
        self.opti = opti
        self.x_bar_var = x_bar
        self.u_bar_var = u_bar
        self.x_meas_par = x_meas

        # 求解器设置
        opts = {'printLevel': 'none'}
        opti.solver('osqp', opts)

    def solve(self, x_current):
        """
        在线求解: 给定当前状态, 返回控制输入

        u_apply = u_bar_0* + K * (x_current - x_bar_0*)
        """
        self.opti.set_value(self.x_meas_par, x_current)

        try:
            sol = self.opti.solve()
            x_bar_0 = sol.value(self.x_bar_var[:, 0])
            u_bar_0 = sol.value(self.u_bar_var[:, 0])

            # Tube 反馈律
            u_apply = u_bar_0 + self.K @ (x_current - x_bar_0)

            return u_apply, {'status': 'optimal', 'x_bar_0': x_bar_0}

        except RuntimeError:
            # 不可行 -> 紧急制动
            return np.zeros(self.nu), {'status': 'infeasible'}

使用示例：

# 四旋翼高度控制
A = np.array([[1, 0.05], [0, 1]])       # 双积分器, Ts=0.05
B = np.array([[0.00125], [0.05]])
Q = np.diag([100, 1])
R = np.array([[10]])
N = 20
x_max = np.array([5.0, 3.0])            # |z| <= 5m, |vz| <= 3m/s
u_max = np.array([5.0])                  # |a| <= 5m/s^2
W_max = np.array([0.005, 0.05])          # 风扰动

mpc = RigidTubeMPC(A, B, Q, R, N, x_max, u_max, W_max)
# 仿真
x = np.array([2.0, 0.5])  # 初始状态
for k in range(200):
    u, info = mpc.solve(x)
    w = np.random.uniform(-W_max, W_max)  # 随机扰动
    x = A @ x + B @ u + w

练习¶

练习 3.13.3.1 ⭐⭐：修改上述代码，把 Rigid Tube 替换为 Homothetic Tube（添加 $\alpha_k$ 标量决策变量和收缩约束 $\alpha_{k+1} \ge \alpha^* \alpha_k + \rho$）。比较两者的可行域大小。

练习 3.13.3.2 ⭐⭐⭐：对上述四旋翼模型，用 MPT3 (MATLAB) 或 pypoman (Python) 精确计算多面体形式的 mRPI，与 box 近似对比体积差异。

练习 3.13.3.3 ⭐⭐⭐（跨章综合）：回顾 3.7 节 ISS 理论。证明 Rigid Tube MPC 闭环系统是 ISS 稳定的：构造 Lyapunov 函数 $V(x) = V_N^*(\bar x^*(x)) + \|x - \bar x^*(x)\|_P^2$（其中 $P$ 满足 $A_K^\top P A_K - P \preceq -\gamma I$），推导 ISS 增益 $\gamma_{\text{ISS}}(\|w\|)$。

桥梁¶

→ 3.7 ISS：tube 收缩性 = $A_K$ 的 ISS-Lyapunov 条件。
→ 3.11.6：本节是 3.11 Tube 章节的**完整深化**。
→ §3.13.8：把 mRPI 用 GP 后验协方差构造，得到 learning-based tube。

§3.13.4 LPV Tube MPC：把非线性嵌入参数变化¶

严格定义¶

LPV（Linear Parameter-Varying）系统： $$x_{k+1}=A(\theta_k)x_k+B(\theta_k)u_k,\qquad \theta_k\in\Theta=\mathrm{conv}\{\theta^1,\dots,\theta^q\}.$$

$A(\theta),B(\theta)$ 对 $\theta$ 仿射。LPV embedding 是把非线性系统 $x_{k+1}=f(x_k,u_k)$ 写成 $x_{k+1}=A(\rho(x_k,u_k))x_k+B(\rho)u_k$，其中 $\rho$ 是内生 scheduling 变量（如速度、姿态角）。

核心方法（Hanema–Tóth–Lazar 2016 CDC, DOI:10.1109/CDC.2016.7798472；Automatica 2017 85:137–144）¶

管道 + 策略序列：$X_{i|k}=z_{i|k}\oplus\alpha_{i|k}\Omega$（homothetic），$\Omega$ 为 periodic contractive 终端集；顶点控制 $\Pi_{i|k}(x,\theta)=\sum_j\lambda_j(\theta)[K_j(x-z_{i|k})+v_{i|k}]$，$\lambda_j$ 为 $\theta$ 的重心坐标。

Anticipation：未来 $\theta$ 的预告集 $\Theta_{i|k}\supseteq\{\theta_{k+i}\}$，可利用位置/速度等已知状态收缩。

在线问题：单个 LP，$\mathcal{O}(N)$ 变量；递归可行 + 渐近稳定。

工程含义¶

LPV-Tube 是**机器人非线性鲁棒 MPC 的实用近道**：四旋翼把倾角 $\phi$ 作 $\theta$、四足把 CoM 高度作 $\theta$、机械臂把关节角作 $\theta$——避免直接解非线性 tube 的双层 NMPC。保守性取决于 $\Theta$ 包络紧致度。

LPV 嵌入的完整构造方法¶

为什么不直接用非线性 Tube MPC？ 非线性 Tube MPC（Mayne 2011）需要解两个 NLP——名义 NMPC + 辅助 NMPC，计算量是名义 NMPC 的 2-3 倍。对于实时性要求高的机器人（如 200 Hz 四旋翼），这不可接受。LPV 嵌入把非线性系统"包裹"在一个参数变化的线性系统族中，使得 Tube MPC 的所有线性工具（mRPI、Pontryagin 差、QP）都直接可用——代价是引入了包络保守性。

四旋翼 LPV 嵌入示例：

四旋翼的非线性动力学（高度-俯仰平面简化）： $$\dot v_x = -\frac{T}{m}\sin\phi, \quad \dot v_z = \frac{T}{m}\cos\phi - g$$

令 scheduling 变量 $\theta = \phi$（俯仰角），定义 $\theta \in [-\phi_{\max}, \phi_{\max}]$（如 $\pm 30°$）。离散化后： $$x_{k+1} = \underbrace{\begin{bmatrix}1 & T_s & 0 & 0\\ 0 & 1 & -\frac{T_0 T_s}{m}\sin\theta_k & 0\\ 0 & 0 & 1 & T_s\\ 0 & 0 & \frac{T_0 T_s}{m}\cos\theta_k & 1\end{bmatrix}}_{A(\theta_k)} x_k + B u_k$$

由于 $\sin\theta$ 和 $\cos\theta$ 都可以被 $\theta \in [-30°, 30°]$ 的凸包覆盖（$\sin\theta \in [-0.5, 0.5]$, $\cos\theta \in [0.87, 1]$），$A(\theta)$ 属于**两个顶点矩阵的凸包** $\text{conv}\{A(\theta^1), A(\theta^2)\}$（取 $\theta^1 = -30°$, $\theta^2 = +30°$）。

然后对每个顶点设计独立的反馈 $K_1, K_2$，用重心坐标插值： $$K(\theta) = \lambda_1(\theta) K_1 + \lambda_2(\theta) K_2, \quad \lambda_1 + \lambda_2 = 1, \lambda_i \ge 0$$

反事实推理：如果把整个 $[-30°, 30°]$ 当作单一不确定性集（而非 LPV），则需要一个对所有 $\theta$ 鲁棒的单一 $K$——这等价于把非线性项完全当作加法扰动，保守性极大。LPV 的价值在于**利用了 $\theta$ 是可测量的这一信息**：每步可以根据当前 $\theta_k$ 选择最优的反馈 $K(\theta_k)$。

常见陷阱¶

陷阱	现象	正确做法
$\Theta$ 取太大	LPV 退化为鲁棒 MPC	利用状态信息收缩 $\Theta_{i\\|k}$
忽略 $\theta$ 的变化率	$\theta_{k+1}$ 可能跳出 $\Theta$	加约束 $\\|\theta_{k+1}-\theta_k\\| \le \delta$
顶点太多	顶点 LMI 爆炸	用椭球包络替代多面体

练习¶

练习 3.13.4.1 ⭐⭐：对上述四旋翼 LPV 模型，写出完整的 LPV-Tube MPC 在线 LP。标明每个约束对应的物理含义。

练习 3.13.4.2 ⭐⭐⭐：比较 LPV-Tube MPC 与 SLQ-NMPC（§3.11/3.12 中的方法）在四旋翼大角度机动时的性能。什么时候 LPV 近似开始失效？

§3.13.5 随机 MPC 概述：三大范式鸟瞰¶

严格定义¶

机会约束 (chance constraint)：对每步每条约束 $g_j$， $$\Pr[g_j(x_k,u_k,w_k)\le 0]\ge 1-\varepsilon_j.$$

随机 OCP： $$\min_\pi\mathbb{E}\!\left[\sum_{k=0}^{N-1}\ell(x_k,u_k)+V_f(x_N)\right],\quad \mathrm{s.t.}\ x_{k+1}=f(x_k,u_k,w_k),\,u_k=\pi_k(\cdot),\,\Pr[x_k\in\mathcal{X}]\ge 1-\varepsilon.$$

三大范式¶

1. 解析重构。假设 $w\sim\mathcal{N}(\mu,\Sigma)$ 与约束线性 $g=a^\top x+b$，则 $$\Pr[a^\top x+b\le 0]\ge 1-\varepsilon\iff a^\top\mu+b+\Phi^{-1}(1-\varepsilon)\|\Sigma^{1/2}a\|_2\le 0$$ （SOCP）。非高斯用 Chebyshev–Cantelli： $$a^\top\mu+b+\sqrt{\tfrac{1-\varepsilon}{\varepsilon}}\|\Sigma^{1/2}a\|_2\le 0$$ （分布无关，显著更保守）。

2. Scenario-based。用 $N_s$ 个独立样本替换概率约束（§3.13.6）。

3. Distributionally Robust。$\inf_{\mathbb{Q}\in\mathcal{P}_{\mathrm{amb}}}\Pr_\mathbb{Q}[g\le 0]\ge 1-\varepsilon$（§3.13.7）。

核心对照（Mesbah 2016, IEEE CSM 36(6):30–44, DOI:10.1109/MCS.2016.2602087）¶

范式	前提	保守性	在线复杂度	工程场景
解析（高斯）	$w$ 高斯，$g$ 仿射	紧	SOCP	过程控制、HVAC
解析（Cantelli）	知 $\mu,\Sigma$	保守	SOCP	不知分布但知矩
Scenario	可采样	$(\varepsilon,\beta)$ 可调	QP（规模 $\propto N_s$）	建筑、风机、自动驾驶
DRO	知模糊集	抗分布误判	SOCP/SDP	数据驱动、sim-to-real

Open-loop vs Closed-loop 方差传播¶

Open-loop：$u_k=\bar u_k\Rightarrow\Sigma_{k+1}=A\Sigma_kA^\top+\Sigma_w$（单调增大、保守）。

Closed-loop（仿射反馈 $u_k=\bar u_k+K_k(x_k-\bar x_k)$）：$\Sigma_{k+1}=(A+BK_k)\Sigma_k(A+BK_k)^\top+\Sigma_w$；若 $A+BK_k$ 稳定，$\Sigma_k$ 有界。这是 tube stochastic MPC 的标准技巧。

与 Robust MPC 的关系¶

取 $\varepsilon=0$ 且 $w\in\mathcal{W}$ 有界支撑 ⇒ 退化为 Robust MPC。$\varepsilon>0$ 允许稀有违反，换取可行域扩张与性能提升——对无界高斯扰动这是唯一可行路径。

解析 Chance Constraint 的完整推导¶

问题：约束 $a^\top x + b \le 0$ 需以概率 $\ge 1-\varepsilon$ 满足，其中 $x$ 受随机扰动 $x \sim \mathcal{N}(\mu, \Sigma)$。

推导： $$\Pr[a^\top x + b \le 0] \ge 1 - \varepsilon$$ $$\iff \Pr[a^\top x \le -b] \ge 1 - \varepsilon$$

因为 $a^\top x \sim \mathcal{N}(a^\top \mu, a^\top \Sigma a)$，标准化： $$\Pr\left[\frac{a^\top x - a^\top \mu}{\sqrt{a^\top \Sigma a}} \le \frac{-b - a^\top \mu}{\sqrt{a^\top \Sigma a}}\right] \ge 1-\varepsilon$$

\[\iff \Phi\left(\frac{-b - a^\top \mu}{\sqrt{a^\top \Sigma a}}\right) \ge 1-\varepsilon\]

\[\iff \frac{-b - a^\top \mu}{\sqrt{a^\top \Sigma a}} \ge \Phi^{-1}(1-\varepsilon)\]

\[\iff a^\top \mu + b + \Phi^{-1}(1-\varepsilon)\sqrt{a^\top \Sigma a} \le 0\]

\[\iff a^\top \mu + b + \Phi^{-1}(1-\varepsilon)\|\Sigma^{1/2} a\|_2 \le 0\]

这是一个**二阶锥约束** (SOCP)，可以高效求解。

物理含义：均值必须离约束面至少 $\Phi^{-1}(1-\varepsilon)$ 个标准差的距离。当 $\varepsilon = 0.05$ 时，$\Phi^{-1}(0.95) = 1.645$（约 1.6 个标准差）；$\varepsilon = 0.01$ 时，$\Phi^{-1}(0.99) = 2.326$。

与 Tube MPC 的精确对比：

维度	Tube MPC	Stochastic MPC (Gaussian)
不确定性模型	$w \in \mathcal{W}$（有界集）	$w \sim \mathcal{N}(0, \Sigma_w)$（无界）
约束处理	$x \in \mathcal{X} \ominus Z$（确定性收紧）	$\Pr[x \in \mathcal{X}] \ge 1-\varepsilon$（概率收紧）
保守性	对 $\mathcal{W}$ 内所有扰动鲁棒（100%）	允许 $\varepsilon$ 违反率（< 100%）
可行域	较小（为最坏情况预留）	较大（只为典型情况预留）
安全保证	绝对（$\forall w \in \mathcal{W}$）	概率（$1-\varepsilon$ 概率）

本质洞察：Tube MPC 和 Stochastic MPC 不是两种不同的方法，而是**同一个 trade-off 的两个极端**。Tube 对应 $\varepsilon = 0$（绝对安全），Stochastic 对应 $0 < \varepsilon < 1$（概率安全）。工程师的选择取决于**违反后果的严重程度**：飞机防撞用 Tube（$\varepsilon = 0$），HVAC 舒适度控制用 Stochastic（$\varepsilon = 5\%$ 可接受）。

Tube Stochastic MPC：两者的统一¶

Cannon-Kouvaritakis-Raković (2010, Automatica) 提出了 Stochastic Tube MPC，将两种方法统一：

用线性反馈 $u_k = \bar u_k + K(x_k - \bar x_k)$ 把闭环方差限制在 $\Sigma_k \preceq \bar\Sigma$
在名义问题中用概率约束 $a^\top \bar x_k + \Phi^{-1}(1-\varepsilon)\|(\bar\Sigma)^{1/2} a\|_2 \le g_i$
反馈 $K$ 的设计使 $\Sigma_{k+1} = (A+BK)\Sigma_k(A+BK)^\top + \Sigma_w$ 有界

这是当前机器人高斯扰动场景（如 SLAM 不确定性、传感器噪声）下的标准做法。

常见陷阱¶

陷阱	现象	正确做法
高斯假设用于非高斯扰动	实际违反率远超 $\varepsilon$	用 Chebyshev-Cantelli 或 scenario
忽略方差传播的闭环效应	开环方差单调增大导致过保守	用仿射反馈 $u = \bar u + K(x-\bar x)$
$\varepsilon$ 设得太小	SOCP 约束过紧→不可行	对安全关键约束 $\varepsilon=10^{-3}$，非关键 $\varepsilon=0.05$
多步约束用联合概率	单步 $1-\varepsilon$ 不保证全程 $1-\varepsilon$	用 Bonferroni $\varepsilon_k = \varepsilon/N$ 或 conditional 方法

练习¶

练习 3.13.5.1 ⭐⭐：对二维系统 $x_{k+1} = Ax_k + Bu_k + w_k$，$w \sim \mathcal{N}(0, \Sigma_w)$，$\mathcal{X} = \{x: |x_i| \le 5\}$，$\varepsilon = 0.05$。推导每步的 SOCP 约束具体形式。如果 $\Sigma_w = 0.1 I_2$，计算约束收紧量并与 Tube MPC（$\mathcal{W} = [-0.3, 0.3]^2$，覆盖 3-sigma）对比可行域大小。

练习 3.13.5.2 ⭐⭐⭐：证明 Chebyshev-Cantelli 不等式 $\Pr[X \ge \mu + t] \le \frac{\sigma^2}{\sigma^2 + t^2}$ 并用它推导分布无关的 chance constraint 重写。与高斯假设对比保守性。

§3.13.6 Scenario MPC：把概率约束变成 $N_s$ 条硬约束¶

核心定理（Calafiore–Campi 2006，IEEE TAC 51(5):742–753，DOI:10.1109/TAC.2006.875041）¶

定理 3.13.4（Scenario 样本复杂度）。考虑凸机会约束规划 $$\min_{x\in\mathbb{R}^d}c^\top x\quad\mathrm{s.t.}\ \Pr_\delta\{f(x,\delta)\le 0\}\ge 1-\varepsilon,$$ 其中 $f(\cdot,\delta)$ 对每个 $\delta$ 凸。取 $N$ 个独立样本 $\{\delta^{(i)}\}_{i=1}^N$ 求解 $$\min c^\top x\quad\mathrm{s.t.}\ f(x,\delta^{(i)})\le 0,\ i=1,\dots,N,$$ 若 $$N\ge\left\lceil\frac{2}{\varepsilon}\!\left(\ln\frac{1}{\beta}+d\right)\right\rceil,$$ 则 $\Pr^N\{V(x_N^*)>\varepsilon\}\le\beta$，其中 $V(x)=\Pr\{f(x,\delta)>0\}$。

Campi–Garatti 2008 精确界（SIAM J. Optim. 19(3):1211–1230，DOI:10.1137/07069821X）： $$\Pr^N\{V(x_N^*)>\varepsilon\}\le\sum_{i=0}^{d-1}\binom{N}{i}\varepsilon^i(1-\varepsilon)^{N-i}.$$ Fully-supported 问题取等号。实用中用该精确 Beta 尾分布数值反解 $N$，样本数可比 $\lceil(2/\varepsilon)(\ln(1/\beta)+d)\rceil$ 少一半。

证明骨架：Helly 定理给出 $d$ 维凸问题的 support constraint 数 $\le d$ ⇒ $x_N^*$ 仅由最多 $d$ 个样本决定 ⇒ 组合论证 $\Pr\{V>\varepsilon\}\le\binom{N}{d}(1-\varepsilon)^{N-d}$ ⇒ 反解 $N$。

Scenario MPC（Schildbach–Fagiano–Frei–Morari 2014, Automatica 50(12):3009–3018, DOI:10.1016/j.automatica.2014.10.035）¶

Support rank 技巧：对每条约束 $j$ 的支持秩 $\rho_j\le d$（仅约束某子空间维度），样本数 $$\sum_{i=0}^{\rho_j-1}\binom{N_j}{i}\varepsilon_j^i(1-\varepsilon_j)^{N_j-i}\le\beta_j.$$

闭环 "average-in-time" 解释：把机会约束理解为闭环时间平均违反率 $\le\varepsilon$，在 receding horizon 下用 martingale/ergodic 论证成立。

示例计算：$\varepsilon=0.05$、$\beta=10^{-4}$、$d=20$（典型 $N=10$ MPC with $n_u=2$）： - Calafiore–Campi 界：$N\ge(2/0.05)(\ln 10^4+20)=40\cdot(9.21+20)\approx 1170$。 - Campi–Garatti 精确：数值解约 $N\approx 580$。 - Support rank $\rho=1$（单一约束）：$N\approx 190$。

与蒙特卡洛的本质区别¶

维度	蒙特卡洛	Scenario approach
目的	验证：估计 $V(x)=\Pr\{f>0\}$	合成：求满足概率约束的 $x^*$
使用样本方式	后验评估 $\hat V=\frac{1}{N}\sum\mathbf{1}\{f>0\}$	作为硬约束参与优化
收敛率	$\mathcal{O}(1/\sqrt{N})$	$\mathcal{O}(d/\varepsilon\cdot\log(1/\beta))$
概率保证	大数律 / Hoeffding	VC / support constraint / compression
典型 $N$（$\varepsilon=5\%$）	$10^4$（估尾概率）	$10^2$–$10^3$（合成决策）

工程含义¶

Scenario MPC 是**分布无关（distribution-free）、仅需采样能力的方法——**对机器人 sim-to-real 非常友好：从 domain randomization 的仿真器直接抽 $N_s$ 条轨迹作为 scenario 即可。缺点：在线 QP 规模 $\propto N_s$，对 $\varepsilon<1\%$ 不实用。

Scenario MPC 的完整推导与直觉¶

为什么 Scenario 方法有效？ 这背后是一个深刻的组合-概率论结果。让我们从直觉出发：

类比：想象你要设计一座桥，需要承受"所有可能的车辆荷载"。你无法列举所有可能荷载，但你可以随机采样 $N$ 种典型荷载，设计一座能承受这 $N$ 种荷载的桥。问题是：你需要多少样本 $N$，才能保证这座桥对 $1-\varepsilon$ 比例的真实荷载也安全？

Calafiore-Campi 的天才发现是：对**凸**问题，答案只取决于**决策变量维度 $d$** 和**违反率 $\varepsilon$，而**与扰动分布无关！

为什么凸性是关键？ 凸问题的最优解只由最多 $d$ 个"支撑约束"决定（Helly 定理的推论）。删除任何一个非支撑样本不改变最优解。这意味着解的"复杂度"由 $d$ 而非 $N$ 控制——这就是为什么样本需求 $N = O(d/\varepsilon)$ 而非 $O(1/\varepsilon^2)$（后者是非凸情形的典型需求）。

如果问题非凸会怎样？ 样本复杂度可能指数增长（$N = O(2^d/\varepsilon)$），且没有全局保证。这是 Scenario 方法的根本局限——也是为什么 NMPC 中直接用 scenario 需要额外的凸近似或分段线性化。

Scenario MPC 的 Python 实现骨架¶

import numpy as np
from scipy.optimize import linprog
import cvxpy as cp

def scenario_mpc(x0, A, B, N_horizon, N_scenarios, epsilon, beta,
                 Q, R, x_max, u_max, sampler):
    """
    Scenario MPC 求解器

    Args:
        x0: 当前状态
        A, B: 系统矩阵
        N_horizon: 预测时域
        N_scenarios: scenario 数量 (由 Campi-Garatti 界确定)
        epsilon: 允许违反率
        beta: 置信参数
        Q, R: 代价权重
        x_max, u_max: 约束
        sampler: 扰动采样器 (从仿真器/domain randomization 获取)

    Returns:
        u0_opt: 第一步最优控制
    """
    nx, nu = A.shape[0], B.shape[1]
    d = nu * N_horizon  # 决策变量维度

    # Campi-Garatti 精确界: 计算所需最小样本数
    N_min = compute_campi_garatti_bound(d, epsilon, beta)
    N_s = max(N_scenarios, N_min)

    # 采样 scenarios (从仿真器 domain randomization)
    W_scenarios = [sampler(N_horizon) for _ in range(N_s)]

    # 构建 scenario QP
    u_var = cp.Variable((nu, N_horizon))
    x_var = cp.Variable((nx, N_horizon + 1))

    cost = 0
    constraints = [x_var[:, 0] == x0]

    # 名义动力学 + 代价
    for k in range(N_horizon):
        cost += cp.quad_form(x_var[:, k], Q) + cp.quad_form(u_var[:, k], R)
        constraints += [cp.norm_inf(u_var[:, k]) <= u_max]
    cost += cp.quad_form(x_var[:, N_horizon], Q)

    # 对每个 scenario 添加硬约束
    for i in range(N_s):
        x_scen = cp.Variable((nx, N_horizon + 1))
        constraints += [x_scen[:, 0] == x0]
        for k in range(N_horizon):
            w_k = W_scenarios[i][k]  # 第 i 个 scenario 在时刻 k 的扰动
            constraints += [
                x_scen[:, k+1] == A @ x_scen[:, k] + B @ u_var[:, k] + w_k,
                cp.norm_inf(x_scen[:, k+1]) <= x_max  # 状态约束
            ]

    prob = cp.Problem(cp.Minimize(cost), constraints)
    prob.solve(solver=cp.OSQP, eps_abs=1e-6)

    return u_var[:, 0].value

def compute_campi_garatti_bound(d, epsilon, beta):
    """
    Campi-Garatti 2008 精确界的数值反解

    找最小 N 使得 sum_{i=0}^{d-1} C(N,i) * eps^i * (1-eps)^{N-i} <= beta
    """
    from scipy.special import comb
    for N in range(d, 10000):
        prob_sum = sum(comb(N, i) * epsilon**i * (1-epsilon)**(N-i) 
                      for i in range(d))
        if prob_sum <= beta:
            return N
    return 10000  # 退化为保守界

Scenario MPC 与 Domain Randomization 的精确对应¶

Domain Randomization (RL)	Scenario MPC
物理参数 $\xi \sim \mathcal{U}(\Xi)$	扰动样本 $w^{(i)} \sim \mathbb{P}$
训练环境数量	Scenario 数量 $N_s$
策略对所有环境鲁棒	控制序列对所有 scenario 可行
无理论保证（经验）	$(1-\varepsilon, 1-\beta)$ 概率保证
训练时间长	在线 QP 规模大

关键洞察：Domain randomization 本质上就是 scenario approach 在 RL 中的实现。如果你用 $N$ 个随机环境训练 RL 策略，Campi-Garatti 理论可以给出该策略在新环境中失败概率的上界——前提是策略参数化使问题凸（线性策略或 convex parameterization）。

§3.13.7 分布鲁棒 MPC：Wasserstein 球内最坏情形¶

严格定义¶

Wasserstein-$p$ 距离： $$W_p(\mathbb{P},\mathbb{Q})=\left(\inf_{\gamma\in\Pi(\mathbb{P},\mathbb{Q})}\int\|x-y\|^p\,\mathrm{d}\gamma(x,y)\right)^{1/p}.$$

经验中心模糊集：$\mathbb{B}_\varepsilon(\hat{\mathbb{P}}_N)=\{\mathbb{Q}:W_p(\mathbb{Q},\hat{\mathbb{P}}_N)\le\varepsilon\}$，$\hat{\mathbb{P}}_N=\frac{1}{N}\sum_i\delta_{\hat\xi_i}$。

DRO： $$\hat J_N(\varepsilon)=\inf_x\sup_{\mathbb{Q}\in\mathbb{B}_\varepsilon(\hat{\mathbb{P}}_N)}\mathbb{E}^\mathbb{Q}[h(x,\xi)].$$

核心定理（Mohajerin Esfahani–Kuhn 2018, Math. Program. 171:115–166, DOI:10.1007/s10107-017-1172-1）¶

定理 3.13.5（Wasserstein DRO Tractable Reformulation）。对 Lipschitz 凸 $h(x,\cdot)$、$p=1$： $$\sup_{\mathbb{Q}\in\mathbb{B}_\varepsilon}\mathbb{E}^\mathbb{Q}[h(x,\xi)]=\inf_{\lambda\ge 0,s_i,z_{ik}}\left\{\lambda\varepsilon+\frac{1}{N}\sum_{i=1}^N s_i:[-h]^*(z_{ik}-\hat\xi_i)+\sigma_\Xi(z_{ik})-\langle z_{ik},\hat\xi_i\rangle\le s_i,\|z_{ik}\|_*\le\lambda\right\}.$$

Finite-sample 保证：若 $\mathbb{P}$ light-tail，则 $\Pr\{W_1(\hat{\mathbb{P}}_N,\mathbb{P})>\varepsilon\}\le c_1e^{-c_2N\varepsilon^{\max\{m,2\}}}$；取 $\varepsilon_N(\beta)=(\ln(c_1/\beta)/c_2N)^{1/\max\{m,2\}}$ 给出 $\Pr\{\mathbb{E}^\mathbb{P}[h(\hat x_N,\xi)]\le\hat J_N(\varepsilon_N)\}\ge 1-\beta$。

证明骨架：(i) Kantorovich 对偶把 $\sup_\mathbb{Q}$ 转为 $\inf_\lambda\{\lambda\varepsilon+\mathbb{E}^{\hat{\mathbb{P}}_N}[\sup_\xi(h(x,\xi)-\lambda\|\xi-\hat\xi_i\|)]\}$；(ii) 对凸 $h$ 的内 $\sup$ 用共轭函数 $[-h]^*$ 展开得有限凸约束；(iii) 经验测度集中（Fournier–Guillin）给出 $\varepsilon_N$ 样本保证。

DR-DeePC（Coulson–Lygeros–Dörfler 2022, IEEE TAC 67(7):3289–3304, DOI:10.1109/TAC.2021.3097706）¶

对 Hankel 数据矩阵噪声构造 Wasserstein DRO。核心结果：L2 正则化的 DeePC 等价于 Wasserstein DR-DeePC——正则半径 = Wasserstein 半径 $\varepsilon_N$。即 $$\lambda_g\|g\|_2^2\iff\text{对数据矩阵做椭球 robust LS}.$$

Robust–DR–Scenario 的连续谱¶

\[\underbrace{\text{Robust MPC}}_{\varepsilon\to\infty,\,\mathcal{P}=\text{所有支撑于}\mathcal{W}}\ \longleftrightarrow\ \underbrace{\text{DR-MPC}}_{0<\varepsilon_N<\infty}\ \longleftrightarrow\ \underbrace{\text{Scenario MPC}}_{\varepsilon\to 0,\,\hat{\mathbb{P}}_N}.\]

$\varepsilon\to\infty$：仅保留支撑集信息 → min-max。
中等 $\varepsilon_N$：out-of-sample $(1-\beta)$-保证，保守性可调。
$\varepsilon\to 0$：退化为样本平均近似（SAA）/ scenario。

KL vs Wasserstein 对比¶

维度	KL 球 $\mathrm{KL}(\mathbb{Q}\\|\hat{\mathbb{P}}_N)\le\eta$	Wasserstein 球 $W_p\le\varepsilon$
支撑	$\mathbb{Q}\ll\hat{\mathbb{P}}_N$（只在样本点加权）	样本邻域任意位置，支持连续
连续真分布	KL 到经验 = $+\infty$，保证破坏	$W_1(\hat{\mathbb{P}}_N,\mathbb{P})=O(N^{-1/m})$
重构	Exponential cone	SOCP / LP
控制应用	风险敏感、相对熵约束	DeePC、tube DR、sim-to-real

机器人首选 Wasserstein——几何直觉（扰动相近分布应相近）、连续噪声下非平凡保证、导出正则化解释。

Wasserstein DRO 的完整推导链¶

问题：给定 $N$ 个扰动样本 $\{\hat\xi_1, \ldots, \hat\xi_N\}$，如何求解 $$\inf_x \sup_{\mathbb{Q} \in \mathbb{B}_\varepsilon(\hat{\mathbb{P}}_N)} \mathbb{E}^{\mathbb{Q}}[h(x, \xi)]$$

为什么不直接用样本均值？ 如果只用 $\hat{\mathbb{P}}_N = \frac{1}{N}\sum_i \delta_{\hat\xi_i}$，则 $\inf_x \mathbb{E}^{\hat{\mathbb{P}}_N}[h] = \inf_x \frac{1}{N}\sum_i h(x, \hat\xi_i)$——这就是 Sample Average Approximation (SAA)。SAA 的问题是**过拟合**：它对已观察到的样本最优，但对新样本（out-of-sample）可能很差。DRO 通过在分布空间中添加"余量" $\varepsilon$ 来对抗过拟合。

推导 Step 1：Kantorovich 对偶

对 $W_1$ (Wasserstein-1) 距离，Kantorovich 对偶定理给出： $$\sup_{\mathbb{Q}: W_1(\mathbb{Q}, \hat{\mathbb{P}}_N) \le \varepsilon} \mathbb{E}^{\mathbb{Q}}[h(x, \xi)] = \inf_{\lambda \ge 0} \left\{\lambda \varepsilon + \frac{1}{N}\sum_{i=1}^N \sup_{\xi} [h(x, \xi) - \lambda\|\xi - \hat\xi_i\|]\right\}$$

直觉：$\lambda$ 是 Wasserstein 球约束的 Lagrange 乘子。内层 $\sup_\xi$ 找的是"在距样本 $\hat\xi_i$ 不太远的地方，使损失最大的扰动"——$\lambda$ 惩罚"远离样本"的代价。

推导 Step 2：利用凸共轭简化内层

当 $h(x, \cdot)$ 关于 $\xi$ 凸时，$\sup_\xi[h(x,\xi) - \lambda\|\xi - \hat\xi_i\|]$ 可以通过共轭函数展开为有限维约束： $$= \begin{cases} h(x, \hat\xi_i) & \text{if } \text{Lip}(h) \le \lambda \\ +\infty & \text{otherwise}\end{cases}$$

对 Lipschitz 代价 $h$，当 $\lambda \ge \text{Lip}(h)$ 时内层有限，否则无穷。因此： $$\sup_{\mathbb{Q}} \mathbb{E}^{\mathbb{Q}}[h] = \frac{1}{N}\sum_i h(x, \hat\xi_i) + \varepsilon \cdot \text{Lip}(h)$$

这就是 DRO 的核心结果：最坏情况期望 = 样本均值 + $\varepsilon$ 乘以 Lipschitz 常数。

推导 Step 3：选择半径 $\varepsilon_N$

Fournier-Guillin (2015) 集中不等式给出：对 $m$ 维轻尾分布 $\mathbb{P}$， $$\Pr\{W_1(\hat{\mathbb{P}}_N, \mathbb{P}) > \varepsilon\} \le c_1 e^{-c_2 N \varepsilon^{\max\{m, 2\}}}$$

反解：$\varepsilon_N(\beta) = \left(\frac{\ln(c_1/\beta)}{c_2 N}\right)^{1/\max\{m,2\}}$

取此 $\varepsilon_N$，则以概率 $\ge 1-\beta$，真实分布 $\mathbb{P} \in \mathbb{B}_{\varepsilon_N}(\hat{\mathbb{P}}_N)$，因此 DRO 解具有 out-of-sample 性能保证。

关键性质：$\varepsilon_N \to 0$ 当 $N \to \infty$（数据越多越不保守）；$\varepsilon_N$ 随维度 $m$ 增大而增大（维度灾难）。

双重解读： - 角度 1（鲁棒优化）：DRO 在分布空间中寻找最坏分布，是一种**分布鲁棒性**保证。 - 角度 2（正则化）：DRO 等价于在 SAA 目标上加 Lipschitz 正则项——是一种**结构化正则化**，防止过拟合。

这两个角度的统一是 DRO 最深刻的理论贡献：鲁棒性 = 正则化。

从 DRO 到 DeePC 的正则化等价¶

Coulson-Lygeros-Dörfler (2022 TAC) 的核心定理可以这样理解：

DeePC 的正则化项 $\lambda_g\|g\|_2^2$ 不是一个"调参技巧"——它有精确的鲁棒性含义。当数据矩阵 $H$ 受到噪声扰动 $H + \Delta H$ 时， $$\min_g \|Hg - b\|^2 + \lambda\|g\|^2 = \min_g \max_{\|\Delta H\| \le \rho} \|(H+\Delta H)g - b\|^2$$

其中 $\rho = \sqrt{\lambda}$。这就是 Tikhonov 正则化的鲁棒最小二乘解释。

将此与 Wasserstein DRO 连接：如果把数据矩阵的每一列看作一个"样本"，则列扰动的范数界 $\rho$ 等价于 Wasserstein 半径 $\varepsilon$。因此 DeePC 的 $\lambda_g$ 直接对应 DRO 的半径参数——$\lambda_g$ 越大，对数据噪声越鲁棒，但也越保守。

工程选取 $\lambda_g$ 的指导原则： - 数据质量好（受控实验、低噪声）：$\lambda_g = 10^{-4}$–$10^{-2}$ - 数据质量中等（开环采集、中等噪声）：$\lambda_g = 10^{-2}$–$1$ - 数据质量差（闭环数据、高噪声）：$\lambda_g = 1$–$100$

§3.13.8 自适应 / 学习型鲁棒 MPC¶

框架综述（Hewing–Wabersich–Menner–Zeilinger 2020, ARCRAS 3:269–296, DOI:10.1146/annurev-control-090419-075625）¶

三条主线： 1. 学习预测模型：$\hat f(x,u)+\hat d(x,u)$，$\hat d$ 由 GP/NN 学残差动力学。 2. 学习 MPC 参数：代价、约束、终端集的可微参数化（Amos 2018 differentiable MPC）。 3. 学习安全集：用数据扩展 RPI / 终端集（Wabersich 2021 PSF）。

GP-MPC（Koller–Berkenkamp–Turchetta–Krause 2018, CDC, arXiv:1803.08287）¶

构造：$x_{k+1}=f(x_k,u_k)+g(x_k,u_k)$，$g\sim\mathcal{GP}(0,k(\cdot,\cdot))$，后验 $\mathcal{N}(\mu_g(x,u),\sigma_g^2(x,u))$。用 Chowdhury–Gopalan 界 $\beta_\delta$（RKHS 正则性假设），高概率 tube： $$\Pr\{\|g(x,u)-\mu_g(x,u)\|\le\beta_\delta\sigma_g(x,u)\}\ge 1-\delta.$$

收紧约束 $h(x)+\beta_\delta\sigma(x)\le 0$；用椭球前向传播多步可达集，终端回到 backup 安全集 $\mathcal{X}_{\mathrm{safe}}$。

工程含义：用 GP 方差直接驱动约束收紧——不确定大的地方 tube 变粗，确定的地方变细，自适应保守性。

GP-MPC 的完整推导链¶

为什么用 GP 而非神经网络学残差？ GP 的核心优势不在于预测精度（NN 往往更好），而在于**提供校准的不确定性估计**。对于安全关键的机器人控制，知道"模型不确定在哪里"比知道"模型预测值是多少"更重要。GP 的后验方差 $\sigma^2(x,u)$ 直接量化了"在 $(x,u)$ 处模型有多不确定"——这就是 tube 半径的自然来源。

GP 残差动力学的完整设置：

已知标称模型（如 URDF 动力学）：$x_{k+1} = f_{\text{nom}}(x_k, u_k)$

真实动力学包含未建模项：$x_{k+1} = f_{\text{nom}}(x_k, u_k) + g(x_k, u_k)$

其中 $g: \mathbb{R}^{n_x+n_u} \to \mathbb{R}^{n_x}$ 是未知残差（气动力、摩擦、柔性等）。

GP 建模：$g \sim \mathcal{GP}(\mathbf{0}, k(\cdot, \cdot))$，核函数通常选 Squared Exponential (SE)： $$k(\mathbf{z}, \mathbf{z}') = \sigma_f^2 \exp\left(-\frac{1}{2}(\mathbf{z}-\mathbf{z}')^\top \Lambda^{-1} (\mathbf{z}-\mathbf{z}')\right)$$

其中 $\mathbf{z} = [x^\top, u^\top]^\top$，$\Lambda = \text{diag}(l_1^2, \ldots, l_{n_x+n_u}^2)$ 是长度尺度矩阵。

后验预测（给定 $N_D$ 个数据点 $\{(\mathbf{z}_i, g_i)\}_{i=1}^{N_D}$）： $$\mu_g(\mathbf{z}^*) = k_*^\top (K + \sigma_n^2 I)^{-1} \mathbf{g}$$ $$\sigma_g^2(\mathbf{z}^*) = k(\mathbf{z}^*, \mathbf{z}^*) - k_*^\top (K + \sigma_n^2 I)^{-1} k_*$$

约束收紧公式：利用 Chowdhury-Gopalan (2017) 的 RKHS 浓度界： $$\Pr\{\|g(\mathbf{z}) - \mu_g(\mathbf{z})\| \le \beta_\delta \sigma_g(\mathbf{z}), \forall \mathbf{z}\} \ge 1-\delta$$

其中 $\beta_\delta = B + 4\sigma_f\sqrt{\gamma_{N_D} + 1 + \ln(1/\delta)}$（$B$ 是 RKHS 范数上界，$\gamma_{N_D}$ 是信息增益）。

因此安全约束变为： $$h(x) + \beta_\delta \sigma_g(x, u) \le 0$$

直觉：在数据稠密的区域，$\sigma_g$ 小，约束几乎不收紧（模型可信）；在数据稀疏的区域，$\sigma_g$ 大，约束大幅收紧（模型不可信，保守行驶）。这正是"uncertainty-aware"控制的本质。

GP-MPC 在线字典管理¶

GP 推断的复杂度是 $O(N_D^3)$（Cholesky 分解）。对于在线应用（>100 Hz），$N_D$ 必须严格限制。

信息增益筛选（Kabzan 2019 AMZ）：新数据点 $\mathbf{z}_{\text{new}}$ 只有当其信息增益超过阈值时才加入字典： $$\Delta I = \frac{1}{2}\log(1 + \sigma_n^{-2}\sigma_g^2(\mathbf{z}_{\text{new}})) > \epsilon_I$$

物理含义：只保留"在模型不确定的区域"采集的数据——在已经建模准确的区域（$\sigma_g$ 小）采集的数据不增加信息量，丢弃不影响预测质量。

典型参数（四旋翼高速飞行，Torrente 2021）： - 字典大小：$N_D = 300$–$500$ 点 - GP 推断时间：$< 1$ ms（预计算 Cholesky） - 信息增益阈值：$\epsilon_I = 0.01$ - 核长度尺度：$l \sim [0.5, 0.5, 0.1, 0.1]$（位置方向宽、速度方向窄）

从 GP-MPC 到安全 RL 的桥梁¶

GP-MPC 可以看作是 Model-Based RL (MBRL) 在安全控制中的特例：

MBRL	GP-MPC
学习世界模型	学习残差动力学 $g$
模型不确定性驱动探索	模型不确定性驱动约束收紧
最大化奖励	最小化代价 + 保证约束
无安全保证	高概率安全保证（RKHS 界）

Predictive Safety Filter (PSF) 是两者的桥梁：让 RL 策略自由优化性能，GP-MPC 安全层只在必要时介入——最小化对 RL 策略的修正，同时保证安全约束不被违反。

Predictive Safety Filter (Wabersich–Zeilinger 2021, Automatica 129:109597, arXiv:1812.05506)¶

给定任意 RL 策略 $u_L$（不安全），在线解 $$u^*=\arg\min_u\|u-u_L\|^2\quad\mathrm{s.t.}\ \text{名义动力学,收紧约束,终端安全集}.$$

不可行时回退 backup。RL 落地的即插即用安全模块——与 RL 策略完全解耦，策略可任意改而安全性保持。

LBMPC (Aswani–Gonzalez–Sastry–Tomlin 2013, Automatica 49(5):1216–1226, DOI:10.1016/j.automatica.2013.02.003)¶

双模型：名义线性 $\hat f$ 保证约束 + Robust 稳定（用 min-RPI 紧化）；性能模型 $\hat f+g(x,u)$（学习项）仅用于优化 cost。安全-性能解耦。

机器人实证¶

四旋翼高速：Torrente–Kaufmann–Föhn–Scaramuzza 2021 RA-L，GP 残差补偿气动力，高速跟踪误差降 70%。
自动驾驶赛车：Kabzan–Hewing–Liniger–Zeilinger 2019 RA-L，MPCC+GP，在线字典 300 点（信息增益筛选），AMZ 无人赛车实车圈速显著提升。
越野机器人：Ostafew–Schoellig–Barfoot 2016 IJRR，RC-LB-NMPC，Grizzly 900 kg 越野 2 m/s 实验。

与 domain randomization/RARL 的关系¶

Domain Randomization（Tobin 2017, arXiv:1703.06907）：训练时随机化物理参数 $\xi\in[\xi^{\mathrm{lo}},\xi^{\mathrm{hi}}]$——本质是**有限样本 scenario training**；等价于把 $\Theta$ 建模为 $\xi$ 的凸包，策略对所有样本鲁棒。
RARL（Pinto 2017, arXiv:1703.02702）：零和博弈 $\max_\mu\min_\nu\mathbb{E}[\sum r]$——对标 closed-loop min-max 的 RL 版本。
推荐组合：训练用 RARL/DR 得到策略 $\pi_\theta$，部署时 MPC safety filter 兜底——sim-to-real 的 2026 年最稳范式。

常见陷阱¶

陷阱	现象	处理
GP 字典爆炸	实时性破坏	信息增益阈值筛选，字典大小固定
$\beta_\delta$ 常量选择	理论值过大不可用	实用 $\beta\in[1,3]$（经验 3σ），放弃严格理论
Learning 项进入约束	约束违反率失控	Aswani LBMPC 原则：学习只入 cost，不入约束
Safety filter 震荡	$u^*$ 跳变	加平滑正则 $\\|u-u_{\mathrm{prev}}\\|^2$

桥梁¶

→ 第六批 RL / Safe RL / MBRL：PSF 直接作为 RL safety shield；GP-MPC 与 MBRL 的 GP 世界模型共享框架。
→ 第五批 SLAM：用估计协方差 $\Sigma_{\hat x}$ 驱动 tube 收紧，实现 SLAM-MPC 耦合（Muntwiler 2022 differentiable MHE）。

§3.13.A Constraint Tightening 的精确计算（选学）¶

多面体 Minkowski 差的算法¶

对 $A=\{x:Fx\le g\}$，$B=\{x:Hx\le k\}$（V-多面体 $\mathrm{conv}\{b_1,\dots,b_p\}$）：$A\ominus B=\{x:Fx\le g-h_B(F^\top)\}$，$h_B(F_i^\top)=\max_j F_ib_j$。复杂度 $O(mp)$，其中 $m$ 为约束行数、$p$ 为顶点数。

Zonotope 近似（CORA toolbox）¶

Zonotope：$Z=c\oplus\bigoplus_{i=1}^p[-g_i,g_i]=\{c+\sum_i\alpha_ig_i:\alpha_i\in[-1,1]\}$，由中心 $c$ 与生成子 $g_i$ 表示。Minkowski 和极快：$(c_1,G_1)\oplus(c_2,G_2)=(c_1+c_2,[G_1\,G_2])$；Minkowski 差通过 interval enclosure 近似。

优势：表示复杂度 $O(p)$ 而非 V/H-多面体的 $O(2^n)$；线性变换 $AZ=(Ac,AG)$ 保 zonotope 性质。适合高维（$n\ge 10$）可达集计算。

椭球近似¶

$E=\{x:(x-c)^\top P^{-1}(x-c)\le 1\}$。Minkowski 和 $E_1\oplus E_2\subseteq\tilde E$：存在 $p\in[0,1]$ 使 $\tilde P=(1+p)P_1+(1+1/p)P_2$（minimum-volume outer ellipsoid）。LMI 优化求 $p^*$。

工具对照¶

工具	表示	适合问题	语言
MPT3	H/V-多面体	精确 mRPI、显式 MPC	MATLAB
CORA	Zonotope、多项式 zonotope、Taylor 模型	非线性可达、验证	MATLAB
Polyhedra.jl	H/V-多面体	Julia 生态	Julia
pypoman	H/V-多面体	Python	Python

§3.13.B 非线性 Tube MPC：Contraction/Incremental Lyapunov¶

Köhler–Müller–Allgöwer 2020（IEEE TAC 65(7):3576–3583, DOI:10.1109/TAC.2019.2949350, arXiv:1909.12765）¶

增量可稳假设：存在度量 $V_\delta(x,z)$（incremental Lyapunov）与局部反馈 $\kappa_\delta(x,z,v)$ 使 $$V_\delta(f(x,\kappa_\delta(x,z,v)),f(z,v))\le\rho V_\delta(x,z),\quad\rho\in[0,1).$$

等价于存在 Control Contraction Metric (CCM)（Manchester–Slotine 2017 TAC）$M(x)$ 使 $\dot M+A^\top M+MA\preceq -2\lambda M$。

核心贡献：参数化 terminal ingredients——对任意可达稳态 $(x^r,u^r)$ 族，离线构造 $V_f(x;x^r)=(x-x^r)^\top P(x^r)(x-x^r)$ 满足 $$V_f(f(x,\kappa_f(x;x^r));x^{r+})-V_f(x;x^r)\le -\ell(x,\kappa_f(x;x^r);x^r,u^r).$$

在线 NMPC 复杂度与名义同阶（单层），不需对每条参考轨迹重做离线设计。

与 §3.7.E contraction theory 对接¶

Contraction 度量 $M(x)$ 生成 tube 半径函数 $\rho(x,w)$；机器人实现：Singh–Majumdar–Slotine–Pavone 2017 "Robust Online Motion Planning via Contraction Theory and CCM"（RSS）。

工程含义¶

非线性 tube MPC 在 2020 年后实用化的关键——不再需要 Mayne 2011 双层 NMPC。机器人中 CCM 已用于四旋翼 agile flight、操纵机械臂轨迹跟踪、腿足质心轨迹。

Contraction Theory 与 Tube 半径的精确关系¶

回顾 contraction（收缩性）：系统 $\dot x = f(x,t)$ 在度量 $M(x)$ 下收缩，如果**任意两条轨迹之间的距离以指数率衰减**： $$\frac{d}{dt} d_M(x_1(t), x_2(t)) \le -\lambda \cdot d_M(x_1(t), x_2(t))$$

其中 $d_M$ 是 Riemannian 距离，$\lambda > 0$ 是收缩率。

与 tube 半径的联系：如果标称轨迹 $\bar x(t)$ 满足 $\dot{\bar x} = f(\bar x, u_{\text{nom}})$，真实轨迹 $x(t)$ 满足 $\dot x = f(x, u) + w$（$\|w\| \le w_{\max}$），则在收缩度量下：

\[d_M(x(t), \bar x(t)) \le e^{-\lambda t} d_M(x(0), \bar x(0)) + \frac{w_{\max}}{\lambda}\]

极限 tube 半径 = $w_{\max}/\lambda$——完全类似于线性系统中 $\|Z_\infty\| \sim \|w_{\max}\| / (1-\|\text{eig}(A_K)\|)$！

这就是 Contraction Tube 的核心数学： 1. 离线用 SOS/SDP 求解 CCM $M(x)$ 和收缩率 $\lambda$ 2. Tube 半径 $= w_{\max}/\lambda$（在 $M(x)$ 度量下的球） 3. 在线只解一个 NMPC（与名义问题同阶），约束在 $M(x)$ 度量下收紧

与线性 Tube 的完美类比：

线性 Tube	Contraction Tube
$A_K$ Schur（特征值 < 1）	收缩率 $\lambda > 0$
mRPI $Z_\infty = \bigoplus A_K^i \mathcal{W}$	Tube 半径 $= w_{\max}/\lambda$
约束收紧 $\mathcal{X} \ominus Z$	约束在 $M(x)$-度量下收紧
反馈 $K$（固定线性矩阵）	CCM 反馈（状态依赖，来自 geodesic）
在线解一个 QP	在线解一个 NMPC

本质洞察：所有 Tube MPC（线性/Homothetic/Elastic/Contraction）都在做同一件事：用反馈把误差压缩到有界集合内，然后在收紧的约束内做名义优化。它们的区别只在于"如何参数化反馈"和"如何计算有界集合"。理解了这个统一视角，整个 Tube MPC 家族就通透了。

练习¶

练习 3.13.B.1 ⭐⭐⭐：对一维非线性系统 $\dot x = -x^3 + u + w$（$|w| \le 0.1$），用 Lyapunov 方法找收缩率 $\lambda$（提示：取 $V = x^2/2$，分析 $\dot V$）。计算 tube 半径。

练习 3.13.B.2 ⭐⭐⭐⭐：对 2D Dubins car $\dot x = v\cos\theta$, $\dot y = v\sin\theta$, $\dot\theta = \omega + w$（$|w| \le 0.1$），设计 contraction-based tube MPC。讨论收缩度量 $M(x)$ 的选择对 tube 形状的影响。

§3.13.C 多阶段鲁棒 MPC：介于 min-max 与 tube 之间¶

核心方法（Lucia–Finkler–Engell 2013, J. Process Control 23(9):1306–1319, DOI:10.1016/j.jprocont.2013.08.008）¶

Scenario tree：参数 $d\in\{d^{(1)},\dots,d^{(s)}\}$（顶点集合），每条路径 $i$ 独立决策 $u_k^{(i)}$；robust horizon $N_r$ 后分支冻结，总路径数 $\le s^{N_r}$。

OCP： $$\min_{\{x^{(i)},u^{(i)}\}}\sum_{i=1}^{s^{N_r}}\omega_i\sum_{k=0}^{N-1}\ell(x_k^{(i)},u_k^{(i)})\quad\mathrm{s.t.}$$ - 动力学：$x_{k+1}^{(i)}=f(x_k^{(i)},u_k^{(i)},d_k^{(i)})$； - 路径约束：$g(x_k^{(i)},u_k^{(i)})\le 0$； - 非预期性 (non-anticipativity)：$u_k^{(i)}=u_k^{(j)}$ 若路径 $i,j$ 在时刻 $k$ 前历史相同。

与 min-max / tube 对比¶

方法	recourse	保守性	在线复杂度
Open-loop min-max	无	最高	$\mathcal{O}(1)$ NLP（但指数约束）
Tube MPC	线性反馈 $K$	中高	$\mathcal{O}(N)$ QP
Multi-stage	完整 recourse（$N_r$ 内）	中	$\mathcal{O}(s^{N_r}\cdot N)$ NLP
Closed-loop min-max	完整 recourse	最低	$\mathcal{O}(q^N)$

典型 $N_r=1$–$2$，$s=2$–$3$——3–9 条路径，NLP 可实时。

Multi-stage 的完整构造与物理直觉¶

核心创新：Multi-stage 是 min-max 的"工程化版本"——只在前 $N_r$ 步允许分支（recourse），之后冻结到单条路径。物理直觉：短期内（$N_r$ 步 = 0.1-0.5 秒），扰动的不确定性最大（因为刚发生），控制器需要为每种可能性准备响应；长期（$N_r$ 步之后），扰动的影响被反馈逐渐消解，分支的价值递减——所以冻结后的保守性损失不大。

Scenario tree 的构造：

以 $s = 3$ 个顶点（如质量 $m \in \{m_{\text{low}}, m_{\text{nom}}, m_{\text{high}}\}$）、$N_r = 2$ 为例：

t=0         t=1         t=2         t=3 ... t=N
 │
 ├──m_lo──┬──m_lo───────────────────────── 路径1
 │         ├──m_nom──────────────────────── 路径2
 │         └──m_hi───────────────────────── 路径3
 │
 ├──m_nom─┬──m_lo───────────────────────── 路径4
 │         ├──m_nom──────────────────────── 路径5
 │         └──m_hi───────────────────────── 路径6
 │
 └──m_hi──┬──m_lo───────────────────────── 路径7
           ├──m_nom──────────────────────── 路径8
           └──m_hi───────────────────────── 路径9

总路径数 = $s^{N_r} = 3^2 = 9$。对比完整 min-max：$3^N = 3^{20} \approx 10^9$。

非预期性约束的含义：在 $t=0$ 时刻，所有 9 条路径共享同一个 $u_0$（因为此时尚不知道扰动是哪个值）。在 $t=1$ 时刻，路径 1-3 共享一个 $u_1^{(1)}$，路径 4-6 共享 $u_1^{(2)}$，路径 7-9 共享 $u_1^{(3)}$（因为只有第一步的扰动已实现）。这就是"信息可用时再做决策"的精确数学化。

跨领域类比：Multi-stage MPC 的 scenario tree 与金融中的"多期随机规划"完全同构——投资组合经理在不确定的市场状态（利率涨/跌/平）下逐步做决策，也是前几期分支、后几期冻结。控制工程和金融工程在这里共用同一套数学。

do-mpc 实现（Fiedler et al. 2023 CEP 140:105676）¶

import do_mpc
import numpy as np

# 1. 定义模型
model_type = 'continuous'
model = do_mpc.model.Model(model_type)

# 状态和输入
x = model.set_variable('_x', 'x', shape=(4,1))  # [位置, 速度, 角度, 角速度]
u = model.set_variable('_u', 'u', shape=(1,1))   # 推力

# 不确定参数
mass = model.set_variable('_p', 'mass')

# 动力学 (四旋翼简化)
dx = ca.vertcat(x[1], u/mass - 9.81, x[3], 0)
model.set_rhs('x', dx)
model.setup()

# 2. 配置 multi-stage MPC
mpc = do_mpc.controller.MPC(model)
mpc.set_param(n_horizon=20, n_robust=2, t_step=0.05)

# 不确定性的顶点值 (名义 ± 20%)
mpc.set_uncertainty_values(mass=np.array([0.8, 1.0, 1.2]))

# 代价函数
mterm = (x[0]-1)**2 + x[1]**2  # 终端代价
lterm = (x[0]-1)**2 + x[1]**2 + 0.01*u**2  # 阶段代价
mpc.set_objective(mterm=mterm, lterm=lterm)

# 约束
mpc.bounds['lower','_x','x'] = np.array([[-10], [-5], [-1], [-5]])
mpc.bounds['upper','_x','x'] = np.array([[10], [5], [1], [5]])
mpc.bounds['lower','_u','u'] = np.array([[-15]])
mpc.bounds['upper','_u','u'] = np.array([[15]])

mpc.setup()

后端 IPOPT + MA27/MUMPS；稳定性分析见 Lucia-Subramanian-Limon-Engell 2020 SCL。

计算量实测（do-mpc + IPOPT，Intel i7）：

配置	NLP 变量数	求解时间
$N=20, N_r=0$ (名义)	100	2 ms
$N=20, N_r=1, s=3$	300	8 ms
$N=20, N_r=2, s=3$	900	35 ms
$N=20, N_r=2, s=5$	2500	150 ms

工程经验：$N_r = 1$-$2$，$s = 2$-$3$ 是实时可行的上限（对 50 Hz MPC）。

与 Tube MPC 的精确权衡¶

维度	Rigid Tube	Multi-stage ($N_r=2, s=3$)
保守性	中（固定 tube）	低（2步 recourse）
在线复杂度	同名义 QP	9× 名义 NLP
离线设计	需要计算 mRPI + K	只需指定顶点
非线性支持	需要 LPV 或双层	直接支持
实时性	极好 ($< 1$ ms)	一般 (10-100 ms)
适用频率	200+ Hz	10-50 Hz

选型建议： - 高频控制（> 100 Hz）：Rigid Tube（不可替代） - 中频决策（10-50 Hz）+ 少量顶点：Multi-stage - 需要同时处理非线性 + 不确定性 + 长预测：Multi-stage + warm-start

练习¶

练习 3.13.C.1 ⭐⭐：对一维系统 $x_{k+1} = a x_k + b u_k$，$a \in \{0.8, 1.0, 1.2\}$，$N=5$，$N_r=2$，写出完整的 multi-stage NLP（列出所有路径的动力学约束和非预期性约束）。

练习 3.13.C.2 ⭐⭐⭐：用 do-mpc 实现上述四旋翼 multi-stage MPC，与名义 MPC 对比在质量突变（$m$ 从 1.0 跳到 1.3）时的跟踪性能。

§3.13.D Data-Driven MPC：DeePC 与 Willems 基本引理¶

Willems 基本引理（Willems–Rapisarda–Markovsky–De Moor 2005, SCL 54(4):325–329, DOI:10.1016/j.sysconle.2004.09.003）¶

定理 3.13.6。设 $\mathcal{B}$ 为可控 LTI 系统、阶 $n$、复杂度 $\ell$。若 $u^d\in\mathbb{R}^{mT}$ 持续激励阶 $L+n$（即 $H_{L+n}(u^d)$ 行满秩），则 $$\mathrm{colspan}\,H_L\!\left(\!\begin{bmatrix}u^d\\y^d\end{bmatrix}\!\right)=\mathcal{B}|_L:=\{\text{$\mathcal{B}$ 的长度 $L$ 轨迹}}.$$

即**任意长度 $L$ 的系统轨迹可表示为 $w=H_L(w^d)g$**。

DeePC（Coulson–Lygeros–Dörfler 2019, ECC, DOI:10.23919/ECC.2019.8795639, arXiv:1811.05890）¶

分块 Hankel $\begin{bmatrix}U_p\\U_f\\Y_p\\Y_f\end{bmatrix}=\begin{bmatrix}H_{T_{\mathrm{ini}}+N}(u^d)\\H_{T_{\mathrm{ini}}+N}(y^d)\end{bmatrix}$。优化： $$\min_{g,u,y,\sigma_y}\sum_{k=0}^{N-1}\|y_k-y_k^{\mathrm{ref}}\|_Q^2+\|u_k\|_R^2+\lambda_g\|g\|_2^2+\lambda_\sigma\|\sigma_y\|_1$$ s.t. $$\begin{bmatrix}U_p\\Y_p\\U_f\\Y_f\end{bmatrix}g=\begin{bmatrix}u_{\mathrm{ini}}\\y_{\mathrm{ini}}+\sigma_y\\u\\y\end{bmatrix},\ u\in\mathcal{U},\,y\in\mathcal{Y}.$$

假设：(i) LTI 可控；(ii) 持续激励阶 $\ge T_{\mathrm{ini}}+N+n$；(iii) $T_{\mathrm{ini}}\ge\ell$（lag）。

确定性等价：DeePC ≡ state-space MPC（由基本引理）。

正则化 = 隐式 DRO（Dörfler–Coulson–Markovsky 2023, IEEE TAC 68(2):883–897, DOI:10.1109/TAC.2022.3148374）¶

定理：$\lambda_g\|g\|_2^2$ 等价于对数据矩阵扰动 $\|\Delta\|_F\le\rho$ 的椭球 robust LS；Lipschitz 目标下 $$\sup_{\mathbb{Q}\in\mathbb{B}_\varepsilon^{W_1}(\hat{\mathbb{P}}_N)}\mathbb{E}^\mathbb{Q}[\ell]=\mathbb{E}^{\hat{\mathbb{P}}_N}[\ell]+\varepsilon\cdot\mathrm{Lip}(\ell)\cdot\|g\|_*.$$

即正则化 DeePC 就是一个 Wasserstein DR-DeePC。

工程含义¶

无需系统辨识直接 I/O 数据做 RHC；机器人应用：ETH 四旋翼（Elokda 2021）、赛车（Huang 2022）、同步电机（Huang 2023）、电网变换器。**适合模型难建但数据好采**的场景（四旋翼、电机），不适合接触富/强非线性（腿足、柔体）。

Willems 基本引理的深层含义¶

为什么这个结果如此惊人？ 传统控制设计的流程是：数据 → 辨识模型 → 用模型设计控制器。Willems 引理说：你可以完全跳过辨识——只要数据足够丰富（持续激励），Hankel 矩阵本身就包含了系统的所有信息。

跨领域类比：这类似于机器学习中的 kernel trick。传统方法先学特征表示再做回归，kernel 方法直接在数据空间中做回归（通过 Gram 矩阵）。DeePC 的 Hankel 矩阵就是控制系统的"Gram 矩阵"。

持续激励条件的物理含义：输入信号 $u^d$ 必须"覆盖"系统的所有模态。形式化地，$H_{L+n}(u^d)$ 行满秩意味着输入的 Hankel 矩阵包含了足够多的线性无关行——每个行对应一种"频率组合"。实践中，**白噪声输入**自然满足持续激励（以概率 1），但需要足够长的数据 ($T \ge (m+1)(L+n) - 1$)。

如果持续激励不满足会怎样？ Hankel 矩阵列空间退化——不再等于所有可行轨迹。DeePC 的约束 $Hg = [u_{\text{ini}}; y_{\text{ini}}; u; y]$ 可能无解或解不唯一，导致控制器漂移或失稳。这就是为什么 $\sigma_y$ 松弛项和 $\lambda_g$ 正则化在实际中不可或缺。

DeePC 的局限性与适用边界¶

系统类型	DeePC 适用性	原因
LTI（线性时不变）	完美适用	基本引理精确成立
缓慢时变 LTV	适用（加遗忘）	用滑动窗口更新 Hankel
弱非线性	可用（正则化）	$\lambda_g$ 隐式处理非线性残差
强非线性（腿足接触）	不适用	线性假设完全破坏
混杂系统（切换）	不适用	单个 Hankel 无法表示多模态
高维系统（$n > 20$）	困难	所需数据量 $T$ 与维度立方增长

与系统辨识 + MPC 的对比：

维度	SysID + MPC	DeePC
工程步骤	辨识 → 验证 → 设计 → 调参	采数据 → 调 $\lambda_g, \lambda_\sigma$
模型可解释性	有（状态空间矩阵）	无（黑箱 Hankel）
数据需求	可以分步采集	需一次性持续激励
鲁棒性	取决于辨识误差	内建正则化鲁棒性
非线性扩展	自然（NMPC）	困难（需特殊处理）
适合谁	有模型经验的控制工程师	想快速出结果的应用工程师

从 DeePC 到 Koopman MPC 的演进（2022-2025 前沿）¶

DeePC 的核心局限是线性假设。2022 年后的一个重要方向是将 Koopman 算子与数据驱动控制结合：

Koopman 算子的核心思想（Brunton 2022 综述）：任何非线性系统 $x_{k+1} = f(x_k, u_k)$ 在无穷维"可观测函数空间"中都可以表示为线性算子。如果能找到有限维近似（通过字典函数 $\psi(x) = [\psi_1(x), \ldots, \psi_L(x)]^\top$），则：

\[\psi(x_{k+1}) \approx A_K \psi(x_k) + B_K u_k\]

这将非线性系统"提升"到高维线性空间，然后可以直接应用 Willems 基本引理！

Koopman-DeePC（Lian et al. 2021, L4DC）：在提升空间中构建 Hankel 矩阵，用 $z_k = \psi(x_k)$ 替代原始 $y_k$，得到非线性系统的数据驱动预测控制。

工程价值：腿足机器人的接触动力学高度非线性且难建模，但如果能找到合适的 Koopman 字典（如径向基函数 + 多项式），就可以在不建模的情况下实现 MPC——这是 2024-2025 年的活跃研究方向。

练习¶

练习 3.13.D.1 ⭐⭐：对一维系统 $x_{k+1} = 0.9 x_k + 0.1 u_k + w_k$（$w \sim \mathcal{N}(0, 0.01)$），用 Python 实现 DeePC：采集 $T = 200$ 步持续激励数据，构建 Hankel 矩阵，求解一步 DeePC 并与 MPC（已知模型）对比跟踪性能。

练习 3.13.D.2 ⭐⭐⭐：在上述系统中，逐渐增大噪声方差 $\sigma_w^2 = 0.01, 0.05, 0.1, 0.5$，画出 DeePC 性能随噪声的退化曲线。找到使 DeePC 开始不如朴素 PID 的临界噪声水平。讨论 $\lambda_g$ 如何随噪声调整。

练习 3.13.D.3 ⭐⭐⭐⭐（跨章综合）：对 2-DOF 机械臂（§3.15 辨识专题中的模型），比较三种控制策略：(a) 辨识后 MPC；(b) DeePC；(c) GP-MPC。从数据效率、计算量、鲁棒性三个维度做定量对比。

§3.13.E Robust MPC 与 RL 的融合¶

三条融合路径¶

1. MPC-as-safety-shield + free RL：RL 策略 $\pi_\theta$ 自由探索，MPC PSF（§3.13.8）在线兜底。代表：Wabersich 2021 PSF、ETH safe RL 系列。

2. MPC-in-the-loop RL：MPC 作为可微策略（Amos 2018 OptNet, Agrawal 2019 CVXPY Layer），通过 MPC 反传梯度训练代价参数。代表：Gros–Zanon 2020 "Data-Driven Economic NMPC using RL"。

3. RL-as-warmstart for MPC：RL 策略作为 warmstart 给 MPC 求解器；或 RL 学终端值函数 $V_f$ 替代 Economic MPC 中的 turnpike。代表：Lowrey 2019 POLO、Bhardwaj 2022 blending MPC and RL。

与 sim-to-real 的连接¶

Domain randomization 在仿真采样 $\xi^{(i)}\sim\mathcal{U}(\Xi)$，训练策略对所有 $\xi^{(i)}$ 鲁棒——本质是 scenario approach，给出 $(1-\varepsilon,1-\beta)$ 概率保证（Campi–Garatti 样本界）。

Tube 思想的 RL 版本：把 RL 策略分解为 $\pi(x)=\pi_{\mathrm{nom}}(\bar x)+K(x-\bar x)$，$\bar x$ 由名义 RL 规划生成，$K$ 由额外的 tracking 网络学习——直接对标 tube MPC 的 $u=\bar u+K(x-\bar x)$。

RARL vs Closed-loop min-max MPC¶

RARL 的 $\min\max$ 优化对应离散时间零和博弈——与 closed-loop min-max MPC 同构，但用策略梯度近似求解（避免 $q^N$ 爆炸）。RARL 是 min-max MPC 在 RL 领域的"可扩展近似"。

Robust Bayesian Optimization¶

Koller GP-MPC 的 GP UCB 界与 Robust BO（Bogunovic 2018 "Adversarially Robust Optimization"）的 RKHS 界同源——安全 RL / MBRL 底层数学。

MPC + RL 融合的 2024-2025 工程实践¶

当前最成熟的架构（"RL 做规划、MPC 做兜底"）：

                  ┌──────────────────────────────────────┐
                  │           训练阶段 (仿真)              │
                  │                                      │
                  │  Domain Randomization + PPO/SAC      │
                  │  → π_θ(observation → action)          │
                  │                                      │
                  └────────────────┬─────────────────────┘
                                   │ 部署
                  ┌────────────────▼─────────────────────┐
                  │           部署阶段 (实机)              │
                  │                                      │
                  │  观测 → π_θ → u_RL (期望动作)          │
                  │                    │                  │
                  │                    ▼                  │
                  │  ┌────────────────────────────┐      │
                  │  │ MPC Safety Filter           │      │
                  │  │                            │      │
                  │  │ min ||u - u_RL||²          │      │
                  │  │ s.t. 动力学可行性           │      │
                  │  │      约束满足(tube/CBF)     │      │
                  │  │      终端安全集可达         │      │
                  │  │                            │      │
                  │  │ → u_safe (安全修正后)       │      │
                  │  └────────────────────────────┘      │
                  │                    │                  │
                  │                    ▼                  │
                  │              执行器                    │
                  └──────────────────────────────────────┘

安全层介入频率的实证数据（来自 ETH RSL 实验）：

场景	RL 策略质量	Safety filter 介入率	性能损失
训练后直接部署（无 DR）	差	30-50%	高
中等 DR 训练	中	5-15%	低
充分 DR + 良好 reward shaping	好	< 2%	极低
人为恶意输入	极差	80-100%	-

关键观察：如果 RL 策略本身足够好（充分训练），safety filter 几乎不介入——它只在"黑天鹅"事件（异常扰动、传感器失效、分布外状态）时起保护作用。这意味着 safety filter 的计算开销可以认为是"保险费"——平时不用但必须有。

Safe RL 的数学保证总结¶

方法	保证类型	适用场景	参考
PSF (Wabersich 2021)	确定性鲁棒	有界扰动 + 已知动力学	Automatica
GP-MPC Filter	高概率 ($1-\delta$)	未知残差 + RKHS 假设	Koller 2018
RARL (Pinto 2017)	经验鲁棒	对手训练	ICML
Constrained RL	期望违反率 ≤ ε	Lagrangian 方法	Tessler 2018
Formal verification	数学证明	低维 + 简单网络	SMT/SOS

2025 年前沿方向： 1. 可微 MPC 作为策略层（Amos 2018, Agrawal 2019）：把 MPC 的 KKT 条件微分，通过 MPC 反传梯度训练代价参数——让 RL 直接"学会如何设置 MPC 参数"。 2. Hamilton-Jacobi 安全值函数（§3.18 HJ 可达性）：离线计算安全值函数 $V(x)$，在线当 $V(x) \le 0$ 时切换到最优安全控制——与 CBF-QP 互补但计算更重。 3. Diffusion Policy + Safety Layer：2024 年的扩散策略生成连续动作分布，安全层对分布做截断投影。

练习¶

练习 3.13.E.1 ⭐⭐：实现一个最简单的 PSF：给定 cart-pole 的名义动力学和约束（$|x| \le 2.4$, $|\theta| \le 12°$），当 RL 策略输出的 $u_{\text{RL}}$ 会导致下一步违反约束时，用 QP 投影到安全控制。

练习 3.13.E.2 ⭐⭐⭐：设计一个实验，量化 safety filter 对 RL 策略性能的影响。在 MuJoCo 中训练一个四足 locomotion 策略，分别测量有/无 safety filter 时的速度、稳定性和关节限位违反次数。

练习 3.13.E.3 ⭐⭐⭐⭐（跨章综合）：结合 §3.13.3（Tube MPC）、§3.13.8（GP-MPC）和 §3.18（HJ 可达性），设计一个三层安全架构： - 外层：HJ 值函数判断当前状态是否安全 - 中层：Tube MPC 做中期规划 - 内层：CBF-QP 做紧急安全修正说明每层的时间尺度、计算频率和相互关系。

前置自测¶

在开始本章之前，请确认你能回答以下问题：

#	自测题	前置
Q1	MPC 四条件（终端集、终端代价、阶段代价、时域长度）分别起什么作用？	3.11
Q2	什么是递归可行性？为什么它对 MPC 闭环稳定性至关重要？	3.11
Q3	什么是 ISS 稳定性？它与 Lyapunov 稳定性的区别是什么？	3.7
Q4	QP 的标准形式是什么？如何把 MPC 写成 QP？	3.12
Q5	什么是 CBF？它如何通过单步 QP 保证安全集前向不变？	3.8

如果 Q1-Q3 答不出 ≥ 2 题，先回 3.7/3.11 复习。

§3.13.F 鲁棒 MPC 在机器人的部署案例¶

四足：ETH ANYmal Perceptive Locomotion¶

Grandia et al. "Perceptive Locomotion through Nonlinear MPC", IEEE T-RO 39(5):3402–3421, 2023, DOI:10.1109/TRO.2023.3275384。SLQ-MPC（连续时间 NMPC）+ elevation map 感知 + CBF 安全层。接触位置不确定用 tube 处理，地形高度用 SLAM 协方差收紧约束。

无人机：UZH RPG 高速飞行¶

Torrente–Kaufmann–Föhn–Scaramuzza 2021 RA-L "Data-Driven MPC for Quadrotors" (arXiv:2102.05773)。GP 残差补偿气动力，高速跟踪误差降 70%。Bauersfeld et al. RSS 2021 NeuroBEM——BEM 物理 + NN 残差混合。

自动驾驶：DR-MPC + DeePC¶

Kabzan–Hewing–Liniger–Zeilinger 2019 AMZ 赛车 MPCC+GP；Huang et al. 2022 赛车 DeePC 实证。

机械臂接触：Carron et al. 2019¶

Carron–Arcari–Wermelinger–Hewing–Hutter–Zeilinger 2019 RA-L "Data-Driven MPC for Trajectory Tracking with a Robotic Arm"——GP 学接触残差 + tube 收紧。Sleiman–Farshidian–Hutter 2021 RA-L 统一 loco-manipulation NMPC。

SLAM-MPC 耦合¶

Muntwiler et al. L4DC 2022 "Learning-based MHE through Differentiable Convex Optimization Layers"——MHE 估计作为可微层嵌入 MPC，估计协方差驱动约束收紧。

部署案例的深度分析：ETH ANYmal Perceptive Locomotion¶

系统架构（Grandia et al. 2023 T-RO）：

LiDAR + IMU → 地形高程图 → elevation mapping
     │                              │
     ▼                              ▼
VIO 状态估计 → x_hat ────→ NMPC (SLQ) ────→ 摆腿轨迹
                                │                    │
                                ▼                    ▼
                           全身控制 (WBC) ─────→ 关节力矩 → 执行器

鲁棒性在每层的体现：

层	不确定性来源	处理方法	对应本章内容
状态估计	传感器噪声、VIO 漂移	EKF 协方差 → 约束收紧	§3.13.5（概率约束）
地形感知	高程图不确定性	地形方差 → 足端安全距离	§3.13.8（GP 类似）
质心 MPC	接触位置/时序误差	Tube 思想 + 保守足底力	§3.13.3（Tube）
全身控制	模型误差 + 延迟	WBC 鲁棒 QP + 力矩限幅	§3.13.7（CBF-QP）

关键设计决策：为什么选 SLQ-NMPC（连续时间 DDP）而非 Tube MPC？ 1. 四足的非线性动力学（接触切换、浮基）使线性 Tube 保守性过大 2. SLQ-NMPC 的 RTI 策略天然提供一步 feedback（类似 tube 反馈 $K$） 3. 地形不确定性通过**收紧落脚点约束**处理（非显式 tube）

实机性能指标： - MPC 频率：50 Hz（SLQ + HPIPM） - WBC 频率：400 Hz（hierarchical QP） - 扰动恢复：侧向 60 N 推力，0.3 秒内恢复 - 地形适应：±15 cm 高度误差下稳定行走

工程经验：实际机器人部署中，鲁棒性的 80% 来自反馈结构（高频 WBC + 状态估计），只有 20% 来自显式鲁棒 MPC 公式化。一个好的名义 NMPC + 高频反馈控制器，往往比一个理论完美但频率低的鲁棒 MPC 更实用。这是"理论指导工程"而非"理论替代工程"的典范。

部署案例：自动驾驶赛车 DeePC (AMZ Racing 2022)¶

场景：ETH AMZ 无人赛车，时速 100+ km/h，低摩擦赛道，轮胎-地面模型高度非线性且难以精确建模。

为什么选 DeePC 而非 model-based MPC？ 1. 轮胎模型（Pacejka 魔术公式）参数在赛道不同位置不同——辨识跟不上变化 2. 实际侧向力与 Pacejka 预测差 20-30%——模型不够用 3. DeePC 只需要历史 I/O 数据——每圈自动更新 Hankel 矩阵

关键实现细节： - 数据长度：$T = 500$（50 Hz 采样，10 秒数据） - $T_{\text{ini}} = 5$（0.1 秒初始化窗口） - $N = 30$（0.6 秒预测） - $\lambda_g = 0.1$（中等正则化——数据质量中等） - $\lambda_\sigma = 100$（强松弛——噪声较大） - 求解器：OSQP，warm-start，求解时间 < 5 ms

圈速对比： | 控制器 | 最快圈速 | 安全性（完赛率） | |--------|---------|----------------| | MPC + Pacejka 模型 | 100% | 85% | | DeePC（本方案） | 97% | 95% | | 人类最快车手 | 95% | 90% |

DeePC 虽然圈速略慢（因为正则化引入保守性），但完赛率最高——鲁棒性的实际价值。

练习¶

练习 3.13.F.1 ⭐⭐：查阅 Grandia et al. 2023 的论文，画出 ANYmal 的完整控制架构图。标注每层的频率、输入输出和对应的理论方法。

练习 3.13.F.2 ⭐⭐⭐（跨章综合）：设计一个四旋翼的"风扰动拒绝"实验： - 使用 Tube MPC（§3.13.3）处理风的加法扰动 - 使用 GP-MPC（§3.13.8）在线学习风场模型 - 比较两种方法在阵风（瞬时 5 m/s）和恒定侧风（持续 3 m/s）下的表现差异

练习 3.13.F.3 ⭐⭐⭐⭐：设计一个完整的 sim-to-real MPC 部署流程： 1. 仿真阶段：用 domain randomization 训练 + scenario 理论验证 2. 桥接阶段：用 DeePC 从真实数据快速建立数据驱动控制器 3. 精细阶段：用 GP-MPC 在线学习残差并收紧约束写出每个阶段的数学保证和过渡条件。

§3.13.G 2024-2025 年鲁棒 MPC 前沿进展¶

Varying-Tube MPC（Wang et al. IEEE TAC 2024）¶

传统 Rigid Tube 用**固定截面** $Z$ 约束所有时步。2024 年的核心突破是允许管道截面**随预测时步变化**：

\[e_k \in Z_k, \quad Z_{k+1} \supseteq A_K Z_k \oplus \mathcal{W}\]

关键数学结果：设 $Z_0 = \{e_0\}$（初始误差已知为 0），则 $Z_k = \bigoplus_{i=0}^{k-1} A_K^i \mathcal{W}$——前几步 $Z_k$ 比 $Z_\infty$ 小得多！约束收紧量递增：$\mathcal{X} \ominus Z_k$ 在 $k=0$ 时最宽松，逐步收紧到 $\mathcal{X} \ominus Z_\infty$。

可行域扩张：varying-tube 的初始可行域严格包含 rigid-tube 的初始可行域（$Z_0 = \{0\} \subset Z_\infty$）。对于短时域问题（$N$ 小），可行域扩张可达 30-50%。

Data-Driven Robust MPC（Buerger et al. 2025）¶

将在线数据与鲁棒 MPC 结合： 1. 在线收集扰动数据 $\{w_0, w_1, \ldots, w_K\}$ 2. 用 Wasserstein DRO 构造数据驱动的模糊集 $\mathcal{P}_{\text{amb}}$ 3. 在 Tube MPC 中用 $\mathcal{P}_{\text{amb}}$ 替代固定的 $\mathcal{W}$

核心创新：$\mathcal{W}$ 不再是固定的先验设计——它随数据积累而自动缩小（数据越多越不保守）。数学上，$\varepsilon_N \propto N^{-1/m}$（$N$ 为数据量），所以 tube 半径以 $O(N^{-1/m})$ 速率收缩。

Distributed Robust MPC（多机器人协同）¶

当多个机器人需要协同完成任务时，每个机器人的扰动通过通信延迟和模型差异传播到其他机器人。Distributed Tube MPC 把每个智能体的 tube 独立设计，通过**邻居 tube 的叠加**约束自身：

\[x_i \in \mathcal{X}_i \ominus Z_i \ominus \sum_{j \in \mathcal{N}(i)} Z_j^{\text{coupling}}\]

应用：无人机编队飞行（每架无人机独立 MPC + 邻居扰动估计）、多机械臂协作操作。

练习¶

练习 3.13.G.1 ⭐⭐⭐：对比 Rigid Tube 和 Varying-Tube 的可行域大小。对二维双积分器，画出两种方法的初始可行集（在相平面上），计算面积比。

教材对照表¶

教材	章节	本专题对应节
Rawlings–Mayne–Diehl 2017 第 2 版 Ch.3	3.2 Inherent Robustness; 3.3 Min-Max; 3.4 Tube-based (linear); 3.5 Tube-based (nonlinear); 3.6 Stochastic MPC	§3.13.1–3, §3.13.5, §3.13.B
Kouvaritakis–Cannon 2016（Springer）	Ch.3–5 Robust; Ch.6–10 Stochastic/Scenario	全专题核心参考
Borrelli–Bemporad–Morari 2017 Ch.15	Robust Linear MPC (mp-LP, min-max)	§3.13.2
Mesbah 2016 IEEE CSM 综述	Stochastic MPC overview	§3.13.5
Hewing et al. 2020 ARCRAS	Learning-based MPC	§3.13.8

关键定理清单（10 个）¶

Scokaert–Mayne 1998：反馈 min-max MPC ISS 稳定（IEEE TAC 43(8):1136, DOI:10.1109/9.704989）。
Mayne–Seron–Raković 2005 Rigid Tube：鲁棒指数稳定 + 约束满足（Automatica 41(2):219, DOI:10.1016/j.automatica.2004.08.019）。
Raković et al. 2005 mRPI：$(s,\alpha)$ 外近似 $Z_\infty\subseteq(1-\alpha)^{-1}F_s$（IEEE TAC 50(3):406, DOI:10.1109/TAC.2005.843854）。
Raković 2012 Homothetic Tube：参数化管道 + 递归可行（IEEE TAC 57(11):2746, DOI:10.1109/TAC.2012.2191174）。
Mayne–Kerrigan–van Wyk–Falugi 2011 Nonlinear Tube：双层 MPC 鲁棒稳定（IJRNC 21(11):1341, DOI:10.1002/rnc.1758）。
Köhler–Müller–Allgöwer 2020 Ref-Generic Terminal：参数化 terminal ingredients + 指数稳定（IEEE TAC 65(7):3576, arXiv:1909.12765）。
Calafiore–Campi 2006 Scenario：$N\ge\lceil(2/\varepsilon)(\ln(1/\beta)+d)\rceil$ 概率可行（IEEE TAC 51(5):742, DOI:10.1109/TAC.2006.875041）。
Campi–Garatti 2008 精确界：$\Pr\{V>\varepsilon\}\le\sum_{i=0}^{d-1}\binom{N}{i}\varepsilon^i(1-\varepsilon)^{N-i}$（SIAM J Opt 19(3):1211, DOI:10.1137/07069821X）。
Mohajerin Esfahani–Kuhn 2018 Wasserstein DRO：凸 tractable reformulation + $\varepsilon_N$ out-of-sample 保证（Math Prog 171:115, DOI:10.1007/s10107-017-1172-1）。
Willems 2005 Fundamental Lemma：持续激励 → Hankel 列空间 = 所有轨迹（SCL 54(4):325, DOI:10.1016/j.sysconle.2004.09.003）。

陷阱速查（12 条）¶

#	陷阱	场景	正确做法
1	开环 min-max 在非平凡 $\mathcal{W}$ 下一定不可行	$\mathcal{W}$ 含较大加法扰动	用 feedback min-max 或直接 tube
2	Tube $K$ 用过软 LQR 权重	$A_K$ 接近临界稳定 → $Z$ 巨大	$Q,R$ 使 $A_K$ 特征模 $\le 0.7$
3	mRPI $(s,\alpha)$ 近似精度不足	管道过大	精度 $\varepsilon\le 1\%\\|\mathcal{X}\\|$，$s\in[10,30]$
4	Scenario MPC 用 $\varepsilon<0.1\%$	$N_s$ 爆炸，QP 超实时	用 DRO 或 analytic Cantelli 替代
5	解析 chance constraint 用高斯假设处理非高斯	违反率远超 $\varepsilon$	Chebyshev–Cantelli 或 scenario
6	Wasserstein $\varepsilon$ 固定不随 $N$ 调	样本多了仍保守	$\varepsilon_N\propto N^{-1/m}$
7	GP-MPC 字典无上限	实时性破坏	信息增益阈值筛选
8	GP tube $\beta_\delta$ 用严格 RKHS 界	$\beta\gg 10$ 不可用	实用 $\beta\in[1,3]$
9	Learning 项进入约束	约束违反率失控	LBMPC 原则：学习只入 cost
10	Multi-stage $N_r$ 取过大	$s^{N_r}$ 爆炸	$N_r\in\{1,2\}$，$s\le 3$
11	DeePC Hankel $T_{\mathrm{ini}}<\ell$	初始态不唯一 → 漂移	$T_{\mathrm{ini}}\ge$ 系统 lag 上界
12	DeePC $\lambda_\sigma$ 用 L2 而非 L1	对离群样本敏感	$\\|\sigma_y\\|_1$（sparse slack）

自测题（5 道）¶

Q1（Tube MPC 与 mRPI）。给定 $A=\begin{bmatrix}1&1\\0&1\end{bmatrix}$，$B=\begin{bmatrix}0\\1\end{bmatrix}$，$\mathcal{W}=[-0.1,0.1]^2$，取 $K$ 使 $A_K$ 特征值在 $\{0.5,0.7\}$。(a) 写出 $Z_\infty$ 的定义；(b) 用 $(s,\alpha)=(5,0.1)$ 估计 $Z_\infty$ 的外近似；(c) 若 $\mathcal{X}=[-5,5]^2$，计算 $\mathcal{X}\ominus F(5,0.1)$ 的近似表达。

Q2（Scenario 样本复杂度）。对 $d=30$、$\varepsilon=2\%$、$\beta=10^{-5}$，用 Calafiore–Campi 估计 $N$；再用 Campi–Garatti 精确 Beta 尾分布数值反解 $N$，比较两者差距。

Q3（Wasserstein DR）。给定 $N=100$ 样本、$m=3$、$\beta=10^{-3}$、light-tail 常数 $c_1=c_2=1$，计算 $\varepsilon_N$。若 $\ell$ 的 Lipschitz 常数 $=2$，给出 out-of-sample cost 上界的表达式。

Q4（DeePC）。系统 $x_{k+1}=0.9x_k+0.1u_k+w_k$，$w\sim\mathcal{N}(0,0.01^2)$。采集 $T=200$ 条 I.I.D. 持续激励数据，$T_{\mathrm{ini}}=3$，$N=10$。(a) Hankel 矩阵维度；(b) 持续激励需要的最小阶；(c) 若 $\lambda_g=0.1$，讨论其隐式 DRO 半径。

Q5（Safety filter）。给定 RL 策略 $u_L=\pi_\theta(x)$ 和名义动力学 $\hat f(x,u)$，写出 predictive safety filter 的 OCP；证明在名义动力学 + 有限扰动 $\mathcal{W}$ 下递归可行性。

时间预算¶

节	内容	预计时间
§3.13.1	不确定性分类	3 h
§3.13.2	Min-max MPC	5 h
§3.13.3	Tube MPC 完整家族	12 h
§3.13.4	LPV Tube MPC	4 h
§3.13.5	随机 MPC 概述	4 h
§3.13.6	Scenario MPC	6 h
§3.13.7	DR-MPC	8 h
§3.13.8	Learning-based Robust MPC	8 h
档位 3 合计		50 h
§3.13.A	Tightening 精确算法	5 h
§3.13.B	非线性 tube (contraction)	6 h
§3.13.C	Multi-stage robust	5 h
§3.13.D	DeePC 数据驱动	8 h
§3.13.E	Robust MPC + RL	6 h
§3.13.F	机器人部署案例	5 h
档位 4 合计		35 h

资源附录¶

论文（按出现顺序）¶

Scokaert & Mayne 1998, IEEE TAC 43(8):1136, DOI:10.1109/9.704989
Mayne, Seron, Raković 2005, Automatica 41(2):219, DOI:10.1016/j.automatica.2004.08.019
Raković, Kerrigan, Kouramas, Mayne 2005, IEEE TAC 50(3):406, DOI:10.1109/TAC.2005.843854
Raković et al. 2012 Parameterized Tube, IEEE TAC 57(11):2746, DOI:10.1109/TAC.2012.2191174
Raković, Levine, Açıkmeşe 2016 Elastic Tube, ACC, DOI:10.1109/ACC.2016.7525471
Mayne, Kerrigan, van Wyk, Falugi 2011, IJRNC 21(11):1341, DOI:10.1002/rnc.1758
Köhler, Müller, Allgöwer 2020, IEEE TAC 65(7):3576, arXiv:1909.12765
Hanema, Tóth, Lazar 2016 CDC, DOI:10.1109/CDC.2016.7798472; Automatica 2017 85:137
Calafiore & Campi 2006, IEEE TAC 51(5):742, DOI:10.1109/TAC.2006.875041
Campi & Garatti 2008, SIAM J Opt 19(3):1211, DOI:10.1137/07069821X
Schildbach et al. 2014, Automatica 50(12):3009, DOI:10.1016/j.automatica.2014.10.035
Lucia, Finkler, Engell 2013, J Proc Control 23(9):1306, DOI:10.1016/j.jprocont.2013.08.008
Mesbah 2016, IEEE CSM 36(6):30, DOI:10.1109/MCS.2016.2602087
Mohajerin Esfahani & Kuhn 2018, Math Prog 171:115, arXiv:1505.05116
Coulson, Lygeros, Dörfler 2019 ECC DeePC, arXiv:1811.05890
Coulson, Lygeros, Dörfler 2022, IEEE TAC 67(7):3289, DOI:10.1109/TAC.2021.3097706
Willems, Rapisarda, Markovsky, De Moor 2005, SCL 54(4):325
Van Parys, Kuhn, Goulart, Morari 2016, IEEE TAC 61(2):430
Dörfler, Coulson, Markovsky 2023, IEEE TAC 68(2):883, arXiv:2101.01273
Hewing, Wabersich, Menner, Zeilinger 2020, ARCRAS 3:269
Koller, Berkenkamp, Turchetta, Krause 2018 CDC, arXiv:1803.08287
Wabersich & Zeilinger 2021, Automatica 129:109597, arXiv:1812.05506
Aswani, Gonzalez, Sastry, Tomlin 2013, Automatica 49(5):1216
Kabzan, Hewing, Liniger, Zeilinger 2019, IEEE RA-L 4(4):3363
Torrente, Kaufmann, Föhn, Scaramuzza 2021, IEEE RA-L 6(2):3769, arXiv:2102.05773
Ostafew, Schoellig, Barfoot 2016, IJRR 35(13):1547
Tobin et al. 2017 Domain Randomization, arXiv:1703.06907
Pinto et al. 2017 RARL, arXiv:1703.02702
Grandia et al. 2023 ANYmal Perceptive Loco, IEEE T-RO 39(5):3402

工具链¶

acados (github.com/acados/acados)：SQP_RTI + HPIPM，嵌入式 NMPC 主力；支持 soft constraint 近似 chance constraint。
MPT3 (mpt3.org)：MATLAB 多面体工具箱，mRPI/RPI/显式 MPC，Tube MPC 设计标配。
CORA (tumcps.github.io/CORA)：Zonotope/polynomial zonotope 可达集，非线性鲁棒分析。
do-mpc (www.do-mpc.com)：Python + CasADi，multi-stage NMPC/MHE，Lucia 组开发。
CasADi (web.casadi.org)：符号 AD + NLP 建模。
Drake (drake.mit.edu)：机器人动力学 + 优化 + 规划（MIT/TRI）。
DeePC 代码：github 上 deepc、DeePC-HUNT、Dörfler 组的 deepc-workshop。

教材¶

Rawlings, Mayne, Diehl, Model Predictive Control: Theory, Computation, and Design, 2nd ed., Nob Hill 2017/2020, ISBN 978-0-9759377-8-5。
Kouvaritakis & Cannon, Model Predictive Control: Classical, Robust and Stochastic, Springer 2016, DOI:10.1007/978-3-319-24853-0。
Borrelli, Bemporad, Morari, Predictive Control for Linear and Hybrid Systems, Cambridge 2017。

中文资源¶

知乎：搜索"模型预测控制学习笔记"、"DeePC 数据驱动预测控制"、"CasADi 入门"、"Tube MPC"、"分布鲁棒优化 DRO"；作者如"大饼博士"、"CyberSail"、"李柏"常有系统整理。
B 站：关键词"MPC 教程"、"鲁棒 MPC"、"DeePC"、"CasADi MATLAB"、"Sébastien Gros RL-MPC"（带中字）；DR_CAN 系列现代控制论基础。
CSDN：搜"do-mpc 教程"、"acados Python 安装"、"GP-MPC 无人机"——多为实现导向源码博客。
中国学者 MPC 课程：上交席裕庚《预测控制》（最早中文 MPC 教材）；港科大沈劭劼组四旋翼 MPC；西湖大学/浙大相关控制课程；北航段海滨无人系统控制；哈工大高会军非线性控制。

示例代码骨架¶

Rigid Tube MPC（Python + CasADi + MPT3-like mRPI）：

import casadi as ca, numpy as np
# 1. 离线：mRPI 外近似 Z = F(s,alpha)
def compute_mrpi(Ak, W_vertices, s_max=20, eps=1e-3):
    # 迭代 alpha^o(s) = max_j h_W((Ak^s)^T f_j)/g_j
    # 返回 (s, alpha, Z_vertices)
    ...
# 2. 收紧约束 X⊖Z, U⊖KZ
X_tight = pontryagin_diff(X, Z)
U_tight = pontryagin_diff(U, K @ Z)
# 3. 在线 QP（CasADi opti）
opti = ca.Opti('conic')
xbar = opti.variable(nx, N+1); ubar = opti.variable(nu, N)
for i in range(N):
    opti.subject_to(xbar[:,i+1] == A@xbar[:,i] + B@ubar[:,i])
    opti.subject_to(H_X_tight @ xbar[:,i] <= g_X_tight)
    opti.subject_to(H_U_tight @ ubar[:,i] <= g_U_tight)
opti.subject_to(H_Xf @ xbar[:,N] <= g_Xf)
opti.subject_to(x - xbar[:,0] in Z)  # 初始管道约束
opti.minimize(sum(stage_cost) + terminal_cost)
# 4. 控制律：u = ubar[:,0] + K @ (x - xbar[:,0])

DeePC（Python + cvxpy）：

import cvxpy as cp, numpy as np
Up, Yp, Uf, Yf = build_hankel(u_data, y_data, T_ini, N)
g = cp.Variable(Up.shape[1]); u = cp.Variable(N*nu); y = cp.Variable(N*ny)
sigma = cp.Variable(T_ini*ny)
cons = [Up@g == u_ini, Yp@g == y_ini + sigma, Uf@g == u, Yf@g == y]
obj = cp.sum_squares(y - y_ref) + R*cp.sum_squares(u) + lam_g*cp.sum_squares(g) + lam_s*cp.norm(sigma,1)
prob = cp.Problem(cp.Minimize(obj), cons); prob.solve()

桥梁总览（承前启后）¶

→ 3.7 ISS/Lyapunov：Tube MPC 的稳定性直接依赖 $A_K$ 的 ISS 增益；非线性 tube（§3.13.B）用 incremental ISS。
→ 3.8 CBF：Predictive Safety Filter 是 CBF 的预测视界扩展；DCBF = discrete CBF 是 MPC 约束的正则化视角。
→ 3.11 MPC 稳定性：本专题的 tube/stochastic 架构完全建立在 3.11 四条件之上，只是把名义模型换成"名义 + 管道"。
→ 3.12 数值求解：acados RTI + HPIPM 直接支持 tube MPC（只是更大的 QP）；scenario MPC 需特殊稀疏结构利用。
→ 第四批刚体动力学：$M(q)\ddot q+\cdots=\tau+\tau_{\mathrm{ext}}$，$\tau_{\mathrm{ext}}$ 是典型 $\mathcal{W}$；参数不确定性（$m,I$）直接对接 §3.13.4 LPV。
→ 第五批 SLAM 估计：SLAM 协方差 $\Sigma_{\hat x}$ 驱动 tube/chance constraint；differentiable MHE 把估计嵌入 MPC（Muntwiler 2022）。
→ 第六批 RL/Safe RL/MBRL：PSF = safe RL 的硬保障；DeePC = 非参数 MBRL 的控制侧对应物；domain randomization = scenario approach 的 RL 化身；RARL = closed-loop min-max 的策略梯度近似。

总结：Phase E MPC 的闭环¶

**第三批 13 个子专题（Phase E MPC）**从 3.1 LQR 基础、3.2–3.5 凸优化与 QP、3.6 KKT 与灵敏度、3.7 ISS 与 Lyapunov、3.8 CBF 安全、3.9–3.10 非线性规划与 SQP、3.11 MPC 四条件与稳定性、3.12 数值求解实时性，直至本专题 3.13 鲁棒与随机 MPC——完成了一条从**数学基础 → 名义理论 → 数值实现 → 真实世界**的完整路径。

核心收获：机器人工程师在 2026 年面对的**不是要不要用 MPC，而是用哪一层 MPC**。底层实时跟踪（1 kHz）用 3.12 的 RTI + tube（§3.13.3）；中层任务规划（100 Hz）用 NMPC + learning-based residual（§3.13.8）；上层行为决策（10 Hz）用 scenario / DR MPC（§3.13.6–7）甚至 DeePC（§3.13.D）。RL 策略作为性能加速器、MPC 作为安全保障——"RL 做规划、MPC 做兜底" 是后 ChatGPT 时代具身智能最成熟的控制栈。

本专题的三条主线再强调： 1. Tube 是工程主力——Rigid/Homothetic/Elastic/LPV/Contraction 的 tube 族覆盖从线性到非线性、从加法扰动到参数不确定的主要机器人场景；在线复杂度与名义 MPC 同阶。 2. 机会约束通过 scenario / DRO 解决——Calafiore–Campi 给出分布无关的 $O(d/\varepsilon\log(1/\beta))$ 样本复杂度；Wasserstein DRO 在数据驱动场景（DeePC）下等价于 L2 正则。 3. 学习与鲁棒融合——GP-MPC 用后验方差驱动 tube 收紧；PSF 使 RL 策略即插即用安全；domain randomization ≈ scenario approach；RARL ≈ closed-loop min-max 的可扩展近似。

最核心的认知跃迁：从 3.11 的"MPC 稳定性是名义数学性质"到 3.13 的"MPC 稳定性必须以某种概率/鲁棒意义定义在不确定性类上"。名义稳定只是初始条件可行性的代名词——真实的鲁棒稳定永远是相对于某个明确的不确定性集合（$\mathcal{W}$、$\Theta$、$\mathcal{P}_{\mathrm{amb}}$）而言的。在机器人工程中，谁把不确定性建得最准，谁的 MPC 就最鲁棒——这就是 Phase E 的最终定论。

至此，第三批 Phase E MPC 13 个子专题完整闭环。下一批（第四批刚体动力学）会把"$f(x,u)$"的结构从"任意光滑函数"精细化到"Lie 群上的 $M(q)\ddot q+\cdots=\tau$"，为全身 MPC 与接触动力学打基础；第五批 SLAM 与第六批 RL/MBRL 则分别从"状态哪里来"和"策略如何学"两个方向，与本专题的 tube/chance/learning 架构无缝对接。

Phase E 结束。下一站：第四批。

故障排查手册¶

症状	可能原因	排查步骤	相关节
Tube MPC 不可行（QP infeasible）	$Z$ 过大导致收紧后约束集为空	1. 打印 $\mathcal{X} \ominus Z$ 的体积；2. 检查 $A_K$ 特征值是否接近 1；3. 减小 $\mathcal{W}$ 或换更激进的 $K$	§3.13.3
Tube MPC 可行但性能很差	$Z$ 过大导致名义 MPC 过于保守	1. 比较名义 MPC 与 Tube MPC 的代价差异；2. 尝试 Homothetic Tube 替代 Rigid Tube；3. 检查 $\mathcal{W}$ 是否过估	§3.13.3.3
Scenario MPC 在线太慢	$N_s$ 过大或 QP 结构未利用	1. 用 Campi-Garatti 精确界重新计算 $N_s$；2. 利用 support rank 降低 $N_s$；3. 切到 DRO（SOCP 不随 $N$ 增长）	§3.13.6
DeePC 控制信号振荡	$\lambda_g$ 太小导致过拟合噪声	1. 增大 $\lambda_g$；2. 检查数据是否持续激励；3. 增加 $T_{\text{ini}}$	§3.13.D
GP-MPC 实时性不满足	GP 字典太大	1. 用信息增益阈值筛选；2. 固定字典大小（如 300 点）；3. 换 Sparse GP（如 FITC）	§3.13.8
Safety filter 控制跳变	CBF/tube 约束突然 active/deactive	1. 加松弛变量 $\delta$ 平滑过渡；2. 加 $\\|u - u_{\text{prev}}\\|^2$ 正则项；3. 检查 $h(x)$ 的 Lipschitz 常数	§3.13.8

跨章综合练习¶

综合练习 1（Tube MPC + ISS + CBF） ⭐⭐⭐

考虑一个四旋翼高度控制问题，同时受风扰动 $w \in [-0.5, 0.5]$ m/s$^2$ 和地面约束 $z \ge 0.2$ m。

(a) 设计 Rigid Tube MPC 保证高度约束在扰动下始终满足。写出 mRPI 集的计算过程。

(b) 在 Tube MPC 基础上，设计一个 CBF 安全层 $h(x) = z - 0.2$，使用 HOCBF（相对度 2）。写出约束 $\dot\psi_1 + \alpha_2(\psi_1) \ge 0$ 的具体形式。

(c) 比较 Tube MPC 单独使用和 Tube MPC + CBF 双层架构的性能差异。什么时候需要双层？

(d) 回顾 3.7 节 ISS 理论：证明 Tube MPC 闭环系统是 ISS 稳定的，给出 ISS 增益的具体表达式。

综合练习 2（Scenario MPC + Domain Randomization + RL） ⭐⭐⭐

(a) 一个四足机器人的模拟器支持随机化 8 个参数（4 条腿的摩擦系数和质量）。使用 Campi-Garatti 界，计算 $\varepsilon = 3\%$、$\beta = 10^{-4}$ 时需要多少 scenario。

(b) 将 scenario approach 理论应用于 domain randomization 训练的 RL 策略。如果训练用了 1000 个随机环境，能得到什么样的 out-of-sample 保证？

(c) 设计一个 MPC safety filter + RL 策略的联合架构。写出 safety filter 的 OCP，说明当 RL 策略输出不安全动作时，filter 如何最小化修正量。

本章常见误解汇总¶

误解	正确理解
Tube MPC 中的管道（tube）是固定大小的	Rigid Tube MPC 的管道确实固定（mRPI 集），但 Homothetic/Elastic Tube 允许管道大小随时间变化——可变管道在初始偏差大时提供更大可行域
Min-max MPC 是最佳的鲁棒 MPC 方法	Min-max MPC 是理论最优但计算量指数增长（对所有扰动序列穷举），工程中不可行；Tube MPC 通过解耦名义优化和误差界控制，用少量保守性换取计算可行性
随机 MPC（Chance-constrained）允许约束被违反所以不安全	随机 MPC 提供概率保证（如 99.5% 满足约束），在很多场景下比确定性鲁棒 MPC 更实用——后者对最坏情况的保守设计可能导致性能严重退化或不可行
DeePC 完全不需要系统模型	DeePC 用数据直接预测，但隐含假设系统是线性时不变（LTI）或通过正则化近似处理弱非线性；强非线性系统下 DeePC 的 Hankel 矩阵列空间无法覆盖真实行为
GP-MPC 可以处理任意未知动力学	GP 模型的预测能力仅限于训练数据覆盖的状态空间区域（局部性），外推能力差；且 GP 推断的 $O(N_\text{data}^3)$ 复杂度限制了在线数据量
Scenario MPC 的场景数越多越好	场景数 $N_s$ 增加提高概率保证但线性增加 QP 规模；Campi-Garatti 理论给出精确的 $N_s$-$\varepsilon$-$\beta$ 关系，允许以最小场景数达到目标置信水平
鲁棒 MPC 可以完全替代自适应控制	鲁棒 MPC 假设不确定集已知且有界；当不确定性的结构或范围未知时，需要在线辨识（自适应）来缩小不确定集，两者是互补而非替代关系

本章小结¶

方法	核心思想	保守性	计算量	适用场景
Min-max MPC	对所有扰动最优	最低（开环最高）	指数	理论基准
Rigid Tube	固定管道 + 约束收紧	中	同名义 MPC	机器人工程主力
Homothetic Tube	可变管道缩放	中低	多 $N$ 标量	可行域需扩张时
Elastic Tube	可变形状 + 反馈	低	$O(N^2)$	性能敏感场景
LPV Tube	非线性嵌入 LPV	中	LP	非线性系统近似
Scenario MPC	随机样本作硬约束	可调 $(1-\varepsilon)$	$O(N_s)$ QP	无模型/sim-to-real
DR-MPC	Wasserstein 球内最坏	可调 $\varepsilon_N$	SOCP	数据驱动
GP-MPC	GP 方差驱动 tube	自适应	GP 推断	学习型系统
DeePC	数据直接预测控制	$\lambda_g$ 控制	QP	模型难建的系统
Safety Filter	投影到安全控制	最小侵入	QP	RL 安全部署

本章目标¶

学完本专题后，你应能：

**区分**加法扰动、参数不确定和分布不确定三类不确定性，并为给定机器人选择合适的建模方式
完整推导 Rigid Tube MPC 的全部步骤：$K$ 设计 → mRPI 计算 → 约束收紧 → 在线 QP → 反馈律
实现 Tube MPC 的 Python/CasADi 代码，能在仿真中验证鲁棒性
使用 Scenario approach 和 Wasserstein DRO 处理概率性不确定性
理解 GP-MPC 和 Safety Filter 如何为 RL 策略提供安全保障
设计完整的"RL 策略 + MPC 安全层"部署架构
**判断**在给定工程约束下（频率、维度、不确定性类型），应选择哪种鲁棒 MPC 方法

延伸阅读¶

资料	内容	难度	推荐理由
Rawlings-Mayne-Diehl 2017 Ch.3	Robust MPC 完整理论	⭐⭐⭐	最权威的教科书
Kouvaritakis-Cannon 2016 全书	鲁棒/随机 MPC 专著	⭐⭐⭐	最全面的专题书
Hewing et al. 2020 ARCRAS	Learning-based MPC 综述	⭐⭐	10 页即掌握全貌
Mesbah 2016 IEEE CSM	Stochastic MPC 综述	⭐⭐	三大范式清晰对比
Dörfler-Coulson-Markovsky 2023 TAC	DeePC 正则化 = DRO	⭐⭐⭐⭐	最深刻的 DeePC 论文
Brunton et al. 2022 SIAM Review	Koopman 算子综述	⭐⭐⭐	数据驱动非线性
do-mpc 官方教程	Multi-stage NMPC 实现	⭐⭐	可运行代码
acados + CasADi tutorial	NMPC C++/Python 部署	⭐⭐	机器人实战

中文资源： - 知乎"大饼博士"DRO 系列（Wasserstein DRO 推导最清楚的中文资料） - 知乎"CyberSail"Tube MPC 系列（含 MATLAB 代码） - B 站 DR_CAN MPC 系列（基础推导） - CSDN "do-mpc 入门"（Python 实现导向）

累积项目：本章新增模块¶

在第三批累积项目"从零构建四旋翼 NMPC 控制器"中，本章新增**鲁棒性模块**：

模块	对应节	实现内容	预计代码量
Tube 层	§3.13.3	离线 mRPI + 在线约束收紧	300 行 Python
安全层	§3.13.8	CBF-QP + Predictive Safety Filter	200 行 C++
数据驱动层	§3.13.D	DeePC 备份控制器	400 行 Python
评估模块	全部	Monte Carlo 约束违反率统计	200 行 Python

项目代码位于 projects/quadrotor_mpc/robust/。

术语对照表¶

英文术语	中文	首次出现节
Tube MPC	管道模型预测控制	§3.13.3
Robust Positively Invariant (RPI)	鲁棒正不变集	§3.13.3.1
minimal RPI (mRPI)	最小鲁棒正不变集	§3.13.3.1
Pontryagin Difference	Pontryagin 差/Minkowski 差	§3.13.3.1
Support Function	支持函数	§3.13.3.1
Constraint Tightening	约束收紧	§3.13.3.2
Chance Constraint	机会约束/概率约束	§3.13.5
Scenario Approach	场景方法	§3.13.6
Distributionally Robust	分布鲁棒	§3.13.7
Wasserstein Distance	Wasserstein 距离	§3.13.7
Ambiguity Set	模糊集	§3.13.1
Persistent Excitation	持续激励	§3.13.D
Hankel Matrix	Hankel 矩阵	§3.13.D
Predictive Safety Filter	预测安全滤波器	§3.13.8
Control Contraction Metric	控制收缩度量	§3.13.B
Data-Driven MPC	数据驱动 MPC	§3.13.D
Fundamental Lemma (Willems)	基本引理	§3.13.D
Regularization	正则化	§3.13.D
Moment Ambiguity Set	矩模糊集	§3.13.7

⚠️ 教学专栏：鲁棒 MPC 三大深层陷阱¶

⚠️ 概念误区：认为 Tube MPC 的保守性仅来自 mRPI 近似

新手想法："只要用更多项近似 mRPI（增大 $s$），保守性就可以任意降低。"

现象/后果：即使 $s \to \infty$ 使 $Z \to Z_\infty$ 精确，Rigid Tube 仍然比理想 feedback min-max 保守——因为 Rigid Tube 强制所有扰动实现共享同一个 nominal trajectory $\bar{x}$。在不同扰动场景下，最优 nominal 路径本应不同。

根本原因：保守性有两个来源：(1) $Z$ 的外近似精度（可通过增大 $s$ 改善）；(2) 结构性保守——fixed tube shape 无法适应不同扰动序列。Homothetic tube（允许 $\alpha_k$ 伸缩管道半径）部分缓解 (2)，而 System Level Synthesis (SLS) tube 彻底消除 (2)，但计算量增大 $O(N^2 n_w)$。

正确认知：选择 tube 形状（rigid / homothetic / SLS / nonlinear contraction）本质上是**保守性 vs 在线计算量**的 Pareto 权衡，没有"最好"的选择。

🧠 思维陷阱：认为 Wasserstein DRO 半径 $\varepsilon$ 越小越好

新手想法："$\varepsilon$ 越小，模糊集越紧，优化结果越不保守，性能越好。"

现象/后果：$\varepsilon$ 太小时，out-of-sample 性能剧烈恶化。在 $N=50$ 样本下取 $\varepsilon = 10^{-6}$，模糊集退化为对经验分布 $\hat{\mathbb{P}}_N$ 的点估计——等价于 SAA（Sample Average Approximation），完全没有鲁棒性。

根本原因：$\varepsilon$ 的正确含义不是"用户可调的保守性旋钮"，而是"真实分布 $\mathbb{P}^*$ 以高概率落在模糊集内"的统计保证半径。正确的 $\varepsilon_N$ 由置信度 $\beta$、样本量 $N$ 和分布维度 $m$ 共同决定（如 $\varepsilon_N \sim N^{-1/\max(m,2)}$）。

正确做法：(1) 用 finite-sample 定理确定 $\varepsilon_N$；(2) 如果 $\varepsilon_N$ 太大导致过保守，应增加样本量 $N$ 而非人为压低 $\varepsilon$；(3) 可通过交叉验证选择 $\varepsilon$，但必须保留 out-of-sample test set 验证。

💡 概念误区：将 Scenario MPC 等同于"用 Monte Carlo 采样约束"

新手想法："Scenario MPC 就是随机生成 $N$ 个扰动样本，把每个样本对应的约束全加上。样本越多越好。"

现象/后果：盲目增大 $N$ 导致 QP 规模 $O(N \cdot n_{\text{constraints}})$ 爆炸，求解超时。而减少 $N$ 又无法满足 Calafiore-Campi 的覆盖率保证 $\Pr[\text{violation}] \le \varepsilon$。

根本原因：Scenario Approach 的数学基础不是大数定律（"越多越准"），而是**VC 理论的有限样本泛化界**。$N$ 的下界由 $(d, \varepsilon, \beta)$ 严格确定（$d$ = 决策变量维度）。盲目增减 $N$ 都会破坏理论保证。

正确做法：(1) 先用公式 $N \ge \frac{2}{\varepsilon}(\ln\frac{1}{\beta} + d)$ 计算最小 $N$；(2) 若 $N$ 太大，用 scenario removal / importance sampling 技术减小有效样本数；(3) 考虑 Garatti 2022 的 wait-and-judge 范式——先求解再用 a posteriori 界验证。

综合练习题¶

练习 13.A ⭐⭐（Tube MPC 完整闭环验证）：对双积分器 $x_{k+1} = \begin{bmatrix}1&0.1\\0&1\end{bmatrix}x_k + \begin{bmatrix}0.005\\0.1\end{bmatrix}u_k + w_k$，$w_k \in [-0.05, 0.05]^2$。(a) 设计 $K$ 使 $A_K$ 稳定且 $Z_\infty$ 体积合理；(b) 计算 mRPI $(s=20, \alpha=0.05)$ 外近似；(c) 收紧约束后构建 QP；(d) 仿真 200 步，画出名义轨迹 $\bar{x}$ 和实际轨迹 $x$，验证 $x - \bar{x} \in Z$ 始终成立。

练习 13.B ⭐⭐⭐（Scenario vs Analytic Chance Constraint 对比）：对相同的双积分器系统加上高斯扰动 $w \sim \mathcal{N}(0, 0.02^2 I)$，分别实现：(a) analytic chance constraint（假设高斯，用 $\Phi^{-1}(1-\varepsilon)$ 收紧）；(b) Scenario Approach（$\varepsilon = 5\%, \beta = 10^{-4}$，计算所需 $N$）。比较两者的可行域大小和在线计算时间。讨论：当真实分布是重尾的 Student-t 时，哪种方法更可靠？

练习 13.C ⭐⭐⭐（跨章综合 — 从 ISS 到 Tube 的完整链路）：回顾 §3.7（ISS 理论）和 §3.11（MPC 四条件）。对 Rigid Tube MPC 闭环系统：(a) 构造 ISS Lyapunov 函数 $V(x) = V_N^*(\bar{x}^*(x)) + \|x - \bar{x}^*(x)\|_P^2$；(b) 推导 ISS 增益 $\gamma(\|w\|)$；(c) 解释为什么管道半径 $r(Z)$ 正是 ISS 的 ultimate bound。

练习 13.D ⭐⭐（DeePC 数据需求量估计）：对一个 $n=4, m=2, p=2$ 的系统，$T_{\text{ini}} = 5$，prediction horizon $N = 15$。(a) 计算 Hankel 矩阵的最小列数 $T - T_{\text{ini}} - N + 1$，进而确定最小数据长度 $T$；(b) 持续激励（PE）条件要求 rank 至少为多少？(c) 如果系统是轻度非线性的（$\|f(x,u) - (Ax+Bu)\| \le \delta$），讨论 DeePC 正则化参数 $\lambda_g$ 应如何随 $\delta$ 调整。

练习 13.E ⭐⭐⭐⭐（Safety Filter 设计与验证）：给定一个训练好的 RL locomotion 策略 $\pi_\theta(x)$ 和四足机器人的 SRBD 模型，设计 predictive safety filter：(a) 定义安全集 $\mathcal{C} = \{x : h(x) \ge 0\}$（如质心高度、关节角度限制、ZMP 在支撑域内）；(b) 写出 safety filter 的 OCP 公式（最小化 $\|u - \pi_\theta(x)\|^2$ s.t. 轨迹不离开 $\mathcal{C}$）；(c) 证明在名义动力学 + 有界扰动下的递归可行性；(d) 讨论：如果动力学模型有误差（tube 外的状态），safety filter 还能保证安全吗？如何修复？

\(N\)	\(q=3\)（1D 扰动）	\(q=8\)（3D box 扰动）	\(q=64\)（6D box）
5	243	32,768	\(10^9\)
10	59,049	\(10^9\)	\(10^{18}\)
20	\(3\times 10^9\)	\(10^{18}\)	\(10^{36}\)

陷阱	现象	正确做法
\(K\) 选择过"软"	\(A_K\) 特征值接近单位圆 → \(Z\) 巨大 → 约束收紧后不可行	用 LQR 作初值，调权重矩阵 \(Q,R\) 使 \(Z\) 合理
低估 \(\mathcal{W}\)	仿真测试没出错，实车坏掉	\(\mathcal{W}\) 应包含建模误差 + 延迟 + 量化误差的联合估计
忘记 \(K_f\) 与 \(K\) 的区别	终端集用错误反馈计算	终端局部律 \(K_f\) 可独立于 tube 反馈 \(K\)，但 \(\mathcal{X}_f\oplus Z\subseteq\mathcal{X}\) 必须满足
mRPI 近似精度不够	\((s,\alpha)\) 取得太粗 → 管道过大 → 保守	选 \(\varepsilon\le 1\%\cdot\\|\mathcal{X}\\|\)，\(s\) 通常 10–20

线性 Tube	Contraction Tube
\(A_K\) Schur（特征值 < 1）	收缩率 \(\lambda > 0\)
mRPI \(Z_\infty = \bigoplus A_K^i \mathcal{W}\)	Tube 半径 \(= w_{\max}/\lambda\)
约束收紧 \(\mathcal{X} \ominus Z\)	约束在 \(M(x)\)-度量下收紧
反馈 \(K\)（固定线性矩阵）	CCM 反馈（状态依赖，来自 geodesic）
在线解一个 QP	在线解一个 NMPC

#	陷阱	场景	正确做法
1	开环 min-max 在非平凡 \(\mathcal{W}\) 下一定不可行	\(\mathcal{W}\) 含较大加法扰动	用 feedback min-max 或直接 tube
2	Tube \(K\) 用过软 LQR 权重	\(A_K\) 接近临界稳定 → \(Z\) 巨大	\(Q,R\) 使 \(A_K\) 特征模 \(\le 0.7\)
3	mRPI \((s,\alpha)\) 近似精度不足	管道过大	精度 \(\varepsilon\le 1\%\\|\mathcal{X}\\|\)，\(s\in[10,30]\)
4	Scenario MPC 用 \(\varepsilon<0.1\%\)	\(N_s\) 爆炸，QP 超实时	用 DRO 或 analytic Cantelli 替代
5	解析 chance constraint 用高斯假设处理非高斯	违反率远超 \(\varepsilon\)	Chebyshev–Cantelli 或 scenario
6	Wasserstein \(\varepsilon\) 固定不随 \(N\) 调	样本多了仍保守	\(\varepsilon_N\propto N^{-1/m}\)
7	GP-MPC 字典无上限	实时性破坏	信息增益阈值筛选
8	GP tube \(\beta_\delta\) 用严格 RKHS 界	\(\beta\gg 10\) 不可用	实用 \(\beta\in[1,3]\)
9	Learning 项进入约束	约束违反率失控	LBMPC 原则：学习只入 cost
10	Multi-stage \(N_r\) 取过大	\(s^{N_r}\) 爆炸	\(N_r\in\{1,2\}\)，\(s\le 3\)
11	DeePC Hankel \(T_{\mathrm{ini}}<\ell\)	初始态不唯一 → 漂移	\(T_{\mathrm{ini}}\ge\) 系统 lag 上界
12	DeePC \(\lambda_\sigma\) 用 L2 而非 L1	对离群样本敏感	\(\\|\sigma_y\\|_1\)（sparse slack）

不确定来源	数学建模	典型量级（机器人）
风扰（四旋翼）	加法 \(w\in\mathcal{W}\)，\(\\|w\\|_\infty\le 5\,\mathrm{m/s}\) 对应加速度扰动	3–8 m/s²
地面摩擦 \(\mu\)	参数 \(\theta\in[0.4,1.0]\)	±50%
负载质量 \(m\)	参数 \(\theta=m\in[m_0,m_0+\Delta m]\)	±20%
SLAM 位姿误差	加法 \(v_k\) 或分布 \(\mathbb{P}\in\mathcal{P}_{\mathrm{amb}}\)	位置 1–5 cm、yaw 0.5–2°
接触位置不确定	参数（接触点坐标）+ 加法（法向力扰动）	接触点 ±1 cm
执行器延迟	参数（延迟 \(\tau\)）	5–50 ms

性质	表达式	含义
正齐次	\(h_B(\lambda f) = \lambda h_B(f), \lambda \ge 0\)	方向缩放
次可加	\(h_B(f_1 + f_2) \le h_B(f_1) + h_B(f_2)\)	三角不等式
Minkowski 和	\(h_{A \oplus B}(f) = h_A(f) + h_B(f)\)	支持函数线性叠加
线性变换	\(h_{MB}(f) = h_B(M^\top f)\)	矩阵变换

陷阱	现象	正确做法
\(\Theta\) 取太大	LPV 退化为鲁棒 MPC	利用状态信息收缩 \(\Theta_{i\\|k}\)
忽略 \(\theta\) 的变化率	\(\theta_{k+1}\) 可能跳出 \(\Theta\)	加约束 \(\\|\theta_{k+1}-\theta_k\\| \le \delta\)
顶点太多	顶点 LMI 爆炸	用椭球包络替代多面体

范式	前提	保守性	在线复杂度	工程场景
解析（高斯）	\(w\) 高斯，\(g\) 仿射	紧	SOCP	过程控制、HVAC
解析（Cantelli）	知 \(\mu,\Sigma\)	保守	SOCP	不知分布但知矩
Scenario	可采样	\((\varepsilon,\beta)\) 可调	QP（规模 \(\propto N_s\)）	建筑、风机、自动驾驶
DRO	知模糊集	抗分布误判	SOCP/SDP	数据驱动、sim-to-real

维度	Tube MPC	Stochastic MPC (Gaussian)
不确定性模型	\(w \in \mathcal{W}\)（有界集）	\(w \sim \mathcal{N}(0, \Sigma_w)\)（无界）
约束处理	\(x \in \mathcal{X} \ominus Z\)（确定性收紧）	\(\Pr[x \in \mathcal{X}] \ge 1-\varepsilon\)（概率收紧）
保守性	对 \(\mathcal{W}\) 内所有扰动鲁棒（100%）	允许 \(\varepsilon\) 违反率（< 100%）
可行域	较小（为最坏情况预留）	较大（只为典型情况预留）
安全保证	绝对（\(\forall w \in \mathcal{W}\)）	概率（\(1-\varepsilon\) 概率）

陷阱	现象	正确做法
高斯假设用于非高斯扰动	实际违反率远超 \(\varepsilon\)	用 Chebyshev-Cantelli 或 scenario
忽略方差传播的闭环效应	开环方差单调增大导致过保守	用仿射反馈 \(u = \bar u + K(x-\bar x)\)
\(\varepsilon\) 设得太小	SOCP 约束过紧→不可行	对安全关键约束 \(\varepsilon=10^{-3}\)，非关键 \(\varepsilon=0.05\)
多步约束用联合概率	单步 \(1-\varepsilon\) 不保证全程 \(1-\varepsilon\)	用 Bonferroni \(\varepsilon_k = \varepsilon/N\) 或 conditional 方法

维度	蒙特卡洛	Scenario approach
目的	验证：估计 \(V(x)=\Pr\{f>0\}\)	合成：求满足概率约束的 \(x^*\)
使用样本方式	后验评估 \(\hat V=\frac{1}{N}\sum\mathbf{1}\{f>0\}\)	作为硬约束参与优化
收敛率	\(\mathcal{O}(1/\sqrt{N})\)	\(\mathcal{O}(d/\varepsilon\cdot\log(1/\beta))\)
概率保证	大数律 / Hoeffding	VC / support constraint / compression
典型 \(N\)（\(\varepsilon=5\%\)）	\(10^4\)（估尾概率）	\(10^2\)–\(10^3\)（合成决策）

Domain Randomization (RL)	Scenario MPC
物理参数 \(\xi \sim \mathcal{U}(\Xi)\)	扰动样本 \(w^{(i)} \sim \mathbb{P}\)
训练环境数量	Scenario 数量 \(N_s\)
策略对所有环境鲁棒	控制序列对所有 scenario 可行
无理论保证（经验）	\((1-\varepsilon, 1-\beta)\) 概率保证
训练时间长	在线 QP 规模大

维度	KL 球 \(\mathrm{KL}(\mathbb{Q}\\|\hat{\mathbb{P}}_N)\le\eta\)	Wasserstein 球 \(W_p\le\varepsilon\)
支撑	\(\mathbb{Q}\ll\hat{\mathbb{P}}_N\)（只在样本点加权）	样本邻域任意位置，支持连续
连续真分布	KL 到经验 = \(+\infty\)，保证破坏	\(W_1(\hat{\mathbb{P}}_N,\mathbb{P})=O(N^{-1/m})\)
重构	Exponential cone	SOCP / LP
控制应用	风险敏感、相对熵约束	DeePC、tube DR、sim-to-real

MBRL	GP-MPC
学习世界模型	学习残差动力学 \(g\)
模型不确定性驱动探索	模型不确定性驱动约束收紧
最大化奖励	最小化代价 + 保证约束
无安全保证	高概率安全保证（RKHS 界）

方法	recourse	保守性	在线复杂度
Open-loop min-max	无	最高	\(\mathcal{O}(1)\) NLP（但指数约束）
Tube MPC	线性反馈 \(K\)	中高	\(\mathcal{O}(N)\) QP
Multi-stage	完整 recourse（\(N_r\) 内）	中	\(\mathcal{O}(s^{N_r}\cdot N)\) NLP
Closed-loop min-max	完整 recourse	最低	\(\mathcal{O}(q^N)\)

配置	NLP 变量数	求解时间
\(N=20, N_r=0\) (名义)	100	2 ms
\(N=20, N_r=1, s=3\)	300	8 ms
\(N=20, N_r=2, s=3\)	900	35 ms
\(N=20, N_r=2, s=5\)	2500	150 ms

维度	SysID + MPC	DeePC
工程步骤	辨识 → 验证 → 设计 → 调参	采数据 → 调 \(\lambda_g, \lambda_\sigma\)
模型可解释性	有（状态空间矩阵）	无（黑箱 Hankel）
数据需求	可以分步采集	需一次性持续激励
鲁棒性	取决于辨识误差	内建正则化鲁棒性
非线性扩展	自然（NMPC）	困难（需特殊处理）
适合谁	有模型经验的控制工程师	想快速出结果的应用工程师

症状	可能原因	排查步骤	相关节
Tube MPC 不可行（QP infeasible）	\(Z\) 过大导致收紧后约束集为空	1. 打印 \(\mathcal{X} \ominus Z\) 的体积；2. 检查 \(A_K\) 特征值是否接近 1；3. 减小 \(\mathcal{W}\) 或换更激进的 \(K\)	§3.13.3
Tube MPC 可行但性能很差	\(Z\) 过大导致名义 MPC 过于保守	1. 比较名义 MPC 与 Tube MPC 的代价差异；2. 尝试 Homothetic Tube 替代 Rigid Tube；3. 检查 \(\mathcal{W}\) 是否过估	§3.13.3.3
Scenario MPC 在线太慢	\(N_s\) 过大或 QP 结构未利用	1. 用 Campi-Garatti 精确界重新计算 \(N_s\)；2. 利用 support rank 降低 \(N_s\)；3. 切到 DRO（SOCP 不随 \(N\) 增长）	§3.13.6
DeePC 控制信号振荡	\(\lambda_g\) 太小导致过拟合噪声	1. 增大 \(\lambda_g\)；2. 检查数据是否持续激励；3. 增加 \(T_{\text{ini}}\)	§3.13.D
GP-MPC 实时性不满足	GP 字典太大	1. 用信息增益阈值筛选；2. 固定字典大小（如 300 点）；3. 换 Sparse GP（如 FITC）	§3.13.8
Safety filter 控制跳变	CBF/tube 约束突然 active/deactive	1. 加松弛变量 \(\delta\) 平滑过渡；2. 加 \(\\|u - u_{\text{prev}}\\|^2\) 正则项；3. 检查 \(h(x)\) 的 Lipschitz 常数	§3.13.8

130_鲁棒与随机MPC

专题 3.13：鲁棒与随机 MPC — 从名义模型走向真实世界¶

摘要与本专题定位¶

§3.13.1 不确定性的分类与建模¶

定义（三类不确定性）¶

机器人物理含义映射表¶

不确定性量化的工程方法论¶

从名义 MPC 到鲁棒 MPC 的决策流程¶

练习¶

桥梁¶

§3.13.2 Min-max MPC：数学上的顶峰，工程上的陷阱¶

严格定义¶

核心定理（Scokaert–Mayne 1998，IEEE TAC 43(8):1136–1142，DOI:10.1109/9.704989）¶

计算复杂度爆炸¶

Bemporad–Borrelli–Morari 2003（mp-LP 方法）¶

工程含义¶

为什么 min-max 在工程上不实用——深层原因分析¶

Bemporad–Borrelli–Morari 2003 的工程启示¶

桥梁¶

练习¶

§3.13.3 Tube MPC 完整理论：机器人鲁棒控制的主力¶

3.13.3.1 预备：Pontryagin 差与 mRPI 集¶

不变集理论的来龙去脉¶

Pontryagin（Minkowski）差的完整推导¶

最小鲁棒正不变集 (mRPI) 的完整理论¶

定理 3.13.2（Raković–Kerrigan–Kouramas–Mayne 2005, IEEE TAC 50(3):406–410, DOI:10.1109/TAC.2005.843854）¶

3.13.3.2 Rigid Tube MPC（Mayne–Seron–Raković 2005）¶

Rigid Tube MPC 完整设计流程（从零到部署）¶

Tube MPC 时变截面的 2024 前沿（Sieber et al. 2024, Wang et al. 2024 TAC）¶

3.13.3.3 Homothetic Tube MPC（Raković 2012）¶

3.13.3.4 Elastic Tube MPC（Raković–Levine–Açıkmeşe 2016, ACC）¶

3.13.3.5 Nonlinear Tube MPC（Mayne–Kerrigan–van Wyk–Falugi 2011）¶

常见陷阱¶

Tube MPC 的 Python/CasADi 完整实现¶

练习¶

桥梁¶

§3.13.4 LPV Tube MPC：把非线性嵌入参数变化¶

严格定义¶

核心方法（Hanema–Tóth–Lazar 2016 CDC, DOI:10.1109/CDC.2016.7798472；Automatica 2017 85:137–144）¶

工程含义¶

LPV 嵌入的完整构造方法¶

常见陷阱¶

练习¶

§3.13.5 随机 MPC 概述：三大范式鸟瞰¶

严格定义¶

三大范式¶

核心对照（Mesbah 2016, IEEE CSM 36(6):30–44, DOI:10.1109/MCS.2016.2602087）¶

Open-loop vs Closed-loop 方差传播¶

与 Robust MPC 的关系¶

解析 Chance Constraint 的完整推导¶

Tube Stochastic MPC：两者的统一¶

常见陷阱¶

练习¶

§3.13.6 Scenario MPC：把概率约束变成 \(N_s\) 条硬约束¶

核心定理（Calafiore–Campi 2006，IEEE TAC 51(5):742–753，DOI:10.1109/TAC.2006.875041）¶

Scenario MPC（Schildbach–Fagiano–Frei–Morari 2014, Automatica 50(12):3009–3018, DOI:10.1016/j.automatica.2014.10.035）¶

与蒙特卡洛的本质区别¶

工程含义¶

Scenario MPC 的完整推导与直觉¶

Scenario MPC 的 Python 实现骨架¶

Scenario MPC 与 Domain Randomization 的精确对应¶

§3.13.7 分布鲁棒 MPC：Wasserstein 球内最坏情形¶

严格定义¶

核心定理（Mohajerin Esfahani–Kuhn 2018, Math. Program. 171:115–166, DOI:10.1007/s10107-017-1172-1）¶

DR-DeePC（Coulson–Lygeros–Dörfler 2022, IEEE TAC 67(7):3289–3304, DOI:10.1109/TAC.2021.3097706）¶

Robust–DR–Scenario 的连续谱¶

KL vs Wasserstein 对比¶

Wasserstein DRO 的完整推导链¶

从 DRO 到 DeePC 的正则化等价¶

§3.13.8 自适应 / 学习型鲁棒 MPC¶

框架综述（Hewing–Wabersich–Menner–Zeilinger 2020, ARCRAS 3:269–296, DOI:10.1146/annurev-control-090419-075625）¶

GP-MPC（Koller–Berkenkamp–Turchetta–Krause 2018, CDC, arXiv:1803.08287）¶

GP-MPC 的完整推导链¶

GP-MPC 在线字典管理¶

从 GP-MPC 到安全 RL 的桥梁¶

Predictive Safety Filter (Wabersich–Zeilinger 2021, Automatica 129:109597, arXiv:1812.05506)¶

LBMPC (Aswani–Gonzalez–Sastry–Tomlin 2013, Automatica 49(5):1216–1226, DOI:10.1016/j.automatica.2013.02.003)¶

机器人实证¶

与 domain randomization/RARL 的关系¶

常见陷阱¶