1. 项目起点被“模型”逼出来的无模型思路前阵子同事抱着一堆现场数据来找我说某个非线性对象的预测控制做不动了。对象机理写不出来参数辨识的结果勉强凑合工况一变就崩。他问了我一句很直接的话“能不能只用数据就把预测控制给跑了”我当时就把 MFAC 这套思路给他了只提了两条路一条是 MFAPC一条是 MFAILC。也正是这个需求让我把这套无模型自适应预测控制与迭代学习控制的数值验证仿真程序重新完整梳理了一遍。这篇就聊聊这套仿真程序的设计思路、算法细节、实现过程以及我调试过程中实实在在想记录下来的坑。先给还不熟悉的读者一个定位MFAPCModel-Free Adaptive Predictive Control无模型自适应预测控制和 MFAILCModel-Free Adaptive Iterative Learning Control无模型自适应迭代学习控制同属无模型自适应控制家族核心特征是不需要被控对象的机理数学模型直接基于系统的输入输出数据在线设计控制器。MFAPC 适合当前时刻需要对未来一段输出做优化、系统有非线性或者滞后特征的场景MFAILC 适合系统在固定时间区间上反复执行同一参考轨迹、希望批次间误差不断压缩的场景。仿真程序的目的就是在真实对象无法建模的条件下先通过数值验证把这两种算法的收敛性、参数范围和抗干扰表现摸清楚再决定是否拿到现场用。这个方向适合谁参考我的判断是三类人第一类是正在做控制器设计但被建模问题卡住的工程师第二类是自动化、控制理论方向的研究生需要跑通 MFAPC 和 MFAILC 代码并理解背后的调整逻辑第三类是准备把无模型控制方法落地到实际设备的工艺工程师想先搞清楚仿真阶段要准备哪些测试用例和评价指标。接下来我按自己做项目时的顺序展开。1.1 传统预测控制绕不开的“建模硬伤”工业现场做预测控制最花时间的往往不是控制器本身而是那个被当成“已知条件”的对象模型。你做一个化工反应釜的温度控制反应机理涉及传热、传质、反应动力学参数还随催化剂活性变化你做一个伺服系统的位置控制摩擦力、负载惯量、弹性形变这些因素又混在一起。机理建模要么成本极高要么精度不够辨识得来的参数一旦遇到工况迁移控制器效果立刻掉下来。我见过太多项目死在“建模”这一环节。有人花了一个月做阶跃响应测试得到的传递函数在某个工作点匹配得不错换一个工作点却直接发散。也有人试过用神经网络做非线性辨识离线精度很高在线更新却困难重重。而 MFAC 的思路恰好绕开了这条路不需要机理模型不依赖离线辨识只利用实时采集的输入输出数据在控制器运行过程中同步估计系统的局部动态特征再基于这个局部特征设计控制律。预测控制和迭代学习控制是两条不同的工程技术路线。预测控制强调当前时刻对未来的优化天然擅长处理约束和滞后迭代学习控制强调利用历史批次数据修正当前批次天然擅长重复性任务。把“无模型自适应”这个框架分别与它们结合就得到了 MFAPC 和 MFAILC。这两者的关系不是替代而是互补所以在数值验证时我把它们放在同一个程序框架里跑。1.2 从 MFAC 到 MFAPC、MFAILC 的演化逻辑要理解 MFAPC 和 MFAILC得先理解 MFAC 的地基也就是“动态线性化”。你不需要真的知道系统长什么样你只需要知道在当前这个工作点附近系统的输出变化量和控制输入变化量之间存在一个时变的等价增益关系。这个等价增益就是伪偏导数Pseudo Partial DerivativePPD。MFAC 的基本形式是单步控制每一拍算出控制增量让当前输出逼近参考值。但实际工程里你往往不满足于单步控制系统有滞后你希望控制动作能提前考虑到未来几步的参考轨迹系统有约束你希望控制输入不要猛打猛收系统有大幅值扰动你希望控制动作是平缓的而不是一棒子打过去。预测控制的滚动优化思想正好能补上这块于是就产生了 MFAPC。另一类需求是重复性运动。注塑机的合模过程、点焊机器人的运动轨迹、物料搬运的往复节拍这类过程每个批次开始时的初始条件基本一致参考轨迹也是固定的。传统的反馈控制每一批都是从零开始不会积累经验而迭代学习控制会在每一批结束后把输出误差折算成控制修正量叠加到下一批的控制序列上相当于一台设备在学习自己前一次的错误。把这个思路与伪偏导数的在线估计结合就是 MFAILC。在我的仿真程序里MFAPC 被设计为“时间轴上滚动运行”的控制器MFAILC 被设计为“时间轴内层循环 迭代轴外层循环”的控制器。两者共用同一个被控对象共用同一套数据采集逻辑只是控制器内部结构不同这样就方便对比同一场景下两种策略的差异。1.3 数值仿真要回答的三个核心问题我在设计这套仿真程序的时候并没有一上来就堆代码而是先定义清楚仿真的目标。对我来说数值验证至少要回答三个问题。第一个问题是收敛性给定一个非线性、强耦合、甚至带噪声的被控对象MFAPC 的输出能不能稳定跟踪参考轨迹MFAILC 的批次误差能不能随着迭代次数单调下降。收敛性不是靠理论推导能完全预判的仿真能给出直观的曲线证据。第二个问题是参数敏感性伪偏导数的初始值、步长因子、惩罚因子、预测时域这些参数到底怎么取它们之间是怎么互相影响的哪些参数稍微动一下就发散哪些参数怎么调都不太敏感这些信息只有亲自动手仿真才能积累出直觉。第三个问题是抗扰性现场数据不可能干净测量噪声、对象突变、前批次遗留扰动都会影响控制效果。仿真时我会有意加入不同强度的噪声和时变扰动观察算法失效的边界在哪里。把这三个问题跑清楚一个仿真程序才算真正完成了数值验证的任务而不是仅仅画了一条跟踪曲线。2. 算法核心细节拆解从原理公式到可编程实现仿真程序的基础是算法本身。你如果只是照着论文抄公式跑出来的结果往往是差强人意你如果理解了每个公式背后的代价函数和约束条件才知道参数该怎么调发散时该往哪个方向改。这一节我把 MFAPC 和 MFAILC 的核心公式拆开讲尽量落到“可以直接写代码”的程度。2.1 紧格式动态线性化伪偏导数 PPD 的含义考虑一个单输入单输出的非线性离散系统机理未知我们只知道它每个采样周期的输入输出。MFAC 的紧格式动态线性化Compact Form Dynamic LinearizationCFDL给出这样一个结论在控制输入变化量 Δu(k) 不为零的前提下系统输出增量可以被表示成如下等价形式[ y(k1) - y(k) \phi_c(k) \cdot \Delta u(k) ]其中 φ_c(k) 就是伪偏导数它是一个时变的、有界的标量。这句话的工程含义是我不需要知道系统内部的非线性函数具体是什么只需要在当前工作点附近用一个等效局部增益关系把输入变化和输出变化联系起来。你可以把这个 PPD 理解成“系统在当前工作点上的局部增益”的在线估计值就像路面上某一点的坡度走一步测一步坡度随位置变化但你不需要整条路的地形图。PPD 的在线估计采用投影算法代价函数是估计误差平方加上估计值变化量的惩罚[ J(\phi_c(k)) \left| y(k) - y(k-1) - \phi_c(k)\Delta u(k-1) \right|^2 \mu \left| \phi_c(k) - \hat{\phi}_c(k-1) \right|^2 ]对 φ_c(k) 求极值得到估计更新公式[ \hat{\phi}_c(k) \hat{\phi}_c(k-1) \frac{\eta \Delta u(k-1) \left( \Delta y(k) - \hat{\phi}_c(k-1)\Delta u(k-1) \right)}{\mu \Delta u(k-1)^2} ]公式里的 η 是步长因子μ 是防止分母为零的正则因子。如果 Δu(k-1) 等于零或者非常小分母会拖垮整个估计所以 μ 必须取一个正的小值我通常会在 0.5 到 5 之间选。还有一步必须做重置机制。当估计出来的 PPD 绝对值太小、符号发生变化或者 Δu 接近零时把 φ_c(k) 重置为初始值。这是保证系统不发散的关键后面第 4 节会详细展开。2.2 MFAPC滚动优化与预测时域MFAPC 在动态线性化基础上引入预测控制的核心要素也就是滚动时域优化。控制目标不是只看下一拍而是让未来 N_p 步的预测输出向参考轨迹靠拢同时抑制控制增量的剧烈变化因此目标函数写成[ J \sum_{i1}^{N_p} \left( y^*(ki) - \hat{y}(ki) \right)^2 \lambda \sum_{i1}^{N_u} \Delta u(ki-1)^2 ]其中 y*(ki) 是未来第 i 步的参考值ŷ(ki) 是基于 PPD 数据模型预测的输出λ 是控制增量惩罚因子N_u 是控制时域。预测输出怎么来在预测时域内如果认为 PPD 在当前时刻附近变化不大就可以用当前估计值 φ_c(k) 递推未来输出变化[ \hat{y}(ki) y(k) \hat{\phi}c(k) \sum{j1}^{i} \Delta u(kj-1) ]对目标函数求极值在控制时域取 1 的情况下可以推出 MFAPC 的控制律[ \Delta u(k) \frac{\hat{\phi}_c(k)}{\lambda N_p \cdot \hat{\phi}c(k)^2} \sum{i1}^{N_p} \left( y^*(ki) - y(k) \right) ]这个公式是 MFAPC 仿真里最核心的一行代码。它的结构很有启发预测时域 N_p 越大参考轨迹的未来信息被利用得越多但分母也会随之变大控制动作会更加保守λ 越大控制增量受到的惩罚越强系统响应越平缓φ_c(k) 可以被理解为当前工作点的增益增益越大控制增量应该给得越足。MFAPC 与传统 MPC 相比最大的区别在于预测模型。MPC 用的是状态空间模型或者传递函数模型MFAPC 用的是 PPD 描述的增量数据模型。你用不着花大力气做机理建模代价是 PPD 的估计质量直接决定整个预测精度所以 PPD 估计环节是整个 MFAPC 编程实现里最需要下功夫的地方。2.3 MFAILC迭代轴上的学习律MFAILC 解决的是另一类问题。系统在有限时间区间 [0, T] 上重复运行每个批次的参考轨迹相同、初始条件接近目标是通过学习让输出误差在批次之间逐步减小。第 i 次迭代的输出误差定义为[ e_i(k) y^*(k) - y_i(k) ]MFAILC 的学习律基于伪偏导数数据模型把上一批次的误差换算成当前批次的控制输入修正量[ u_{i1}(k) u_i(k) \frac{\rho \hat{\phi}_c(k)}{\lambda \hat{\phi}_c(k)^2} e_i(k1) ]这个形式与 MFAPC 的控制律有几分相似但含义完全不同。MFAPC 是时间轴上的滚动修正每一拍根据当前估计的未来预测决定控制增量。MFAILC 是迭代轴上的批次修正当前批次完全沿用上一批次的控制序列只在上一批次误差较大的时刻叠加一个修正量。ρ 是迭代学习增益决定每次修正的幅度同一公式里做归一化是为了防止 PPD 过小时修正量失控。打个生活化的比方MFAPC 像开车看导航每开一段路重新规划接下来几秒的路径MFAILC 像练投篮投完一球对比落点与篮筐的偏差下一次出手时调整力度和角度。你不需要知道投篮的完整空气动力学模型只需要根据偏差反复修正。实现 MFAILC 时要注意一层嵌套关系外层是迭代轴内层是时间轴。每次外层迭代跑完整段时间轴才能拿到一整条的误差序列供下一次外层更新 PPD 和控制序列。这个结构在代码里就是两个 for 循环但很多人第一次写会把迭代轴写成步进式的那就失去了“批次学习”的意义。3. 仿真程序设计参数、代码骨架与结果指标算法原理清楚了接下来就是仿真的工程实现。这一节我会给出一个可以直接扩展的仿真程序骨架包含被控对象、控制器主循环、评价指标和出图逻辑并解释每一步的作用。3.1 被控对象选取与参考轨迹设定仿真程序要可信被控对象不能太简单也不能太难。太简单的线性对象无法体现无模型控制算法的价值太复杂的对象又难以定位问题出在算法还是出在仿真环境。我选了一个带状态耦合和非线性项的二阶差分系统[ y(k1) 0.6 \sin(y(k)) 0.4 y(k) 1.2 u(k) - 0.1 u(k)^2 ]这个对象有三个值得测试的特性第一sin 项让系统输出在较大摆动时表现出明显非线性第二输出项系数 0.4 保证系统在没有控制输入时是稳定的方便观察控制器从零开始建立跟踪第三u(k)^2 项引入输入相关的非线性这会让基于线性增量假设的 PPD 模型面临真实挑战仿真结果更接近现场感受。参考轨迹我建议先用幅值交替变化的方波再试正弦波。方波考验的是控制器的快速跟踪能力和超调抑制能力正弦波考验的是连续跟踪能力和平滑性。两条轨迹都要跑分别记录结果。仿真时长我取 200 步采样周期设为 1这个长度足够观察方波跳变后的过渡过程和正弦波的周期跟踪效果。3.2 MFAPC 主循环代码骨架MFAPC 主循环可以分成几个函数被控对象函数、PPD 估计函数、控制器计算函数。我用 Python 风格的伪代码来表达核心逻辑可以直接迁移到 MATLAB 或 C 代码。import numpy as np import matplotlib.pyplot as plt # 被控对象非线性离散系统 def plant(y_k, u_k, u_prev): return 0.6 * np.sin(y_k) 0.4 * y_k 1.2 * u_k - 0.1 * u_prev # 伪偏导数估计 def ppd_estimate(phi_prev, dy, du_prev, mu, eta): denom mu du_prev ** 2 phi_hat phi_prev eta * du_prev * (dy - phi_prev * du_prev) / denom # 重置机制 if abs(phi_hat) 1e-4 or np.sign(phi_hat) ! np.sign(phi_prev): phi_hat 0.5 return phi_hat # 参数设置 T 200 # 仿真步数 N 5 # 预测时域 lam 0.1 # 控制增量惩罚 mu 1.0 # PPD 估计正则因子 eta 1.0 # PPD 估计步长 phi 0.5 # PPD 初始值 y np.zeros(T) u np.zeros(T) y_ref np.zeros(T) y_ref[50:150] 1.0 # 方波参考轨迹 y_ref[150:] -0.5 u_prev 0.0 phi_hist np.zeros(T) for k in range(T - 1): # 1. 计算控制增量 err_sum np.sum(y_ref[k1 : k1N] - y[k]) delta_u phi * err_sum / (lam N * phi ** 2) u[k] u_prev delta_u # 2. 施加控制得到新输出 y_next plant(y[k], u[k], u_prev) y[k1] y_next # 3. 在线更新 PPD 估计 dy y_next - y[k] du_prev u[k] - u_prev phi ppd_estimate(phi, dy, du_prev, mu, eta) phi_hist[k] phi u_prev u[k]这段代码的逻辑顺序是固定的先根据当前 PPD 估计值计算控制增量并施加再采集新的输出再用新的输入输出数据更新 PPD。这个“先控后估”的顺序是由因果性决定的——当前时刻的控制器不可能用到当前时刻尚未测得的输出值。3.3 MFAILC 迭代循环代码骨架MFAILC 的代码结构是双层循环外层迭代轴、内层时间轴。需要注意的是外层迭代的 PPD 估计是基于整条上一批次的数据离线完成的还是基于上一批次逐点完成的这个要明确。我选择在每一轮外层迭代结束后用整条数据估计一条 PPD 曲线供下一批次使用。# 参数设置 I 30 # 迭代批次数量 T 100 # 每批次时间步数 rho 0.9 # 迭代学习增益 lam_ilc 0.1 # 归一化因子 u np.zeros((I, T)) y np.zeros((I, T)) e np.zeros((I, T)) # 初始批次先给一个常值控制让系统进入合理工作范围 u[0, :] 0.2 for k in range(T - 1): y[0, k1] plant(y[0, k], u[0, k], u[0, k-1] if k 0 else 0.0) e[0, k1] y_ref[k1] - y[0, k1] # 迭代学习主循环 for i in range(I - 1): # 基于上一批次数据估计 PPD 曲线 phi_seq np.ones(T) for k in range(1, T - 1): du u[i, k] - u[i, k-1] dy y[i, k] - y[i, k-1] phi_seq[k] ppd_estimate(phi_seq[k-1], dy, du, mu, eta) # 计算当前批次控制输入 for k in range(T - 2): u[i1, k] u[i, k] rho * phi_seq[k] * e[i, k1] / (lam_ilc phi_seq[k] ** 2) # 运行当前批次 for k in range(T - 1): y[i1, k1] plant(y[i1, k], u[i1, k], u[i1, k-1] if k 0 else 0.0) e[i1, k1] y_ref[k1] - y[i1, k1]这个骨架里有一个非常重要但容易被忽视的细节初始批次的控制输入不能随便给零。如果初始批次输出离参考轨迹太远第一次迭代修正量就可能过大导致批次间振荡。我通常会让初始批次先用一个保守的常值控制或者用最基础的 MFAC 先跑一条勉强跟踪的曲线再进入 MFAILC 学习流程。3.4 结果评价指标与可视化仿真跑完不能只看曲线好不好看还要量化。我在这个项目里固定用五个指标指标计算方式意义均方根误差 RMSEsqrt(mean(e(k)^2))整体跟踪精度最大绝对误差 MAEmax(abs(e(k)))最坏情况偏差控制输入方差的平均值mean(var(u))控制动作的平稳性收敛步数首次进入误差带内的时刻响应速度迭代收敛批次首次误差指标低于阈值的批次MFAILC 的学习速度出图方面我建议至少画三张图。第一张是输出 y 和参考轨迹 y_ref 对比一眼看出跟踪效果第二张是控制输入 u 的时间曲线检查控制量有没有剧烈抖振第三张是 PPD 估计值的曲线这个图很多人忽略但它能帮助你定位参数问题——PPD 曲线如果高频振荡或者跳到异常大的数值说明估计环节出了问题这时候就算输出曲线暂时好看控制器也是不健康的。4. 实测中的常见问题与调参心得仿真程序跑通只是开始真正花时间的在于参数调整和异常排查。这一节我把自己实测中遇到的问题和解决方法整理出来给读者减少一些试错成本。4.1 PPD 估计发散症状与对策PPD 发散是最常见的仿真失败模式。症状很明显控制输入曲线出现高频大幅抖振输出像是被什么东西反复拉扯严重时数值直接溢出变成 NaN。我遇到的第一类原因是分母趋零。当 Δu(k-1) 接近零时PPD 估计公式的分母 μ Δu² 接近 μ但分子里的 Δu(k-1) 同样接近零两者相除可能出现数值不稳定。需要确保 μ 不小于 0.1并给 Δu 设置一个下限保护比如 Δu 的绝对值小于 1e-5 时跳过 PPD 更新。第二类原因是重置机制缺失。没有重置或重置条件设置得过松PPD 一旦跑到负值且系统中又存在负的局部增益变化估计值就会越跑越偏。我建议重置条件设为“绝对值小于阈值”或“符号与上一次估计相反”两种情况都触发重置值取初始值 0.5 是一个稳妥的选择。实测中这套重置机制在大多数对象上都能避免发散。第三类原因是步长因子 η 过大。η 越大PPD 跟踪越快但估计值噪声也越大。如果系统本身存在较强非线性或者测量噪声η 超过 2 时 PPD 曲线会明显粗糙。我一般从 1.0 起步如果发现 PPD 曲线高频抖动就降一半。4.2 预测时域与控制时域的调参拉扯MFAPC 里最影响性能的参数就是预测时域 N 和惩罚因子 λ它们之间存在明显的跷跷板效应。N 太小控制器退化成近似单步 MFAC遇到系统滞后时容易出现持续振荡因为控制动作只盯着当前这一拍没有考虑未来一段时间的变化趋势。我在方波参考轨迹上测试过N 取 2 时输出在方波跳变处会超调 20% 以上N 提高到 8 之后超调明显下降。但 N 太大又有另一个问题控制器过于依赖对未来参考轨迹的预测如果参考轨迹本身是突变的或者未来参考值不可预知控制动作会变得保守、迟钝。λ 的作用是抑制控制增量。λ 太小控制量会猛冲猛打尤其是被控对象增益偏大的时候λ 太大输出跟踪变慢方波响应像被什么东西拖着走。我的调参顺序是先把 N 固定在 5把 λ 从 1 开始往下调观察输出出现振荡时把 λ 回调一倍然后再逐步增加 N 观察改进空间。这样调比同时动两个参数要容易定位问题。MFAILC 的学习增益 ρ 同样需要谨慎。ρ 太大批次间修正过猛误差曲线不是单调下降而是来回振荡甚至发散ρ 太小学习进度过慢几十个批次都压不到目标误差。经验上 ρ 在 0.5 到 1.5 之间比较安全先取 0.8观察前五个批次的误差变化趋势再微调。4.3 测量噪声与初始批次的影响现场数据不可能干净仿真时加入测量噪声才能暴露算法的真实鲁棒性。我在对象输出上加均值为 0、标准差为 0.01 的高斯噪声时两种算法都能保持基本跟踪但 PPD 估计曲线会变得粗糙。对策是在 PPD 估计更新时引入一个平滑因子 α将更新公式改成加权形式[ \hat{\phi}_c(k) \hat{\phi}_c(k-1) \alpha \cdot \frac{\eta \Delta u(k-1)(\Delta y(k) - \hat{\phi}_c(k-1)\Delta u(k-1))}{\mu \Delta u(k-1)^2} ]α 取 0.5 到 0.9 之间相当于给估计值加了一阶低通滤波。代价是 PPD 对系统增益变化的响应变慢但换来的是估计曲线稳定得多。在强噪声环境下我宁可选稍大的平滑系数也不愿意让控制器跟着噪声乱动。MFAILC 的初始批次影响很大。如果初始批次完全用零控制输入输出可能落在距离参考轨迹很远的平衡点第一次迭代的误差非常大修正量会直接把控制输入推到不合理的数值。我的经验是先跑一小段时间基础 MFAC让系统先进入参考轨迹附近的区域再把这条控制序列作为 MFAILC 的初始批次。没有这一步MFAILC 的收敛曲线会有一个明显陡峭的第一跳很不好看。4.4 参数速查表与推荐调整方向把我在不同对象上反复测试后得出的参数经验整理成一张速查表方便读者起步时参考参数符号推荐初始值调整方向伪偏导初始值φ(1)0.5过小易发散过大易振荡PPD 步长因子η1.0噪声大时降跟踪慢时升PPD 正则因子μ1.0分母奇异时适当加大PPD 平滑因子α0.8噪声大时增加平滑预测时域N5滞后大时增加响应慢时减少控制惩罚因子λ0.1振荡时增大跟踪慢时减小迭代学习增益ρ0.8收敛慢时增大振荡时减小这张表不保证在所有对象上都最优但作为起点是靠谱的。每个参数在调整时都应该一次只动一个观察输出曲线和 PPD 曲线两个方面的变化再决定下一步动作。同时记录每次实验的 RMSE 和控制方差形成一个小规模参数扫描表比在代码里盲调要高效得多。5. 写在最后的个人实操体会这套仿真程序我前前后后跑了三四轮每次都有新的发现。第一次跑通时我以为算法很顺利结果只是参数恰好选得合适当我故意把参考轨迹改成带突变的复杂波形时PPD 估计就暴露出了问题。所以如果你准备复现或者扩展这套验证流程我最大的建议是不要只盯着输出跟踪曲线多画 PPD 估计曲线多用几个不同特性的参考轨迹多一点耐心做参数扫描。真正让你对 MFAPC 和 MFAILC 建立信心的不是某一次偶然的好结果而是你知道在哪些参数范围内算法能收敛、在哪些条件下它会有怎样的异常表现。这个知识只有亲手跑过仿真代码才会长在自己身上。