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

当你开车时不断根据车道位置纠正方向；恒温器对比室温和设定值后调节加热器；火箭通过摆动喷管让箭体保持垂直。 把硬件全部抽掉，剩下的是同一个想法：测量、比较、动作。控制理论就是研究这个闭环的数学，而它的母语正是常微分方程。

本章把整个 ODE 工具链——拉普拉斯变换（第 4 章）、线性方程组（第 6 章）、特征值稳定性（第 7 章）、非线性稳定性（第 8 章）——熔在一起：它们的目标不再是 描述 动力学，而是 设计 动力学。

本章要点

• 开环 vs 闭环：为什么反馈是最核心的思想
• 传递函数：把 ODE 变成代数
• PID 控制器与 Ziegler-Nichols 整定法
• 根轨迹：闭环极点随增益的变化
• Bode 图、增益裕度、相位裕度
• 状态空间表示与多入多出（MIMO）系统
• 能控性、能观性及其秩条件
• 极点配置和 LQR 最优控制（倒立摆案例）
• Luenberger 观测器与分离原理

前置知识

• 第 4 章 ── 拉普拉斯变换（传递函数的来源）
• 第 6 章 ── 线性方程组与矩阵指数
• 第 7 章 ── 特征值与线性稳定性

1. 开环 vs 闭环

让加热器把房间从 15℃ 升到 22℃。

开环：通电固定功率 $u(t)$、固定时间。如果窗户开着、或室外比预期更冷，最终温度就会偏离目标——而你毫无补救手段。

闭环：测量当前温度，计算误差 $e(t) = r(t) - y(t)$，把控制量 $u(t)$ 设计成误差的函数。系统会 自动纠正 模型误差和外部扰动。

反馈带来的两个核心后果，下文将逐步推出：

• 扰动抑制——一个常值扰动不再会让稳态误差发散。
• 灵敏度降低——被控对象 $G$ 的相对偏差只通过 $S = 1/(1+CGH)$ 的因子（灵敏度函数）影响闭环传递函数 $T = CG/(1+CGH)$。

2. 传递函数

对线性时不变（LTI）系统，输入 $u(t)$、输出 $y(t)$，零初始条件下的拉普拉斯比

把 ODE 中的求导变成乘 $s$、积分变成除 $s$——代数取代了微积分。

一阶系统

阶跃响应：$y(t) = K(1 - e^{-t/\tau})$。时间常数 $\tau$ 决定追踪速度。

二阶系统

阻尼比 $\zeta$ 决定 定性 行为：

这四种情形覆盖了几乎所有 稳定平衡点附近 的线性系统——电路、机械、生物，皆是如此。

3. 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:
    """带反饱和（anti-windup）的离散 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

4. 根轨迹 ── 闭环极点的轨迹

单位反馈下，开环传递为 $K\,L(s)$，闭环特征方程

根轨迹画出闭环极点随 $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}$。

5. Bode 图与稳定裕度

频率响应 $L(j\omega)$ 是同一个传递函数在虚轴上的取值。两个最自然的读数：幅值（dB）和相位（度），都画在对数频率轴上。

• 增益裕度（GM）：在失稳前还能把 $L$ 放大多少倍——读取 $\omega_p$ 处的幅值。
• 相位裕度（PM）：在 $\omega_g$ 处还能容忍多少额外相位滞后。

工程经验：GM > 6 dB（约两倍）、PM > 30 度（最好 > 45 度）才算稳健。右侧 Nyquist 图是几何视角：稳定的充要条件是曲线 不 包围临界点 $-1+0j$。

6. 状态空间表示

对多入多出（MIMO）系统，我们不再用传递函数兜圈子，直接写：

$\mathbf x \in \mathbb R^n$ 是 状态向量——只要给出 $\mathbf u(t)$，它就足以预测未来。$A$ 是动力学矩阵，$B$ 把输入耦合进来，$C$ 抽取测量输出，$D$ 是直通项（通常为 0）。

稳定性、能控性、能观性

三个秩条件回答了"能否设计控制器"这个根本问题：

能控性意味着可以用合适的输入把状态驱动到 任意 目标；能观性意味着可以从输出历史还原全部状态。它们在 $A \leftrightarrow A^T,\,B \leftrightarrow C^T$ 下互为对偶。

7. 极点配置与 LQR

如果 $(A, B)$ 能控，状态反馈 $\mathbf u = -K\mathbf x$ 让闭环变为 $\dot{\mathbf x} = (A - BK)\mathbf x$。可以选 $K$ 使闭环特征值落在 任意 期望位置。

 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))

LQR（线性二次型调节器）通过最小化二次代价

来选择 $K$，在状态偏差和控制能量之间权衡。解来自 代数 Riccati 方程

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)$ 可观测，就 保证 稳定。

8. 观测器与分离原理

状态反馈假设我们能 测量 $\mathbf x$ 的每一个分量；现实中只能看到 $\mathbf y$。Luenberger 观测器 用输出去估计未测量的状态：

估计误差 $\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$ 的特征值的并集。这意味着两部分可以 独立 设计。

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

线性化后状态 $\mathbf x = [x,\; \dot x,\; \theta,\; \dot\theta]^T$，矩阵

观察 $A$ 的特征值，其中一个落在右半平面（开环倒立摆会倒）。我们由 Riccati 解出 $K_\text{LQR}$，然后从 $0.2$ 弧度初始倾角仿真 $\dot{\mathbf x} = (A - BK)\mathbf x$。状态空间图最下排显示：角度平滑衰减、小车回到原点、控制力快速降到一个小的稳态推力——整个设计闭环，不到 30 行 Python。

10. 现代图景

本章其实是把 ODE 这门课带着 目标 重演了一遍。第 7-8 章的稳定性理论关心 系统平衡点能否承受小扰动；控制理论关心 如何工程化地塑造 平衡点和收敛速度。第 4 章的拉普拉斯变换不再只是求显式解，而是直接给出频域鲁棒性指标。第 6 章的矩阵指数则成为状态空间仿真和观测器设计的核心。

向前还有：

• 鲁棒控制 ── $H_\infty$、$\mu$ 综合：在有界模型误差下证明性能。
• 自适应控制 / 模型预测控制（MPC） ── 控制器在线更新，每步求解一个滑动时域优化。
• 非线性控制 ── 反馈线性化、滑模控制、控制 Lyapunov 函数。
• 强化学习 ── 控制器从经验中 学 出来，而不是被设计；与第 18 章 Neural ODE 形成自然衔接。

小结

控制理论 把 ODE 从描述提升为设计：我们不再只是预测系统怎么动，而是告诉它 应该 怎么动。

参考资料

• 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).

系列导航