文章

控制论

发表于 更新于
作者 gctoph
56 分钟阅读
控制论

动态系统

控制论要处理的对象：动态系统。

未来由现在、动作和扰动共同决定

\[x_{t+1} = f(x_t, u_t, w_t)\]

其中：

  • $x_t$ 是系统状态；
  • $u_t$ 是控制动作；
  • $w_t$ 是外部扰动；
  • $x_{t+1}$ 是下一时刻状态。

观测不是状态

如果我们不能直接看到状态，只能看到观测值，则有：

\[y_t = g(x_t, v_t)\]

其中：

  • $y_t$ 是观测；
  • $v_t$ 是观测噪声。

开环与闭环

例子：加热一个房间。

开环控制

可以规定：

每天晚上 8 点到 10 点，暖气固定开到 70%。

系统不管房间现在是冷是热。 今天外面 5°C，暖气开 70%。 明天外面 -10°C，暖气还是开 70%。 房间里已经有很多人，暖气仍然开 70%。

开环控制简单、可预测、不容易被测量噪声扰动。 但没有纠错。

闭环控制

使用恒温器：

1
2
3

目标温度：24°C
如果当前温度低于 24°C，就加热
如果当前温度高于 24°C，就停止加热

结构：

1
2
3

目标温度 → 误差 → 控制器 → 暖气 → 房间温度
    ↑                               |
    └────────── 温度计反馈 ──────────┘

闭环系统可以抵抗扰动。

  • 外面突然变冷，房间温度下降，控制器会自动加热。
  • 太阳出来，房间升温，控制器会减少加热。

反馈让系统有了自我修正能力。

但反馈也带来危险。温度计有延迟，暖气有惯性。

如果控制器设计不合理，房间温度可能在 21°C 和 27°C 之间来回震荡，而非稳定在目标 24°C。

比例控制：为什么“纠错太猛”会出错

考虑最简单的离散系统：

\[y_{t+1} = y_t + u_t\]

目标是把 $y_t$ 控制到 0。

定义误差：

\[e_t = 0 - y_t = -y_t\]

比例控制器是：

\[u_t = K_p e_t = -K_p y_t\]

代入系统：

\[y_{t+1} = y_t - K_p y_t\]

所以：

\[y_{t+1} = (1-K_p)y_t\]

要让 $y_t \to 0$，必须有：

\[|1-K_p| < 1\]

因此：

\[0 < K_p < 2\]

控制不是越强越好。增益太大，系统会振荡，甚至发散。

$K_p$闭环系数 $1-K_p$行为
0.20.8慢慢收敛
0.80.2快速收敛
1.00一步到位
1.5-0.5振荡收敛
2.0-1等幅振荡
2.5-1.5振荡发散

开车例子。车偏左一点，往右打方向。

  • 打得太小，车回不来；
  • 打得合适，车回到车道中央；
  • 打得太猛，车冲到右边，然后又往左猛打，车辆开始摆动。

控制器的任务：让闭环系统稳定地消除误差。

延迟

在控制系统中，我们通常认为负反馈能带来稳定。然而，当系统中存在时间延迟（Delay）时，同样的负反馈机制却可能引发持续的振荡甚至发散。

有一步延迟的系统

现实中，传感器读数或执行器动作往往有延迟。假设引入一步延迟，控制器只能看到上一时刻的状态 $y_{t-1}$：

\[u_t = -K y_{t-1}\]

此时，系统的动态方程变为二阶差分方程：

\[y_{t+1} = y_t - K y_{t-1}\]

取增益 $K=1$，并设定初始值 $y_0 = 1, y_1 = 1$。代入方程 $y_{t+1} = y_t - y_{t-1}$ 进行迭代：

时刻 $t$计算过程状态 $y_t$
0(初始值)$1$
1(初始值)$1$
2$y_2 = y_1 - y_0 = 1 - 1$$0$
3$y_3 = y_2 - y_1 = 0 - 1$$-1$
4$y_4 = y_3 - y_2 = -1 - 0$$-1$
5$y_5 = y_4 - y_3 = -1 - (-1)$$0$
6$y_6 = y_5 - y_4 = 0 - (-1)$$1$
7$y_7 = y_6 - y_5 = 1 - 0$$1$

观察： 系统并没有收敛到 0，而是陷入了周期为 6 的持续震荡。

生活中的类比：洗澡调水温

想象你在调节淋浴水温：

  1. 水太冷，你拧大热水阀门。
  2. 延迟存在：热水需要几秒才能流到喷头，你暂时感觉不到变化。
  3. 你以为没效果，于是继续拧大。
  4. 几秒后，滚烫的热水突然涌出。
  5. 你惊慌地拧大冷水阀门。
  6. 又过了几秒，冰水袭来……

延迟让控制器基于过时的信息行动。原本的“纠错动作”因为滞后，变成了新的“错误来源”。

状态空间

将系统状态扩展为二维向量： \(x_t = \begin{bmatrix} y_t \\ y_{t-1} \end{bmatrix}\)

根据动态方程 $y_{t+1} = y_t - K y_{t-1}$，下一时刻的状态为： \(x_{t+1} = \begin{bmatrix} y_{t+1} \\ y_t \end{bmatrix} = \begin{bmatrix} y_t - K y_{t-1} \\ y_t \end{bmatrix}\)

这可以写成矩阵形式 $x_{t+1} = A_K x_t$，其中系统矩阵 $A_K$ 为： \(A_K = \begin{bmatrix} 1 & -K \\ 1 & 0 \end{bmatrix}\)

系统的稳定性完全由矩阵 $A_K$ 的特征值 $\lambda$ 决定：若所有 $|\lambda| < 1$，系统收敛（稳定）。

案例分析：当 K=1 时

系统矩阵为： \(A = \begin{bmatrix} 1 & -1 \\ 1 & 0 \end{bmatrix}\)

求解特征方程 $\det(\lambda I - A) = 0$：

\[\det \begin{bmatrix} \lambda - 1 & 1 \\ -1 & \lambda \end{bmatrix} = \lambda(\lambda - 1) + 1 = \lambda^2 - \lambda + 1 = 0\]

解得特征值： \(\lambda = \frac{1 \pm i\sqrt{3}}{2}\)

计算模长： \(|\lambda| = \sqrt{\left(\frac{1}{2}\right)^2 + \left(\frac{\sqrt{3}}{2}\right)^2} = \sqrt{\frac{1}{4} + \frac{3}{4}} = 1\)

解释： 由于特征值的模长恰好为 1，且为复数，系统在相平面上做旋转运动而不衰减。这解释了前面观察到的周期性振荡序列。

稳定条件

对于增益 $K$，特征方程为：

\[\lambda^2 - \lambda + K = 0\]

根据离散系统的 Jury 稳定性判据（或直接求解），要使所有特征值落在单位圆内 ($|\lambda| < 1$)，增益 $K$ 必须满足：

\[0 < K < 1\]

对比

系统类型动态方程稳定条件说明
无延迟$y_{t+1} = (1-K)y_t$$0 < K < 2$允许较大的增益
一步延迟$y_{t+1} = y_t - K y_{t-1}$$\mathbf{0 < K < 1}$稳定区间减小

延迟降低了系统能容忍的反馈增益区间。同一个 $K$ 值，在无延迟系统中可能非常稳定，但在有延迟系统中可能导致振荡 ($K=1$) 甚至发散 ($K>1$)。

基础模型：标准二阶闭环系统

弹簧-质量-阻尼系统

考虑一个经典的机械振动系统：

  • $m$：质量
  • $c$：阻尼系数
  • $k$：弹簧刚度
  • $x(t)$：位移
  • $F(t)$：外力输入

根据牛顿第二定律： \(m\ddot{x} + c\dot{x} + kx = F(t)\)

对两边进行拉普拉斯变换（假设初始条件为零）： \(ms^2X(s) + csX(s) + kX(s) = F(s)\)

整理得： \(X(s)(ms^2 + cs + k) = F(s)\)

传递函数为： \(G(s) = \frac{X(s)}{F(s)} = \frac{1}{ms^2 + cs + k}\)

将分母除以 $m$： \(G(s) = \frac{1/m}{s^2 + \frac{c}{m}s + \frac{k}{m}}\)

定义：

  • 自然频率：$\omega_n = \sqrt{\frac{k}{m}}$
  • 阻尼比：$\zeta = \frac{c}{2\sqrt{km}} = \frac{c}{2m\omega_n}$

代入后得到： \(G(s) = \frac{1/m}{s^2 + 2\zeta\omega_n s + \omega_n^2}\)

为了使稳态增益为1（单位阶跃响应最终达到1），分子也乘以 $\omega_n^2$： \(T(s) = \frac{\omega_n^2}{s^2 + 2\zeta\omega_n s + \omega_n^2}\)

传递函数

标准二阶系统的传递函数描述了输入 $R(s)$ 到输出 $Y(s)$ 的动态关系：

\[T(s) = \frac{Y(s)}{R(s)} = \frac{\omega_n^2}{s^2 + 2\zeta\omega_n s + \omega_n^2}\]

它描述了一个具有惯性和阻尼的对象，如何响应外部指令。

参数符号名称物理含义直观理解
自然频率$\omega_n$Natural Frequency系统固有的“运动速度”$\omega_n$ 越大，系统反应越快。
(如：弹簧越硬，振动越快)
阻尼比$\zeta$Damping Ratio系统的“刹车力度”$\zeta$ 决定振荡程度。
(如：减震器越强，晃动越少)

极点向量（令分母为0）：

\[s_{1,2} = -\zeta\omega_n \pm j\omega_n\sqrt{1-\zeta^2}\]

在复平面上，极点向量与负实轴的夹角为 $\theta$。

\[\cos \theta = \frac{\text{邻边}}{\text{斜边}} = \frac{\zeta\omega_n}{\omega_n} = \zeta\] \[s_{1,2} = -\cos \theta \omega_n \pm j \sin \theta \omega_n\]

阻尼的四种状态

阻尼描述系统回到目标时是否“刹得住”。

  • $\zeta = 0$ (无阻尼)：永远震荡。极点在虚轴上 ($\pm j\omega_n$)。
  • $0 < \zeta < 1$ (欠阻尼)：会震荡、最终收敛稳定。极点在左半平面，有虚部。
  • $\zeta \approx 0.7$ (欠阻尼，适中)：小幅震荡、收敛较快。
  • $\zeta = 1$ (临界阻尼，最佳)：最快回到原位且不晃动。极点是两个重合的实根。
  • $\zeta > 1$ (过阻尼)：缓慢收敛。极点是两个不同的实根。

幅频特性

在频域中，我们考察系统对正弦输入 $r(t)=\sin(\omega t)$ 的响应。 输出仍为同频率正弦波，但幅值 $A(\omega)$ 和相位 $\phi(\omega)$ 发生改变：

\[y(t) = A(\omega)\sin(\omega t + \phi(\omega))\]

幅频曲线即绘制 $A(\omega) = |T(j\omega)|$：

\[|T(j\omega)| = \frac{\omega_n^2}{\sqrt{(\omega_n^2 - \omega^2)^2 + (2\zeta\omega_n\omega)^2}}\]

（通常再转换为 dB 表示）

频段条件幅值表现物理含义
低频区$\omega \ll \omega_n$$\approx 1$
$(0 \text{ dB})$		完全跟随
系统有足够时间响应慢变化。
谐振区$\omega \approx \omega_n$峰值
(若 $\zeta$ 小)		共振放大
输入节奏踩中系统固有节奏，能量较大。
高频区$\omega \gg \omega_n$$\to 0$
(急剧衰减)		滤波抑制
惯性导致系统跟不上快速变化。

标准二阶系统是一个低通滤波器。

带宽 ($\omega_{bw}$)

定义：幅值下降到 $-3\text{ dB}$ (即原幅值的 $0.707$ 倍) 时的频率。

\[|T(j\omega_{bw})| = \frac{1}{\sqrt{2}} \approx -3\text{ dB}\]
  • 带宽内：信号基本能无损通过。
  • 带宽外：信号被显著削弱。

相频特性：延迟在频域的体现

幅频曲线回答的是输出变大了还是变小了，相频曲线回答的是输出比输入晚了多少拍？

对于线性稳定系统，若输入为 $r(t)=\sin(\omega t)$，输出为：

\[y(t) = A(\omega)\sin(\omega t + \phi(\omega))\]
  • $\phi(\omega) = \angle T(j\omega)$：即相频特性。
  • $\phi < 0$：相位滞后（Lag）。输出慢于输入（绝大多数物理系统）。
  • $\phi > 0$：相位超前（Lead）。输出快于输入（较少见，通常需特殊控制器产生）。

相位角 $\phi$（弧度或角度）可以直接换算为时间延迟 $\tau$。

\[\tau = \frac{|\phi|}{\omega} \quad (\phi \text{ 为弧度})\]

或者使用角度：

\[\tau = \frac{|\phi_{deg}|}{360^\circ} \times T = \frac{|\phi_{deg}|}{360^\circ} \times \frac{2\pi}{\omega}\]
  • $90^\circ$ 滞后 $\rightarrow$ 晚了 1/4 个周期。
  • $180^\circ$ 滞后 $\rightarrow$ 晚了 半个周期（完全反相）。

例子： 若 $\omega = 10 \text{ rad/s}$，相位滞后 $90^\circ$。 周期 $T \approx 0.628s$，则时间延迟 $\tau = 0.628 / 4 \approx 0.157s$。 意味着：输入达到峰值后，过了 0.157秒，输出才达到峰值。

标准二阶系统的相频轨迹

对于标准二阶系统 $T(s) = \frac{\omega_n^2}{s^2 + 2\zeta\omega_n s + \omega_n^2}$，其相位随频率变化呈现典型的“三段式”特征：

\[\angle T(j\omega) = -\arctan_2(2\zeta\omega_n\omega, \ \omega_n^2 - \omega^2)\]
频率区域条件相位值 $\phi$物理含义
低频区$\omega \ll \omega_n$$\approx 0^\circ$同步跟随
输入信号变化慢，系统完全跟得上，几乎无延迟。
谐振区$\omega = \omega_n$$= -90^\circ$滞后1/4拍
系统处于最敏感状态，能量最大，但响应已明显滞后。
高频区$\omega \gg \omega_n$$\to -180^\circ$严重滞后/反相
惯性太大，系统完全跟不上，输出趋势与输入相反。

纯延迟环节的影响

实际系统中常存在纯时间延迟 $e^{-sT_d}$（如传输延迟、计算耗时）。

  • 幅频特性：$|e^{-j\omega T_d}| = 1$ （不改变幅值）
  • 相频特性：$\angle e^{-j\omega T_d} = -\omega T_d$ （线性增加滞后）

同样的时间延迟，对高频信号的危害远大于低频信号。

例子：延迟 $T_d=0.05s$

  • $\omega=1$: 滞后 $\approx 2.9^\circ$ (忽略不计)
  • $\omega=20$: 滞后 $\approx 57.3^\circ$ (显著影响稳定性)

相位滞后会导致不稳定

在负反馈控制系统中，我们依靠“误差信号”来纠正偏差。

负反馈变正反馈的过程：

  1. 正常情况：检测到误差 $\rightarrow$ 施加纠正力 $\rightarrow$ 误差减小。
  2. 相位滞后过大 ($\approx -180^\circ$)：
    • 检测到误差时，系统施加纠正力。
    • 但由于严重滞后，当纠正力真正起作用时，误差的方向已经变了。
    • 原本的“纠正力”变成了“助推力”，反而放大了误差。
    • 结果：负反馈失效，系统变成正反馈，引发振荡甚至发散。

工程上最危险的状况是：

在增益仍较大（$>1$ 或 $0\text{ dB}$）的频率处，相位接近 $-180^\circ$。

在看 Bode 图的相频部分时，关注三点：

  1. 读滞后量（相位裕度）：
    • 在特定频率 $\omega_0$（幅值增益为0dB对应频率）下，$\phi(\omega_0)$ 距离$-180^\circ$是多少？
    • 例：$\phi=-60^\circ$ 表示该频率成分滞后 60 度。幅值裕度 120 度。

  2. 算时间差：判断机械结构或通信延迟是否在可接受范围内。

  3. 看稳定性边界（幅值裕度）：
    • 观察相位曲线穿过 $-180^\circ$ 时的频率（相位穿越频率）。
    • 检查此时的幅值是否已经衰减到足够小（$<-3\text{ dB}$ 或更低）。
    • 如果幅值还很大而相位已达 $-180^\circ$，系统必振荡。

PID（比例 积分 微分）

比例控制只看当前误差：

\[u_t = K_p e_t\]

但真实系统中，只看当前误差不够。系统可能有长期偏差，可能有惯性，可能正在快速恶化。

PID 把控制动作拆成三项：

连续系统：

\[u(t) = K_p e(t) + K_i \int_0^t e(\tau) d\tau + K_d \frac{de(t)}{dt}\]

离散系统：

\[u_t = K_p e_t + K_i \sum_{k=0}^{t} e_k \Delta t + K_d \frac{e_k - e_{k-1}}{\Delta t}\]

（$\Delta t$ 可以 合并进系数）

\[u_t = K_p e_t + K_i' \sum_{k=0}^{t} e_k + K_d' (e_k - e_{k-1})\]

传递函数形式（频域分析用）在拉普拉斯变换下，PID 控制器的传递函数 $G_c(s)$ 为：

\[G_c(s) = \frac{U(s)}{E(s)} = K_p + \frac{K_i}{s} + K_d s\]

PID 的英文全称是：Proportional-Integral-Derivative

通常称为 PID Controller（比例-积分-微分控制器）。

三个组成部分分别对应：

  • P - Proportional (比例)：响应当前误差。
  • I - Integral (积分)：累积过去误差，消除稳态误差。
  • D - Derivative (微分)：预测未来误差趋势，抑制超调。

P 项：当前误差

P 项解决“现在偏离目标多少”。

温度比目标低 2°C，就加热。 车偏左，就向右打方向。 预算花慢了，就提高投放强度。

I 项：长期偏差

有些系统只用 P 控制会留下稳态误差。

例如房间有一扇窗一直漏风。P 控制会加热，但可能最终停在 22°C，而目标是 24°C。

你已经低于目标很久了。即使当前误差不大，我也要继续累积补偿。

I 项危险点：系统可能严重超调。积分饱和（integral windup）

D 项：变化趋势

D 项看误差变化速度。

如果温度虽然只低 1°C，但仍然在下降，控制器应提前加热。

但 D 项对噪声很敏感。因为计算差分：

\[e_t - e_{t-1}\]

如果观测值本身有噪声，差分会放大噪声。

实际工程中，D 项常常需要滤波，甚至不用 D 项。

PID 设计例子

考虑一个质量为 $m$ 的小车，受到以下三种力的作用：

  1. 弹簧回复力：$-kx$
  2. 阻尼摩擦力：$-b\dot{x}$
  3. 电机输入力：$u$

其中：

  • $x(t)$：小车位置
  • $\dot{x}(t)$：小车速度
  • $u(t)$：控制输入
  • $r(t)$：目标位置
  • 跟踪误差定义为：$e(t) = r(t) - x(t)$

根据牛顿第二定律： \(m\ddot{x} = u - b\dot{x} - kx\) 整理得运动方程： \(m\ddot{x} + b\dot{x} + kx = u\)

取系统参数： \(m=1, \quad b=2, \quad k=20\)

代入运动方程： \(\ddot{x} + 2\dot{x} + 20x = u\)

在零初始条件下进行拉普拉斯变换：

\[s^2X(s) + 2sX(s) + 20X(s) = U(s)\] \[(s^2 + 2s + 20)X(s) = U(s)\]

被控对象的传递函数 $G(s)$ 为：

\[G(s) = \frac{X(s)}{U(s)} = \frac{1}{s^2 + 2s + 20}\]

比例控制 (P) 的局限性

若仅采用比例控制 $u = K_p e = K_p(r-x)$，闭环传递函数为： \(\frac{X(s)}{R(s)} = \frac{K_p G(s)}{1 + K_p G(s)} = \frac{K_p}{s^2 + 2s + 20 + K_p}\)

闭环特征方程为： \(s^2 + 2s + (20 + K_p) = 0\)

对比标准二阶系统 $s^2 + 2\zeta\omega_n s + \omega_n^2 = 0$，可得： \(\omega_n = \sqrt{20 + K_p}, \quad \zeta = \frac{1}{\sqrt{20 + K_p}}\)

分析结论：

  • 增大 $K_p$ $\Rightarrow$ $\omega_n$ 增大（响应变快）。
  • 增大 $K_p$ $\Rightarrow$ $\zeta$ 减小（阻尼降低，振荡增强）。

示例： 取 $K_p = 44$，则 $\omega_n = 8$，$\zeta = 0.125$。系统响应极快但严重欠阻尼。

核心问题：比例控制只能改变刚度项，无法独立调节阻尼项，因此无法同时实现“快速”与“稳定”。

PID 控制器

采用 PID 控制器： \(C(s) = K_p + \frac{K_i}{s} + K_d s = \frac{K_d s^2 + K_p s + K_i}{s}\)

单位负反馈系统的闭环特征方程由 $1 + C(s)G(s) = 0$ 给出： \(1 + \frac{K_d s^2 + K_p s + K_i}{s(s^2 + 2s + 20)} = 0\)

\[s^3 + (2 + K_d)s^2 + (20 + K_p)s + K_i = 0\]

性能指标要求

  • 最大超调量 $M_p < 10\%$
  • 调节时间 $T_s \approx 1.5 \text{ s}$
  • 阶跃响应无稳态误差

极点配置法

选择主导二阶极点参数： \(\zeta = 0.7, \quad \omega_n = 4\)

对应的二阶多项式为： \(s^2 + 2\zeta\omega_n s + \omega_n^2 = s^2 + 5.6s + 16\)

由于引入积分项后系统变为三阶，需选择第三个极点。为确保其不影响主导动态，将其置于远离虚轴的位置，取 $s = -8$，对应因子 $(s+8)$。

期望闭环特征多项式为： \((s^2 + 5.6s + 16)(s + 8)\)

展开计算： \(= s^3 + 8s^2 + 5.6s^2 + 44.8s + 16s + 128\) \(= s^3 + 13.6s^2 + 60.8s + 128\)

将实际特征多项式与期望特征多项式系数逐项匹配：

系数项实际多项式系数期望多项式系数方程解得参数
$s^2$$2 + K_d$$13.6$$2 + K_d = 13.6$$K_d = 11.6$
$s^1$$20 + K_p$$60.8$$20 + K_p = 60.8$$K_p = 40.8$
$s^0$$K_i$$128$$K_i = 128$$K_i = 128$

最终控制器 \(C(s) = 40.8 + \frac{128}{s} + 11.6s\)

闭环传递函数

\(T(s) = \frac{C(s)G(s)}{1 + C(s)G(s)} = \frac{11.6s^2 + 40.8s + 128}{s^3 + 13.6s^2 + 60.8s + 128}\)

稳态误差分析

对于单位阶跃输入 $R(s) = 1/s$，利用终值定理：

\[e_{\infty} = \lim_{s \to 0} s E(s) = \lim_{s \to 0} \frac{s R(s)}{1 + C(s)G(s)} = \lim_{s \to 0} \frac{1}{1 + C(s)G(s)}\]

由于 $C(s)$ 包含积分项 $\frac{K_i}{s}$，当 $s \to 0$ 时，$C(s)G(s) \to \infty$（假设 $K_i > 0$ 且 $G(0)$ 有限）。

\[\therefore e_{\infty} = 0\]

结论：积分项保证了阶跃目标下的零稳态误差。

时域性能估计 (基于主导极点)

基于 $\zeta=0.7, \omega_n=4$ 的二阶近似：

  1. 最大超调量 ($M_p$): \(M_p = e^{-\frac{\zeta\pi}{\sqrt{1-\zeta^2}}} = e^{-\frac{0.7\pi}{\sqrt{0.51}}} \approx 0.046 \quad (4.6\%)\) (满足 $<10\%$ 的要求)

  2. 峰值时间 ($T_p$): \(T_p = \frac{\pi}{\omega_n\sqrt{1-\zeta^2}} = \frac{\pi}{4\sqrt{0.51}} \approx 1.10 \text{ s}\)

  3. 调节时间 ($T_s$, 2%准则): \(T_s \approx \frac{4}{\zeta\omega_n} = \frac{4}{0.7 \times 4} \approx 1.43 \text{ s}\) (与设计目标 $1.5 \text{ s}$ 一致)

注意：实际系统含有第三极点和零点，真实响应会与上述二阶估计略有偏差，但主导动态尺度相符。

频域性能与延迟鲁棒性

开环传递函数 \(L(s) = C(s)G(s) = \frac{11.6s^2 + 40.8s + 128}{s(s^2 + 2s + 20)}\)

基础频域指标：

  • 增益交叉频率：$\omega_c \approx 12.67 \text{ rad/s}$
  • 相位裕度 (PM)：$PM \approx 83.6^\circ$

较大的相位裕度表明系统在无延迟情况下具有良好的鲁棒性。

纯延迟对稳定性的影响

若系统存在纯延迟 $e^{-sT_d}$，开环传递函数变为 $L_d(s) = L(s)e^{-sT_d}$。

  • 相位特性：引入额外滞后 $\phi_d = -\omega T_d \frac{180^\circ}{\pi}$。

在交叉频率 $\omega_c$ 处，剩余相位裕度为： \(PM_{\text{new}} = PM - \omega_c T_d \frac{180^\circ}{\pi}\)

不同延迟下的稳定性评估：

延迟 $T_d$相位损失 $\Delta \phi$剩余相位裕度 $PM_{\text{new}}$状态判断
$0 \text{ ms}$$0^\circ$$83.6^\circ$稳健
$20 \text{ ms}$$14.5^\circ$$69.1^\circ$稳健
$50 \text{ ms}$$36.3^\circ$$47.3^\circ$可用，裕度下降
$100 \text{ ms}$$72.6^\circ$$11.0^\circ$接近振荡边界

令 $PM_{\text{new}} = 0$，估算理论稳定边界：

\[T_{d,\max} = \frac{PM_{\text{rad}}}{\omega_c} = \frac{83.6 \times \frac{\pi}{180}}{12.67} \approx \frac{1.46}{12.67} \approx 0.115 \text{ s}\]

工程建议：此为线性近似下的极限值，实际工程中可接受延迟应显著小于此值以保留安全裕度。

当延迟不可忽略时，高带宽控制器会因“基于过时信息决策”导致振荡。常用处理方法：

  1. 降低增益：减小 $K_p, K_i, K_d$ 以降低交叉频率 $\omega_c$，从而减小 $\omega_c T_d$ 带来的相位损失。代价是响应变慢。
  2. 降低目标带宽：重新设计极点，例如将 $\omega_n$ 从 4 降至 2。牺牲速度以换取对延迟和高频噪声不敏感。
  3. 引入预测控制：如 Smith Predictor。利用模型预测未来状态 $\hat{x}(t+T_d)$，基于预测误差进行控制，补偿物理延迟带来的相位滞后。

状态空间

控制论中更通用的写法是状态空间模型：

\[x_{t+1} = A x_t + B u_t\] \[y_t = C x_t + D u_t\]

其中：

  • $x_t$ 是系统在时刻 $t$ 的状态；
  • $u_t$ 是控制动作；
  • $y_t$ 是我们能观测到的输出；
  • $A$ 描述系统在没有控制输入时如何自己演化；
  • $B$ 描述控制输入如何影响状态；
  • $C$ 描述状态如何映射成观测；
  • $D$ 描述控制输入是否会直接影响输出。

状态反馈

有了状态空间模型之后，一个自然的控制策略是状态反馈：

\[u_t = -Kx_t\]

也就是说，控制动作由当前状态决定。

代入系统方程：

\[x_{t+1} = A x_t + B u_t\]

得到：

\[x_{t+1} = A x_t - B K x_t\]

因此闭环系统变成：

\[x_{t+1} = (A - BK)x_t\]

这说明反馈矩阵 $K$ 会改变闭环系统的动力学。

原来的系统由 $A$ 决定； 加入反馈之后，系统由 $A-BK$ 决定。

所以设计控制器，本质上是在设计闭环系统的动态行为：让它稳定、足够快、不过度振荡、控制动作不过大。

LQR：最简单的最优状态反馈控制

PID 也是一种反馈控制，但它更多依赖工程经验调参。LQR 则把这个问题写成一个明确的优化问题。

LQR 的全称是 Linear Quadratic Regulator，即线性二次型调节器。

它处理的问题是：

在线性系统中，如何选择控制动作，使系统既尽快回到目标，又不要使用过大的控制力？

系统模型为：

\[x_{t+1} = A x_t + B u_t\]

LQR 不直接规定“必须多快回到目标”，而是定义一个长期成本函数：

\[J = \sum_{t=0}^{\infty} \left( x_t^\top Q x_t + u_t^\top R u_t \right)\]

\(x_t^\top Q x_t\) 惩罚状态偏差。它表示：系统偏离目标越远，代价越大。

\(u_t^\top R u_t\) 惩罚控制动作。它表示：控制动作越大，代价越大。

在两件事之间做权衡：状态误差要小，控制动作也不能太大。

例如：

  • 小车要尽快回到目标位置，但电机不能输出无限大的力；
  • 无人机要快速稳定姿态，但不能让电机剧烈抖动；
  • 云服务 autoscaling 要降低排队延迟，但不能频繁扩缩容。

动作本身也有成本。一个控制器不仅要减少误差，也要避免过度控制。

LQR 为什么会得到线性反馈？

LQR 的一个漂亮结果是：虽然我们优化的是整个未来的控制序列，但最终得到的最优控制律非常简单：

\[u_t = -Kx_t\]

也就是说，最优策略仍然是一个线性状态反馈。其中 $K$ 由系统矩阵 $A, B$ 和权重矩阵 $Q, R$ 决定。

系统模型：\(x_{t+1} = Ax_t + Bu_t\)

目标函数：最小化无限时域下的累积代价：

\[J = \sum_{t=0}^{\infty} \left( x_t^\top Q x_t + u_t^\top R u_t \right)\]

因为系统是线性的，代价是二次型，线性变换后的二次型仍是二次型。因此可以设从状态 $x$ 出发的最优未来代价为： \(V(x) = x^\top P x\)

Bellman 方程为： \(V(x) = \min_u \left[ x^\top Q x + u^\top R u + V(Ax + Bu) \right]\)

代入 $V(x) = x^\top P x$： \(V(x) = \min_u \left[ x^\top Q x + u^\top R u + (Ax + Bu)^\top P (Ax + Bu) \right]\)

展开后，只看与 $u$ 有关的项： \(u^\top (R + B^\top P B) u + 2u^\top B^\top P A x\)

对 $u$ 求导，并令其为零： \((R + B^\top P B)u + B^\top P A x = 0\)

解得最优控制： \(u^\star = -(R + B^\top P B)^{-1} B^\top P A x\)

记： \(K = (R + B^\top P B)^{-1} B^\top P A\)

于是得到线性反馈律： \(u^\star = -Kx\)

离散代数 Riccati 方程 (DARE)

矩阵 $P$ 满足离散代数 Riccati 方程：

\[P = Q + A^\top P A - A^\top P B (R + B^\top P B)^{-1} B^\top P A\]

LQR 的计算分两步：

  • 第一步：解 Riccati 方程得到 P；
  • 第二步：由 P 计算反馈增益 K

因此，LQR 不是简单地按当前误差反馈，而是在考虑未来演化后的最优线性反馈。

LQR 的局限

LQR 很优雅，但它也有明确限制。

第一，LQR 假设系统是线性的：

\[x_{t+1}=Ax_t+Bu_t\]

如果系统强非线性，LQR 只能作为局部近似，或者需要使用 iLQR、MPC 等扩展方法。

第二，LQR 的标准形式不直接处理硬约束。

例如：

\[u_{\min} \le u_t \le u_{\max}\] \[x_{\min} \le x_t \le x_{\max}\]

在真实系统中，控制量往往有上限，状态也有安全边界。标准 LQR 会惩罚大动作，但不会天然保证动作不越界。

第三，LQR 需要状态 $x_t$。如果状态不能直接观测，就需要先做状态估计，例如使用 Kalman Filter。

因此在工程架构中，LQR 通常位于这样的层次中：

1
2
3
4
5
6
7
8
9

Raw Observations
↓
State Estimation
↓
LQR / PID / MPC
↓
Guardrails
↓
Execution

LQR 负责给出一个平滑、合理的状态反馈控制动作；但安全边界、输入饱和、业务规则仍然需要额外处理。

MPC：模型预测控制 (Model Predictive Control)

硬约束

LQR 虽然优美，但它不直接处理硬约束。在真实系统中，约束往往至关重要：

  • 执行器限制：小车的电机输出有最大力矩；
  • 安全范围：无人机的倾角和推力都有物理极限；
  • 环境边界：自动驾驶车辆不能越过车道边界；
  • 业务逻辑：云服务 autoscaling 不能一次性扩容过多，也不能低于最小实例数。

此时，仅靠 LQR 中的 $R$ 矩阵惩罚大动作是不够的（那是软约束）。我们需要在优化问题中显式加入硬约束。

MPC 的核心流程

MPC 也被称为滚动时域控制 (Receding Horizon Control)。其基本流程如下：

1
2
3
4
5

1. 估计当前状态 x_t
2. 基于模型，向前预测未来 H 步
3. 求解一个有限时域的约束优化问题，得到最优控制序列
4. 只执行序列中的第一步动作 u_t
5. 下一时刻 t+1，重新估计状态，重复上述过程

形式化

假设系统仍为离散线性系统： \(x_{t+1} = Ax_t + Bu_t\)

在当前时刻 $t$，观测到状态 $x_t$。MPC 向前看 $H$ 步（预测时域/Horizon）。

  • 控制序列：$u_{t|t}, u_{t+1|t}, \dots, u_{t+H-1|t}$
    • $u_{t+k|t}$ 表示：在时刻 $t$ 规划时，预测第 $t+k$ 步使用的控制量。
  • 状态序列：$x_{t|t}, x_{t+1|t}, \dots, x_{t+H|t}$
    • 初始条件：$x_{t|t} = x_t$
    • 动力学约束：$x_{t+k+1|t} = A x_{t+k|t} + B u_{t+k|t}$

MPC 在每个时刻求解以下优化问题：

\[\begin{aligned} \min_{u_{t|t}, \dots, u_{t+H-1|t}} \quad & J = \sum_{k=0}^{H-1} \left( x_{t+k|t}^\top Q x_{t+k|t} + u_{t+k|t}^\top R u_{t+k|t} \right) + x_{t+H|t}^\top P x_{t+H|t} \\ \text{s.t.} \quad & x_{t+k+1|t} = A x_{t+k|t} + B u_{t+k|t}, \quad k=0,\dots,H-1 \\ & u_{\min} \le u_{t+k|t} \le u_{\max} \\ & x_{\min} \le x_{t+k|t} \le x_{\max} \end{aligned}\]

参数含义：

  • $Q$：惩罚状态偏差；
  • $R$：惩罚控制能量；
  • $P$：终端代价矩阵（Terminal Cost），用于近似无限时域后的剩余代价，保证稳定性；
  • $H$：预测时域长度；
  • $u_{\min}/u_{\max}$：输入硬约束；
  • $x_{\min}/x_{\max}$：状态硬约束。

注：如果去掉所有约束，且令 $H \to \infty$，线性二次型 MPC 的解将收敛于 LQR 的解。

示例：带输入上限的小车

考虑一维系统： \(x_{t+1} = x_t + u_t\) 目标：使位置偏差 $x_t \to 0$。

设定：

  • 初始状态：$x_0 = 10$
  • 控制约束：$-2 \le u_t \le 2$
特性LQR (无约束)MPC (有约束)
第一步动作可能计算出 $u_0 = -6.18$受限于约束，计算出 $u_0 = -2$
可行性不可执行 (超出电机能力)可执行 (满足硬件限制)
行为逻辑线性反馈，始终按比例输出远时饱和输出，近时平滑调节

假设预测时域 $H=3$：

  1. 远离目标时 ($x_0=10$)： MPC 发现无论怎么努力，短时间内都无法到达原点，因此会尽可能使用最大允许动作： \(u_{0|0}^\star = -2, \quad u_{1|0}^\star = -2, \quad u_{2|0}^\star = -2\) 执行第一步：$u_0 = -2$，状态变为 $x_1 = 8$。

  2. 接近目标时 ($x_t=1$)： MPC 预测如果使用 $u=-2$，下一步状态将变为 $-1$，产生新的偏差和代价。因此它会选择更温和的动作： \(u_t \approx -1\) 从而平滑地逼近零点，避免过冲。

MPC 的局限性

计算成本高

  • LQR 只需一次矩阵乘法 $u=-Kx$。
  • MPC 需要在每个采样周期内在线求解一个优化问题（通常是二次规划 QP）。
  • 对于高频系统或高维状态空间，实时求解可能成为瓶颈。

依赖模型精度

  • MPC 的核心是“预测”。如果模型 $A, B$ 不准确，预测轨迹将与真实轨迹偏离，导致控制性能下降甚至不稳定。
  • 对策：结合系统辨识、在线参数校准或鲁棒 MPC (Robust MPC)。

参数整定复杂 需要精心设计的参数包括：

  • 预测时域 $H$：太短可能导致短视，太长增加计算量。
  • 权重矩阵 $Q, R, P$：平衡响应速度与能耗。
  • 终端约束/代价：保证闭环稳定性。

可行性问题

  • 如果当前状态已违反约束，或者由于扰动导致未来无论如何控制都无法满足约束，优化问题将无解。
  • 工程对策：软约束 (Soft Constraints) 引入松弛变量 (Slack Variable) $s_k \ge 0$： \(x_{t+k|t} \le x_{\max} + s_k\) 并在目标函数中加入高额惩罚： \(J_{new} = J_{original} + \sum \rho s_k^2\) 这样即使无法完全满足约束，求解器也能找到一个“违反程度最小”的可行解，保证控制器不崩溃。
维度LQRMPC
时域无限时域有限时域 (滚动)
约束无硬约束 (仅软惩罚)显式支持硬约束
计算离线计算 $K$，在线极快在线求解优化问题，计算较重
模型依赖中等高 (预测准确性关键)
适用场景线性、无约束、快速系统有约束、多变量、对安全性要求高的系统

Kalman Filter：从噪声观测到状态估计

在前面的 LQR 和 MPC 讨论中，我们默认控制器可以直接获取完美的系统状态 $x_t$。但在真实系统中，真实状态往往是不可直接观测的。控制器拿到的通常是带噪声的观测值（Noisy Observations/Metrics）

先估计状态 (Estimate State)，再控制状态 (Control State)。

Kalman Filter 处理的是包含两种不确定性的线性高斯系统：

状态演化方程 (Process Model)

\(x_{t+1} = A x_t + B u_t + w_t\)

  • $A$：状态转移矩阵
  • $x_t$：真实状态（隐藏变量）
  • $B$：控制输入矩阵
  • $u_t$：控制动作
  • $w_t$：过程噪声 (Process Noise)，代表模型误差或外部扰动，$w_t \sim \mathcal{N}(0, W)$

观测方程 (Measurement Model)

\(y_t = C x_t + v_t\)

  • $y_t$：观测值（可见变量）
  • $C$：观测矩阵
  • $v_t$：观测噪声 (Measurement Noise)，代表传感器误差，$v_t \sim \mathcal{N}(0, V)$

Kalman Filter 的目标： 根据历史观测序列 ${y_0, \dots, y_t}$ 和控制序列 ${u_0, \dots, u_{t-1}}$，估计当前最可能的真实状态 $\hat{x}_t$。

控制器实际使用的是估计值： \(u_t = \pi(\hat{x}_t) \quad \text{例如 LQR: } u_t = -K\hat{x}_t\)

Kalman Filter不是简单的移动平均。Moving Average 只看历史数据，而 Kalman Filter 结合了：

  1. 系统模型预测 (Model Prediction)
  2. 当前观测值 (Current Measurement)
  3. 模型的不确定性 ($W$)
  4. 观测的不确定性 ($V$)

它在“相信模型”和“相信观测”之间做出最优权衡。

Kalman Filter 的核心算法

每一步包含两个阶段：预测 (Predict) 和 更新 (Update)。

预测步骤

基于上一时刻的估计和控制动作，预测当前状态： \(\hat{x}_{t|t-1} = A \hat{x}_{t-1|t-1} + B u_{t-1}\)

$\hat{x}_{t|t-1}$ 先验估计 (Prior Estimate)，即在看 $y_t$ 之前对 $x_t$ 的猜测。

更新步骤

拿到新观测 $y_t$ 后，修正预测：

  1. 计算新信息 (Innovation)： \(\tilde{y}_t = y_t - C \hat{x}_{t|t-1}\)
    • 含义：观测值与预测值的差异，代表了“新信息”。
  2. 状态更新： \(\hat{x}_{t|t} = \hat{x}_{t|t-1} + L_t (y_t - C \hat{x}_{t|t-1})\)
    • $\hat{x}_{t|t}$：后验估计 (Posterior Estimate)，即最终输出的状态估计。
    • $L_t$：卡尔曼增益 (Kalman Gain)。

$L_t$ 是一个动态调整的权重，决定了我们多相信模型，还是多相信观测。

  • 若观测噪声大 ($V$ 大)：
    • 观测不可靠 $\rightarrow$ 少修正 $\rightarrow$ $L_t$ 较小。
    • 滤波器更依赖模型预测。
  • 若过程噪声大/模型不准 ($W$ 大)：
    • 预测不可靠 $\rightarrow$ 多修正 $\rightarrow$ $L_t$ 较大。
    • 滤波器更依赖当前观测。

示例：估计真实负载

场景：估计服务真实负载 $x_t$，观测值为 noisy QPS $y_t$。

  • 模型：$x_{t+1} = x_t + w_t$ (负载变化缓慢)
  • 观测：$y_t = x_t + v_t$

假设：

  • 模型预测当前负载： $\hat{x}_{t|t-1} = 100$
  • 突然观测到指标飙升：$y_t = 130$
  • 卡尔曼增益：$L_t = 0.3$ (说明观测有一定噪声，不完全可信)

计算： \(\begin{aligned} \hat{x}_{t|t} &= 100 + 0.3 \times (130 - 100) \\ &= 100 + 0.3 \times 30 \\ &= 100 + 9 \\ &= 109 \end{aligned}\)

粒子滤波：Particle Filter

Kalman Filter (KF) 非常优雅且高效，但它有两个核心前提：

  1. 系统近似线性。
  2. 噪声近似高斯分布。

需要一种能表示任意形状概率分布的方法。

Particle Filter (粒子滤波) 解决：当系统非线性、噪声非高斯、状态分布形状复杂时，如何估计隐藏状态？

核心思想：不再维护一个状态估计值，而是维护一群可能的状态假设 (Particles)。从“点估计”到“经验分布”。

Kalman Filter 的表示

\(x_t \sim \mathcal{N}(\hat{x}_t, \Sigma_t)\)

  • 用一个椭圆（均值和协方差）描述不确定性。
  • 假设分布是单峰、对称的高斯分布。

Particle Filter 的表示

\(\{ x_t^{(i)}, w_t^{(i)} \}_{i=1}^N\)

  • $x_t^{(i)}$：第 $i$ 个粒子，代表一个可能的状态假设。
  • $w_t^{(i)}$：第 $i$ 个粒子的权重，代表该假设解释当前观测的可信度 ($\sum w_i = 1$)。

这构成了一个经验分布 (Empirical Distribution)： \(p(x_t | y_{1:t}) \approx \sum_{i=1}^{N} w_t^{(i)} \delta(x_t - x_t^{(i)})\)

  • 如果粒子集中在某处 $\rightarrow$ 高置信度。
  • 如果粒子分成两簇 $\rightarrow$ 存在两种可能的解释（多峰）。
  • 如果粒子分散 $\rightarrow$ 高度不确定。

Particle Filter 不要求线性或高斯假设，通用形式为：

  1. 状态转移方程： \(x_{t+1} = f(x_t, u_t) + w_t\)
  2. 观测方程： \(y_t = h(x_t) + v_t\)

其中 $f(\cdot)$ 和 $h(\cdot)$ 可以是任意非线性函数，$w_t, v_t$ 可以是任意分布的噪声。

算法流程：Predict - Weight - Resample

Particle Filter 每一步迭代包含三个核心步骤：

第一步：预测 (Predict)

对每个粒子，根据系统模型向前推演： \(x_{t|t-1}^{(i)} \sim p(x_t | x_{t-1}^{(i)}, u_{t-1})\)

  • 直观理解：假设上一时刻第 $i$ 个粒子是真的，经过控制 $u_{t-1}$ 和过程噪声后，它现在可能在哪里？
  • 操作：$x_{new} = f(x_{old}, u) + \text{sampled noise}$

第二步：加权 (Weight)

拿到新观测 $y_t$ 后，评估每个粒子与观测的匹配程度： \(w_t^{(i)} \propto p(y_t | x_{t|t-1}^{(i)})\)

  • 似然计算：如果粒子所在的状态能很好地解释观测 $y_t$，则权重增加；否则权重降低。
  • 归一化：确保 $\sum w_t^{(i)} = 1$。

第三步：重采样 (Resample)

经过几轮迭代，大多数粒子的权重会趋近于 0，只有少数粒子权重很大（粒子退化问题）。

  • 操作：根据权重比例，随机复制高权重粒子，丢弃低权重粒子。
  • 目的：将计算资源集中在更可能的状态区域，保持粒子群的多样性及有效性。
  • 触发条件：通常当有效粒子数 $N_{eff} < N/2$ 时触发。

场景 A：机器人定位 (Robot Localization)

场景：机器人有地图，但不知道自己在哪里，需要定位。

  1. 初始化：机器人在未知走廊，在地图撒下一堆均匀分布的粒子。
  2. 移动：机器人向前移动。所有粒子根据运动模型和地图向前移动，并加入随机扰动（模拟轮子打滑）。
  3. 观测：传感器检测到“左墙距离 1m，前墙距离 2m”。
  4. 加权：
    • 位于走廊中间的粒子：预测应看到左右墙距离相等 $\rightarrow$ 与观测不符 $\rightarrow$ 权重降低。
    • 靠近左侧墙壁的粒子：预测符合观测 $\rightarrow$ 权重升高。
  5. 重采样：保留靠近左墙的粒子，复制它们，丢弃其他粒子。
  6. 结果：地图中的粒子群迅速收敛到机器人的真实位置附近。

场景 B：云服务拥塞诊断

观测到 Latency 飙升，可能有多种解释：

  • 粒子群假设：
    • 粒子簇 1：真实负载激增 (High Demand)。
    • 粒子簇 2：下游服务故障 (Downstream Failure)。
    • 粒子簇 3：指标采集异常 (Metric Noise)。
  • 动态更新：
    • 若 QPS 同步上升 $\rightarrow$ 簇 1 权重增加。
    • 若下游 Error Rate 上升 $\rightarrow$ 簇 2 权重增加。
    • 若其他指标正常且很快恢复 $\rightarrow$ 簇 3 权重增加。

控制动作

控制器通常需要确定的状态输入。如何处理粒子群？

  1. 加权平均 (Mean Estimate)： \(\hat{x}_t = \sum_{i=1}^{N} w_t^{(i)} x_t^{(i)}\)
    • 适用于单峰分布，直接送入 LQR/MPC。
  2. 最大后验估计 (MAP)： 选择权重最大的粒子作为状态估计。

局限性与工程挑战

  1. 计算成本高：
    • 需要维护数百甚至数千个粒子。
    • 每个粒子都要运行一次模型预测 $f(\cdot)$ 和观测似然计算 $h(\cdot)$。
  2. 维数灾难 (Curse of Dimensionality)：
    • 随着状态维度增加，覆盖状态空间所需的粒子数呈指数增长。
    • 对策：Rao-Blackwellized Particle Filter (部分状态用 KF 解析求解，部分用 PF)。
  3. 粒子贫化 (Particle Deprivation)：
    • 重采样可能导致粒子多样性丧失，所有粒子聚集在一起，无法应对突发变化。
    • 对策：引入额外噪声 (Jitter)、自适应重采样阈值。
  4. 模型依赖性：
    • 权重计算依赖于准确的观测模型 $p(yx)$。如果模型错误，滤波器会被误导。

其他资料

warning

AI
控制 PID LQR MPC Kalman Filter
本文由作者按照 CC BY 4.0 进行授权