常微分方程（十六）：控制理论基础

开车时，你不断根据车道位置调整方向；恒温器将室温与设定值比较后调节加热器；火箭通过摆动喷管微调推力矢量以保持箭体垂直。 抽掉硬件，核心思想始终如一：测量、比较、动作。控制理论正是研究这一闭环过程的数学框架，而常微分方程（ODE）则是它的自然语言。

本章将此前积累的整个 ODE 工具箱——拉普拉斯变换（第 4 章 ）、线性系统（第 6 章 ）、特征值稳定性（第 7 章 ）以及非线性稳定性（第 8 章 ）——整合为一门统一的学科。其目标不再是被动地 描述 动力学行为，而是主动地 设计 它。

总结#

• 开环 vs. 闭环控制：为何反馈是核心思想
• 传递函数及其如何将 ODE 转化为代数运算
• PID 控制器与 Ziegler-Nichols 整定方法
• 根轨迹：闭环极点如何随控制器增益移动
• Bode 图、增益裕度与相位裕度
• MIMO 系统的状态空间表示
• 能控性、能观性及其秩判据
• 极点配置与 LQR 最优控制在倒立摆上的应用
• Luenberger 观测器与分离原理

前置知识#

• 第 4 章 ── 拉普拉斯变换（传递函数的来源）
• 第 6 章 ── 线性系统与矩阵指数
• 第 7 章 ── 特征值判断线性稳定性的准则

开环 vs. 闭环#

假设我们希望加热器将房间温度从 15°C 升至 22°C。

开环控制：施加固定的功率 $u(t)$ 并持续固定时间。若窗户开着，或室外温度比预期更低，最终温度就会偏离目标，且无法纠正。

闭环控制：实时测量温度，计算参考信号 $r$ 与实测输出 $y$ 之间的误差 $e(t) = r(t) - y(t)$ ，并让控制输入 $u(t)$ 成为误差的函数。系统由此具备自我修正能力，可抵抗建模误差和外部扰动。

我们将严格推导反馈带来的两个关键效果：

• 扰动抑制：阶跃型扰动不再导致稳态误差发散至无穷。
• 灵敏度降低：被控对象 $G$ 的相对变化对闭环传递函数 $T = CG/(1+CGH)$ 的影响被削弱为原来的 $S = 1/(1+CGH)$ 倍，其中 $S$ 称为 灵敏度函数。

传递函数#

微分方程中的微分操作变为乘以 $s$ ，积分变为除以 $s$ ，微积分问题由此转化为代数问题。

一阶系统#

阶跃响应为：$y(t) = K(1 - e^{-t/\tau})$ 。时间常数 $\tau$ 决定了系统跟踪设定值变化的速度。

二阶系统#

阻尼比 $\zeta$ 决定了响应的 定性 形态：

这四种情形广泛存在于电路、机械系统乃至生物系统中，因为任何线性系统在 稳定平衡点附近 的行为都近似于此。

PID 控制器#

每一项都有清晰的物理意义：

• 比例项（$K_p$ ）：如同弹簧，将输出拉向设定值。响应快，但无法消除恒定扰动下的稳态误差。
• 积分项（$K_i$ ）：累积历史误差的“记忆”。能将稳态误差精确归零，但引入一个原点处的极点，可能减慢响应并引发超调。
• 微分项（$K_d$ ）：对 未来 误差趋势作出反应的“阻尼器”。改善暂态性能，但会放大高频测量噪声。

Ziegler-Nichols 整定法（无需模型的起点）#

这是一种实用的经验法则，无需系统模型：

1. 设 $K_i = K_d = 0$ ，逐步增大 $K_p$ 直至闭环响应出现等幅振荡，记录此时的 临界增益 $K_u$ 和 振荡周期 $T_u$ 。
2. 按下表设置 PID 参数：

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22

import numpy as np

class PIDController:
    """带反饱和功能的离散 PID 实现"""
    def __init__(self, Kp, Ki, Kd, u_min=-np.inf, u_max=np.inf):
        self.Kp, self.Ki, self.Kd = Kp, Ki, Kd
        self.u_min, self.u_max = u_min, u_max
        self.integral = 0.0
        self.prev_error = 0.0

    def compute(self, error, dt):
        self.integral += error * dt
        deriv = (error - self.prev_error) / dt if dt > 0 else 0.0
        self.prev_error = error
        u = self.Kp*error + self.Ki*self.integral + self.Kd*deriv
        if u > self.u_max:
            self.integral -= error * dt
            u = self.u_max
        elif u < self.u_min:
            self.integral -= error * dt
            u = self.u_min
        return u

根轨迹 —— 闭环极点的走向#

根轨迹 描绘了当增益 $K$ 从 $0$ 扫描至 $\infty$ 时，闭环极点在复平面上的运动轨迹。它直观地揭示了速度（极点越靠左越快）与阻尼（极点越靠近实轴越平稳）之间的权衡。

需牢记两点：

• 根轨迹有 $n_p - n_z$ 条渐近线，从重心 $\sigma_a = (\sum p_i - \sum z_i)/(n_p - n_z)$ 出发，角度为 $(2k+1)\pi/(n_p - n_z)$ 。
• 根轨迹与虚轴的交点对应闭环极点 $\pm j\omega_c$ ，即 临界稳定 点。Routh-Hurwitz 判据可代数求解对应的临界增益 $K_{\text{crit}}$ 。

Bode 图与稳定裕度#

频率响应 $L(j\omega)$ 即传递函数在虚轴上的取值，通常以对数频率为横轴，分别绘制幅值（dB）和相位（度）。

增益穿越频率 $\omega_g$ 是 $|L|$ 等于 1（0 dB）的频率；相位穿越频率 $\omega_p$ 是相位等于 $-180^\circ$ 的频率。由此定义两个鲁棒性指标：

• 增益裕度（GM）：系统失稳前可将开环增益放大的倍数——在 $\omega_p$ 处读取幅值。
• 相位裕度（PM）：在 $\omega_g$ 处可额外容忍的相位滞后量。

经验法则：稳健设计应满足 GM > 6 dB（约两倍增益）且 PM > 30°（理想情况下 > 45°）。右侧 Nyquist 图提供了几何解释：系统稳定的充要条件是 Nyquist 曲线 不包围 临界点 $-1 + 0j$ 。

状态空间表示#

其中 $\mathbf x \in \mathbb R^n$ 为 状态向量——只要给定输入 $\mathbf u(t)$ ，即可预测系统未来行为。$A$ 为系统矩阵，$B$ 将输入耦合进状态，$C$ 提取可观测输出，$D$ 表示直通项（通常为零）。

稳定性、能控性与能观性#

以下三个秩条件决定了控制器是否可设计：

能控性意味着可通过适当输入将状态驱动至任意目标；能观性意味着可从输出历史中完全重构状态。二者在变换 $A \leftrightarrow A^T,\, B \leftrightarrow C^T$ 下互为对偶。

极点配置与 LQR#

若 $(A, B)$ 能控，则状态反馈律 $\mathbf u = -K\mathbf x$ 可使闭环系统变为 $\dot{\mathbf x} = (A - BK)\mathbf x$ 。此时可任意配置 $A - BK$ 的特征值。

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12

import numpy as np
from scipy.signal import place_poles

A = np.array([[0, 1, 0, 0],
              [0, -0.1, -1.96, 0],
              [0, 0, 0, 1],
              [0, 0.2, 23.5, 0]])  # 倒立摆（线性化）
B = np.array([[0], [1], [0], [-2]])

K_pp = place_poles(A, B, [-1, -2, -3, -4]).gain_matrix
print('极点配置增益 K =', K_pp)
print('闭环特征值 =', np.linalg.eigvals(A - B @ K_pp))

1
2
3
4
5
6

from scipy.linalg import solve_continuous_are

Q = np.diag([10, 1, 100, 1])      # 重罚位置和角度
R = np.array([[1.0]])              # 控制能量惩罚适中
P = solve_continuous_are(A, B, Q, R)
K_lqr = np.linalg.solve(R, B.T @ P)

相比极点配置，LQR 通常以更小的执行器能耗实现同等扰动抑制效果。只要 $(A, B)$ 可镇定且 $(A, \sqrt{Q})$ 可检测，LQR 解就保证闭环稳定。

观测器与分离原理#

估计误差 $\tilde{\mathbf x} = \mathbf x - \hat{\mathbf x}$ 满足 $\dot{\tilde{\mathbf x}} = (A - LC)\tilde{\mathbf x}$ 。若 $(A, C)$ 能观，则可任意配置 $A - LC$ 的特征值——通常将其置于比控制器极点更左的位置，确保估计先于控制动作收敛。

神奇的 分离原理 指出：控制器与观测器组合后的闭环极点，恰好是控制器极点（$A - BK$ 的特征值）与观测器极点（$A - LC$ 的特征值）的并集。因此二者可独立设计。

实战案例：小车上的倒立摆#

观察 $\text{eigvals}(A)$ 可知，其中一个特征值位于右半平面（开环下摆会倒下）。我们通过 Riccati 方程求解 $K_{\text{LQR}}$ ，并从 $0.2$ 弧度初始倾角开始仿真 $\dot{\mathbf x} = (A - BK)\mathbf x$ 。状态空间图最下方显示：摆角平滑衰减至零，小车回到原点，控制力最终降至一个微小的稳态推力——整个设计闭环仅需不到 30 行 Python 代码。

一个具体算例：汽车定速巡航的 PID 手算#

把抽象的 PID 公式落到一个我能手算的小例子上，比看一百行 scipy 代码更能建立直觉。设想一辆车，目标车速 $r = 70$ km/h，当前车速 $y(0) = 60$ km/h。控制器输出 $u$ 直接对应油门开度（百分比）。我取一个简化的一阶动态：每秒加速 $0.4u$ km/h（油门 100% 时一秒加 40 km/h，符合普通家用车从 0 到百公里 7 秒左右的量级），步长 $\Delta t = 1$ s。增益取 $K_p = 0.5$ ，$K_i = 0.1$ ，$K_d = 0.05$ 。

四步迭代手算：

第 0 步：$e_0 = 70 - 60 = 10$ ，积分项 $I_0 = 0$ ，微分项 $D_0 = 0$ （无前一刻误差）。$u_0 = 0.5 \cdot 10 + 0.1 \cdot 0 + 0.05 \cdot 0 = 5$ 。下一步速度 $y_1 = 60 + 0.4 \cdot 5 = 62$ 。

第 1 步：$e_1 = 8$ ，$I_1 = 0 + 10 \cdot 1 = 10$ ，$D_1 = (8 - 10)/1 = -2$ 。$u_1 = 0.5 \cdot 8 + 0.1 \cdot 10 + 0.05 \cdot (-2) = 4 + 1 - 0.1 = 4.9$ 。$y_2 = 62 + 0.4 \cdot 4.9 = 63.96$ 。

第 2 步：$e_2 = 6.04$ ，$I_2 = 10 + 8 = 18$ ，$D_2 = (6.04 - 8)/1 = -1.96$ 。$u_2 = 0.5 \cdot 6.04 + 0.1 \cdot 18 + 0.05 \cdot (-1.96) = 3.02 + 1.8 - 0.098 = 4.722$ 。$y_3 = 63.96 + 0.4 \cdot 4.722 = 65.85$ 。

第 3 步：$e_3 = 4.15$ ，$I_3 = 18 + 6.04 = 24.04$ ，$D_3 = (4.15 - 6.04)/1 = -1.89$ 。$u_3 = 0.5 \cdot 4.15 + 0.1 \cdot 24.04 + 0.05 \cdot (-1.89) = 2.075 + 2.404 - 0.0945 \approx 4.385$ 。$y_4 = 65.85 + 0.4 \cdot 4.385 \approx 67.60$ 。

注意三件事：第一，比例项 $K_p e$ 单调递减（误差缩小自然减小）；第二，积分项 $K_i I$ 反而单调递增——它在累积"还差多少没补上"的历史责任，最终会把稳态误差归零；第三，微分项一直为负，实际上是在主动减速，避免冲过头。如果只用 P，闭环极点是 $1 - 0.4 K_p = 0.8$ ，稳态误差非零；只有 $I$ 项才能消除偏差。

直觉与陷阱：积分饱和（windup）的真实代价#

PID 在课本上看着无害，工程上最常见的事故是积分饱和。还是上面这辆车：假设我突然把目标改成 $r = 200$ km/h（人为的极端值），一秒内 $e = 140$ ，积分项疯狂累积；几秒后即使车速冲到 200，已经"积"了大量正误差，控制器还在猛踩油门。等系统真的回到目标，得花几个周期才能把这部分历史"消化"掉，期间车速会大幅过冲。

我在博客早期写控制器的时候栽过一次：DC-DC 变换器的 PI 控制器在启动时电源电压未稳定，积分项饱和到 50 V 等价指令，输出电容硬生生被冲到额定值的 1.4 倍——电容鼓包。事后加了反饱和（anti-windup）：当输出 $u$ 已经到达执行器上下限时，停止积分累积。本章给出的代码片段里 compute() 函数那段 self.integral -= error * dt 就是这一招，叫条件积分或 back-calculation。

第二个陷阱是微分对噪声的放大。$D$ 项 $\dot e = -\dot y$ ，传感器噪声 $\sigma$ 在频率 $\omega$ 上贡献 $\omega \sigma$ ，高频噪声被直接放大 $K_d$ 倍打进油门。工程实践里几乎没人用纯微分，都是改用一阶低通滤波微分：$D(s) = K_d s / (1 + s/N)$ ，其中 $N \in [10, 30]$ 。Astrom 的 Feedback Systems 第 10 章对这个问题有详细推导。

应用与反例：什么时候 PID 不够#

PID 之所以占领工业界 90% 的回路，是因为它简单、便宜、可整定。但它的局限同样明显：

• 多输入多输出耦合系统：例如四旋翼无人机的姿态-位置耦合，每个 SISO PID 互相干扰，这时必须用 LQR 或 MPC 把整个 $(A, B)$ 一起设计。
• 大延迟过程：化工反应器从加料到温度上升常有 30-300 秒延迟，PID 的微分项失去前瞻能力，只能用 Smith 预估器或模型预测控制。
• 强非线性：磁悬浮系统、机器人在不同关节角度下的惯量矩阵会变 5-10 倍，单一组 PID 增益无法覆盖整个工作区。这时要么 gain scheduling（多组增益按工作点切换），要么用反馈线性化。
• 执行器明显饱和：飞行器的舵面有限位，PID 输出超过物理极限，就需要显式的约束控制——这正是 MPC 的核心优势。

反过来说，绝大多数家用空调、汽车定速巡航、3D 打印机温控、咖啡机、电烙铁恒温——这些一阶或二阶慢动态、单输入单输出、无显著延迟的场合，PID 整定一次就能跑十年不出事。这也是为什么 Astrom 在他那篇著名的综述里说：“PID is here to stay”。

现代图景#

本章实质上是以 目标导向 的方式重述了 ODE 的故事。第 7–8 章的稳定性理论回答的是“平衡点能否承受小扰动”；而控制理论则进一步回答“如何设计平衡点及其收敛速率”。第 4 章 的拉普拉斯变换不仅用于求解析解，如今还直接提供频域鲁棒性指标；第 6 章 的矩阵指数则成为状态空间仿真与观测器设计的核心工具。

前沿方向远不止于此：

• 鲁棒控制 —— $H_\infty$ 、$\mu$ -综合：在有界建模误差下保证性能。
• 自适应控制与模型预测控制（MPC） —— 控制器在线更新，每步求解一个滑动时域优化问题。
• 非线性控制 —— 反馈线性化、滑模控制、控制 Lyapunov 函数。
• 强化学习 —— 控制器通过经验 学习 而非人工设计；这与第 18 章 的 Neural ODE 自然衔接。

总结#

控制理论 将 ODE 从描述工具提升为设计语言。我们不再仅仅预测系统如何运动，而是明确告诉它该去往何处。

下一步#

下一章是 ODE 在物理与工程中的应用大全：经典力学（牛顿方程、Lagrange 方程、Hamilton 方程）、电路分析、热力学过程、流体力学、量子力学——它们的核心都是 ODE。下一章用一个章节连接所有这些应用，让你看到同一个数学工具如何统一描述各种物理现实。

工程实践片段：从 PID 到自动驾驶的横向控制#

说几个我自己上手过的真实案例，看 PID 在不同场景里如何变形。

第一个是无人机姿态环。横滚角 $\phi$ 的目标值由位置环给出，姿态环要在 200 Hz 的频率下跟踪。这里 $D$ 项不能省——四旋翼是欠阻尼系统，没有微分项姿态会震荡。但加速度计噪声大约 $0.1$ °/s，直接微分会把噪声放大成乱码。我们的方案是用陀螺仪角速度反馈代替差分，这等价于把 $D$ 项的滤波带宽精确锁在传感器带宽——本质上是把 PID 改写成 PI + 微分反馈。

第二个是车道保持的横向 PID。误差是车辆中心到车道中心的横向距离，控制量是方向盘转角。坑在哪里？车速变化会改变车辆-道路系统的传递函数：60 km/h 调好的增益放到 120 km/h 就会震荡。工程实践里要么做 gain scheduling（按车速分段切换增益），要么直接换成 LQR——后者在状态向量里把车速放进去，自动适配。

第三个是温控。家用电烙铁、烤箱、3D 打印机喷头——这些都是大热惯性 + 小延迟的一阶系统，PI 就够（不需要 D）。Ziegler-Nichols 整定法在这类系统上特别好用：把 $K_p$ 一点点拉大直到出现稳定振荡，记下临界增益和周期，套表就能直接用。我家的 3D 打印机每年要重新整定一次，因为加热块氧化、风扇灰尘累积都会改变热阻。

这三个例子的共同点：PID 不是孤立的算法，是一个可调节的家族。设计者的真正工作是判断系统结构、选择反馈量、设定带宽——这些都跳出了"$K_p, K_i, K_d$ 三个数"的范畴。但底层的 ODE 直觉（一阶 vs 二阶、阻尼比、稳态误差）对每个场景都适用。

一句话提醒：控制理论的本质不是数学#

我教学生时反复强调一件事：控制理论的难点从来不在拉普拉斯变换、根轨迹、Bode 图——这些是手段。真正的难点是建模。一个好的控制工程师 80% 的时间在思考"我的系统真的是 LTI 吗？““我的传感器在什么频段可信？““执行器有多少时滞和饱和？““扰动的功率谱长什么样？“剩下 20% 才是设计 PID 或 LQR。漂亮的数学方法用错系统比没用还糟糕——一个 LQR 控制器能让稳定的系统变成稳定且高性能的系统，但放在一个被错误辨识的系统上，它会高效率地把你的执行器烧掉。这是 Astrom 在所有著作里反复强调的一点，也是我对学生最严肃的告诫。

参考文献#

• Ogata, Modern Control Engineering, Pearson (2010).
• Franklin, Powell & Emami-Naeini, Feedback Control of Dynamic Systems, Pearson (2015).
• Astrom & Murray, Feedback Systems, Princeton (2008). （免费 PDF）
• Skogestad & Postlethwaite, Multivariable Feedback Control, Wiley (2005).