★ FEATURED ARTICLE 去年我被一套温控对象折腾了整整一周。对象本身并不复杂一阶惯性加点非线性麻烦在于它的特性会随温度段漂移。常温段调好的Kp、Ki、Kd到高温段就震荡改完高温低温又变成响应拖沓。来回试凑了很多组固定参数最终只能在中间地带凑合。这大概是所有用固定参数PID的人都会撞到的墙三参数本质上是对某个工作点附近对象特性的线性近似对象一变近似就失效。要根治得让参数自己跟着对象走。我当时采用的做法是把径向基神经网络RBF接进控制回路做在线辨识用RBF输出对象的Jacobian信息再按梯度下降方向每个控制周期微调PID参数这就是基于RBF的PID参数自适应整定。这篇文章会把原理推导、完整可运行的代码、以及我实测中踩到的坑都放出来给同样跟参数整定较劲的工程师和学生一份能直接参考的方案。1. 固定PID参数的天花板非线性与工况漂移带来的麻烦1.1 常规整定方法到底在做什么Ziegler-Nichols、临界比例度法、响应曲线法这些常规方法说白了就是在一组设定工况下做一个线性化实验。测试对象在某个稳态工作点附近的阶跃响应然后用经验公式映射出一组PID参数。这套方法对线性时不变系统非常成熟也足够高效但有两个前提容易被忽略对象工作点不能大幅移动对象特性不能随时间和环境变化。一旦前提不成立参数就只能在曾经测过的那个点附近可用。控制界有个很形象的说法固定PID参数相当于给人配一副定焦眼镜度数只对一个距离合适换距离就得换眼镜。1.2 我踩过的那类工况漂移到底指什么拿温度对象举例常温下散热系数小、热容大等效增益高高温段散热加剧等效增益下降。也就是说对象传递函数本身就在缓慢变化。类似的现象在电机负载、液压系统、材料挤出工艺里都很常见。更极端一点如果系统里还有死区、饱和、摩擦等非线性环节固定PID就更容易出现在某些区段振荡、某些区段迟钝的症状。我那时候调温控程序发现同一个Kp在目标温度30°C时稳定裕度很好到80°C就出现可观察的振荡这就是对象增益变化最直观的体现。1.3 两条自适应路线的取舍要让PID参数跟着对象动工程上大致两条路线。路线A用神经网络直接把误差、误差变化率映射成三参数本质上是做一个参数规划器由离线训练或在线学习决定映射关系。路线B先用神经网络在线辨识对象的输入输出映射求出对象对控制输入的灵敏度∂y/∂u工程上常把这个灵敏度叫作Jacobian信息再把这个灵敏度接入PID参数的梯度下降更新。路线B物理意义更清晰参数更新方向可以推导、可以调试遇到异常也容易定位所以我最终采用路线B。整篇文章的推导和代码都围绕这条路线展开。2. RBF网络在控制回路里的真实角色在线辨识器与Jacobian计算器2.1 径向基结构为什么适合在线辨识RBF网络的结构不复杂输入层、一个含径向基函数的隐含层、一个线性输出层。隐含层每个节点有一个中心向量c和一个宽度b它对输入x的响应是h_j exp(-‖x-c_j‖² / (2b_j²))。离中心越近节点输出越大离得远输出迅速衰减。这种局部响应的特性非常像一堆小核函数做插值所以它比BP那种全局逼近网络更适合在线快速调整某个区域的输入激活对应的局部节点更新互不干扰不容易出现改一个权值影响全局的问题。对于控制回路里的在线辨识这个局部性很宝贵。2.2 输入向量怎么选影响Jacobian的第一道关口我在代码里取输入x [u(k-1), y(k-1), y(k-2)]也就是上一拍控制量、上一拍输出、上上拍输出。选这三项的考虑是常见一阶非线性对象可以写成y(k) f(y(k-1), u(k-1))再多加一项y(k-2)是为了容纳系统中可能存在的二阶动态又不会让网络输入维度膨胀。必须强调一个经验x的第一维一定要放控制量u因为后面算Jacobian时直接对第一个输入节点求偏导数放其它位置代码会绕。如果对象阶数更高把y(k-3)、y(k-4)甚至r(k-1)一并加上RBF输入维度一般不建议超过5到6否则局部节点数目要跟着涨训练数据不够时很难收敛。2.3 权值、中心和宽度怎么在线更新RBF在线辨识的目标是最小化J 0.5(y(k) - ym(k))²。ym是RBF当前输出y(k)是对象真实输出两者在同一个控制周期里都能拿到不需要等待未来数据。梯度下降更新包含三层参数输出权值w_jΔw_j η_w * err_rbf * h_j中心c_jiΔc_ji η_c * err_rbf * w_j * h_j * (x_i - c_ji) / b_j²宽度b_jΔb_j η_b * err_rbf * w_j * h_j * ‖x - c_j‖² / b_j³这里err_rbf y - ym是网络辨识误差不是控制误差两者必须隔离开。公式看起来多但每个参数的物理意义都很清楚权值负责输出幅度中心负责局部区域定位宽度负责影响范围。三者都是沿着误差减小的方向走一步。到这里先不手推全部公式代码里都有对应实现下面推导控制误差的时候就明白这个Jacobian为什么是整套方案的核心。2.4 控制回路里最关键的∂y/∂u怎么算网络输出对输入向量第i维的偏导数可以沿隐藏层直接算出来对第一维控制量u来说就是∂ym / ∂u Σ_j w_j * h_j * (c_j0 - x_0) / b_j²这个值就是前面说的Jacobian信息工程上常常写成dyu。它回答的问题是如果我让控制量u增加一点点对象输出y会朝哪个方向、以多大强度变化。有了这个灵敏度后面PID参数向哪个方向调就有了依据这也是整套自适应整定方案里信息量最大的一个数字。dyu的值由RBF当前逼近的对象局部特性决定对象本身是非线性的、时变的dyu就会跟着变这正好弥补固定参数PID只认一个工作点的缺陷。3. 从误差异方差到PID增量自适应整定公式的完整推导3.1 控制目标和梯度下降的方向把控制目标写成最小化瞬时误差平方E(k) 0.5 * e(k)²其中e(k) r(k) - y(k)。要对PID三个参数做梯度下降需要计算∂E/∂Kp、∂E/∂Ki、∂E/∂Kd。直接对E求偏导并不困难难点在于对象模型未知∂y/∂Kp这样的传递关系我们手里没有。解决办法是做一个链式展开∂E / ∂Kp -e * ∂y / ∂u * ∂u / ∂Kp其中∂y/∂u就用RBF网络上一节算出来的dyu替换。这样三个参数梯度全部可算这也是自适应这个称谓的由来——每个控制周期都用当前对象的最新灵敏度信息去修正控制器参数而不是用一套固定值走到底。3.2 增量式PID对参数的偏导我选用增量式PID而不是位置式PID一个很实际的原因是增量式输出的是控制量的变化量Δu天然适合限幅也便于在工程上实现抗积分饱和。增量式控制律写成Δu(k) Kp * (e(k) - e(k-1)) Ki * e(k) Kd * (e(k) - 2e(k-1) e(k-2))那么三个偏导分别是∂u / ∂Kp e(k) - e(k-1)∂u / ∂Ki e(k)∂u / ∂Kd e(k) - 2e(k-1) e(k-2)这三项在代码里就是de、e、d2e直接复用计算量几乎可以忽略。3.3 三个参数的更新公式合并链式展开结果并取负梯度方向梯度下降可得ΔKp ηp * e(k) * dyu * (e(k) - e(k-1))ΔKi ηi * e(k) * dyu * e(k)ΔKd ηd * e(k) * dyu * (e(k) - 2e(k-1) e(k-2))Kp(k1) Kp(k) ΔKpKi、Kd同理。这里要特别留意符号方向因为E对Kp的偏导带负号负负得正后正是上面形式如果对象增益为负dyu本身会给出负号算法依然自动沿正确方向调整。这个特性在处理正反向对象时非常方便不需要手动切换符号。3.4 为什么需要预热期和参数限幅公式刚出来的时候最容易让三个参数从第一步就开始猛更新。实际跑下来发现必须加两个保护措施。第一是预热期前若干周期只更新RBF辨识器、不更新PID参数让dyu先落在合理范围内我常用30拍预热。因为初始权重随机、中心随机RBF一开始算出来的dyu可能是噪声用它去调PID参数会把系统带偏。第二是参数限幅Kp、Ki、Kd分别clip在合理区间内我代码里给的是Kp∈[0.01, 5.0]、Ki∈[0.0, 2.0]、Kd∈[0.0, 1.0]。限幅不是为了限制性能而是防止噪声环境下参数漂到离奇数值项目上线时这几乎是保命的设置。4. 能直接跑的Python实现RBF网络类、主循环与可视化4.1 文件结构和替换思路代码我拆成清晰的块按照标题说的详细代码说明文件来组织rbf_network.py放RBF网络类plant_model.py放被控对象main_sim.py放主循环和画图。实际用的时候只需要替换plant_model.py里的对象模型其它部分一般不用动。对象的输入输出接口就是step(u)传入控制量返回新输出。隐藏层节点数默认6对一阶对象足够如果你的对象输入输出关系更复杂先加到8到10试试再往上涨意义不大反而容易在线学习不稳。4.2 RBF网络类完整代码import numpy as np class RBFNetwork: def __init__(self, n_in3, n_hidden6, lr_w0.03, lr_c0.01, lr_b0.01): self.n_in n_in # 输入维度 self.n_hidden n_hidden # 隐含层节点数 self.lr_w lr_w # 输出权值学习率 self.lr_c lr_c # 中心学习率 self.lr_b lr_b # 宽度学习率 # 中心在输入可能范围内均匀分布宽度在 0.6~1.2 之间 self.c np.random.uniform(-1.5, 1.5, (n_hidden, n_in)) self.b np.random.uniform(0.6, 1.2, n_hidden) self.w np.random.uniform(-0.5, 0.5, n_hidden) def _hidden(self, x): h np.zeros(self.n_hidden) for j in range(self.n_hidden): dist2 np.sum((x - self.c[j]) ** 2) h[j] np.exp(-dist2 / (2.0 * self.b[j] * self.b[j])) return h def predict(self, x): return float(np.dot(self._hidden(x), self.w)) def jacobian_u(self, x): # 对第一个输入维控制量u求偏导 h self._hidden(x) dyu 0.0 for j in range(self.n_hidden): dyu self.w[j] * h[j] * (self.c[j, 0] - x[0]) / (self.b[j] ** 2) return dyu def update(self, x, err_rbf): h self._hidden(x) # 输出权值梯度下降 for j in range(self.n_hidden): self.w[j] self.lr_w * err_rbf * h[j] # 中心和宽度的梯度下降 for j in range(self.n_hidden): dist2 np.sum((x - self.c[j]) ** 2) for i in range(self.n_in): # 中心c_ji的梯度 dc self.lr_c * err_rbf * self.w[j] * h[j] * (x[i] - self.c[j, i]) / (self.b[j] ** 2) self.c[j, i] dc # 宽度b_j的梯度 db self.lr_b * err_rbf * self.w[j] * h[j] * dist2 / (self.b[j] ** 3) self.b[j] np.clip(self.b[j] db, 0.2, 3.0)这段代码里有一个容易被忽略的细节中心和宽度的更新都用到了w_j本身如果w_j已经收敛到0附近后续中心和宽度基本不会再变动这其实是RBF局部性的体现——该节点已经退休由其它节点负责当前区域。实际跑的时候不需要特别干预只要隐藏层节点数不是太少系统会自动安排好。4.3 主循环与PID自适应更新完整代码import numpy as np import matplotlib.pyplot as plt class NonlinearPlant: 一阶非线性被控对象y(k1) a*y(k)/(1y(k)^2) u(k) 参数a用于模拟对象增益的离线/在线漂移。 def __init__(self, a1.0): self.a a self.y 0.0 def step(self, u): y_next self.a * self.y / (1.0 self.y ** 2) u self.y y_next return y_next def run_simulation(adaptiveTrue, seed0): np.random.seed(seed) rbf RBFNetwork(n_in3, n_hidden6) plant NonlinearPlant(a1.0) # 初始PID参数故意给一个不算最优的起点 Kp, Ki, Kd 0.8, 0.15, 0.03 # PID参数整定学习率 eta_p, eta_i, eta_d 0.05, 0.02, 0.01 u_prev 0.0 y_prev 0.0 y_prev2 0.0 e_prev 0.0 e_prev2 0.0 sim_steps 1000 warmup 30 # 预热期只更新RBF不更新PID参数 r_hist, y_hist, u_hist [], [], [] kp_hist, ki_hist, kd_hist [], [], [] for k in range(sim_steps): # 参考信号前750拍保持1之后切到0.5 r 1.0 if k 750 else 0.5 # 第350拍模拟对象增益突变 if k 350: plant.a 2.0 y plant.y e r - y de e - e_prev d2e e - 2.0 * e_prev e_prev2 # RBF在线辨识用[u(k-1), y(k-1), y(k-2)]预测当前y(k) x_rbf np.array([u_prev, y_prev, y_prev2]) ym rbf.predict(x_rbf) err_rbf y - ym rbf.update(x_rbf, err_rbf) # 计算Jacobian即输出对控制量的灵敏度 dyu rbf.jacobian_u(x_rbf) # 自适应预热期后按梯度更新PID参数 if adaptive and k warmup: Kp np.clip(Kp eta_p * e * dyu * de, 0.01, 5.0) Ki np.clip(Ki eta_i * e * dyu * e, 0.0, 2.0) Kd np.clip(Kd eta_d * e * dyu * d2e, 0.0, 1.0) # 增量式PID计算控制量 du Kp * de Ki * e Kd * d2e u np.clip(u_prev du, -2.0, 2.0) y_next plant.step(u) r_hist.append(r) y_hist.append(y_next) u_hist.append(u) kp_hist.append(Kp) ki_hist.append(Ki) kd_hist.append(Kd) # 滚动历史变量 u_prev u y_prev2 y_prev y_prev y e_prev2 e_prev e_prev e return r_hist, y_hist, u_hist, kp_hist, ki_hist, kd_hist主循环里的滚动顺序要细心看y_prev保存的是进入第k拍之前的plant.y也就是y(k-1)y_prev2自然就是y(k-2)。x_rbf的三个输入全部来自上一拍之前的历史值而监督信号err_rbf用的是当前拍y(k)这样的时序是因果可实现的不会用到未来数据。实物部署时只需要把plant.step换成真实对象的读取和执行操作其它逻辑完全一致。4.4 绘图与对照实验def plot_results(r_hist, y_hist, u_hist, kp_hist, ki_hist, kd_hist): t np.arange(len(r_hist)) fig, axes plt.subplots(3, 1, figsize(10, 8), sharexTrue) axes[0].plot(t, r_hist, k--, labelreference) axes[0].plot(t, y_hist, labeloutput) axes[0].set_ylabel(y / r) axes[0].legend() axes[0].grid(True) axes[1].plot(t, u_hist, labelcontrol) axes[1].set_ylabel(u) axes[1].legend() axes[1].grid(True) axes[2].plot(t, kp_hist, labelKp) axes[2].plot(t, ki_hist, labelKi) axes[2].plot(t, kd_hist, labelKd) axes[2].set_ylabel(PID gains) axes[2].set_xlabel(step k) axes[2].legend() axes[2].grid(True) plt.tight_layout() plt.show() # 跑自适应 r_hist, y_hist, u_hist, kp_hist, ki_hist, kd_hist run_simulation(adaptiveTrue) plot_results(r_hist, y_hist, u_hist, kp_hist, ki_hist, kd_hist)想看对照组也很简单再调一次run_simulation(adaptiveFalse)此时Kp、Ki、Kd始终固定在初始值。把两组y_hist叠在一起画对象突变后的差距会非常明显。这里adaptiveFalse的那组本质上就是一组没有人工细调过但也能勉强跑的固定PID和自适应组对比体现的是自动走梯度带来的收益。4.5 说明文件建议骨架既然标题明确提到说明文件我建议代码仓库里至少放这几块内容。README部分我一般这么组织每个文件的职责、每个类的初始化参数含义、如何把NonlinearPlant换成自己的对象、学习率调节顺序。还有一个参数速查表把eta_p/eta_i/eta_d、lr_w/lr_c/lr_b、预热期、参数限幅全部列出来标注先调RBF学习率再调整定学习率的顺序。参数之间是联动的不要同时改好几个否则一旦曲线异常你根本不知道是哪一项造成的。我每接手一个新对象第一步永远是重置参数到这套初始值确认曲线形态正常再开始往实际工况上靠。5. 仿真结果怎么读对象突变时自适应到底在干什么5.1 我用的实验场景仿真设置我说明一下方便你对照复现前750拍参考输入为1第750拍切到0.5第350拍让被控对象参数a从1突变到2模拟增益漂移。对照组的固定PID参数没有刻意调优就是初始值0.8/0.15/0.03自适应组也从同一组初始值出发。两组对象模型完全相同唯一的差别是PID参数是否在线调整。5.2 三个容易观察到的现象以我实际跑这个场景的典型结果来看有三个现象值得盯住。现象一预热期后跟踪很快建立。RBF辨识误差迅速下降dyu稳定在一个合理数值输出曲线在前100拍内贴住参考。如果你看到前30拍输出还在乱飘不要急着改代码那是预热期里网络还在找状态正常现象。现象二第350拍对象突变后固定PID的输出开始出现可观察的振荡或拖尾重新稳定需要较长时间自适应组的输出在突变后先有一个尖峰随后Kp、Ki、Kd根据dyu被修正一般几十拍后重新贴回参考。尖峰无法完全避免因为突变发生在对象内部控制器只能在后一拍感知再优秀的自适应也需要一到两拍的响应时间。好的自适应结果是尖峰小、恢复快而不是完全没有尖峰。现象三第750拍参考值从1切到0.5自适应组能平滑过渡说明调整后的参数在新工作点仍然有效。这一点很关键它说明网络学会的不是某一个点的局部映射而是当前工况附近的动态关系。5.3 如果曲线不符合预期按什么顺序查输出发散先看dyu的符号和量级把它单独打印出来。如果dyu明显振荡或者出现过大的数值问题大概率在RBF网络的学习率而不是PID参数的学习率。如果dyu正常但误差不收敛再调大eta_p、eta_i、eta_d每次只调一个。还有一个常见情况是Kp曲线高频抖动这会直接传导到控制量通常加一个小动量项就能压住。这个排查顺序我踩了好多坑才总结出来比瞎调参数高效得多。6. 落地时真正容易踩的坑和三个改进方向6.1 五个必踩的坑先说最容易踩的坑几乎每个照代码自己跑的人都会碰到。第一辨识误差和控制误差混用。err_rbf是y - ym控制误差是r - y。混用的结果是参数更新方向错乱看起来每条公式都对跑起来就是发散。检查方法是打印两个误差的历史曲线它们应该在不同量级和相位。第二中心和宽度随手初始化。中心范围必须覆盖输入的可能区间否则部分节点永远不激活。我的经验是输入范围在[-2, 2]左右时中心取[-1.5, 1.5]均匀分布宽度取输入范围的0.3倍左右。第三学习率一次调很多。RBF的w学习率从0.03、中心0.01、宽度0.01起步PID整定学习率从0.05、0.02、0.01起步。每次只动一个量观察100拍曲线再决定下一步。第四忽略采样周期稳定性。RBF在线辨识默认等间隔采样。实物里如果主循环里有耗时任务导致控制周期忽长忽短dyu会引进额外扰动。用定时器中断或者实时线程跑控制循环别把它丢在普通业务逻辑里。第五参数限幅抄我的不匹配对象。Kp上限5.0、Ki上限2.0是给这个仿真对象用的。对象量纲、比例系数差几个数量级时限幅必须跟着改否则自适应被限幅卡死等于没自适应。6.2 改进方向一给参数更新加动量项动量项公式不复杂ΔKp_new (1 - α) * ΔKp_old α * ηp * e * dyu * deα一般取0.3左右。它能让参数更新的轨迹更平滑尤其在对象噪声比较大的时候能明显减少Kp的高频抖动。代价是响应速度变慢所以α不要超过0.6。实现上就是多维护一个Kp的上一步增量变量每拍做一次加权平均改动量很小。6.3 改进方向二对dyu做滑动平均dyu本身是RBF的局部输出对输入的偏导在线更新时会有波动直接拿它导出的参数自然也跟着波动。我后来在温度控制上给dyu加了一个滑动平均窗口取5到15拍输出曲线平滑很多。更讲究的做法是加带通滤波把低频趋势留下、高频噪声滤掉。注意窗口太大会引入滞后导致参数更新跟不上对象的突变。6.4 改进方向三多工况识别或变学习率更进一步的方案是把对象工作区间分段每个区间各训练一个小RBF辨识器再根据当前y和r切换到对应模型。这比单个大网络更稳缺点是需要预先知道大概的工作区间。另外也可以用自适应学习率比如误差的滑动平均值上升时自动降低学习率避免参数发散。这些属于工程优化核心结构仍然是RBF在线辨识Jacobian加PID梯度整定理解了这层各种扩展都不难。我现在接新的控制项目遇到难整定的对象第一反应已经很少是翻经验公式而是先问一句对象特性稳不稳定固定参数顶不顶得住。如果答案是顶不住我基本就是套这套RBF自适应整定的框架先把结构跑通再根据对象慢慢调学习率。这套方案说不上完美但它是那种我知道每一步在干什么、每条公式为什么长这样的方案对我来说比黑盒调参舒服太多。如果你也正被PID参数整定反复折磨建议先找个一阶非线性对象跑一遍这套代码理解dyu和三个更新公式之后再迁移到自己的对象上。花几个小时把原理走通后面能省下无数个调参的周末。
阅读完成 · 觉得有帮助？
