粒子群算法（PSO）优化PID参数整定：从仿真到工程落地

简介：面向自动控制与智能优化算法学习者，这份基于MATLAB/Simulink环境的压缩包提供了一个使用粒子群优化算法进行PID控制器参数整定的完整示例。资源共包含七个文件，其中两个脚本分别实现粒子群搜索逻辑与误差追踪评估，三个模型用于搭建被控对象和控制系统仿真环境，可模拟不同控制策略、对比优化前后的系统响应，两个辅助文件保存优化过程数据或中间结果。整个压缩包仅26KB，结构清晰，便于快速部署与二次修改。目前已有三百一十八人浏览学习。通过该示例，读者能够理解粒子群算法与PID控制相结合的核心思想，掌握粒子初始化、适应度函数设计、参数迭代更新及仿真验证的完整流程，并可将同一套方法迁移到其他控制器参数优化问题中，适合作为自动控制课程设计、智能优化算法入门以及相关课题研究的参考资料。

1. 调 PID 遇到瓶颈时，PSO 为什么值得试

tunning-PID-by-PSO.rar这类资源在专业论坛和网盘里流传很广，解压出来通常是一段把粒子群算法（PSO）和 PID 调参绑在一起的仿真代码。tunning是tuning的手误，但不妨碍它成为检索关键词。新手上路时，第一反应往往是抄里面的适应度函数和 PSO 循环，改改模型就“跑通”了，可一旦拿到真实电机或者温控箱上，参数又完全不听话。

PID 调参之所以让人头疼，在于误差曲面不是光滑的单峰函数。同一个 Kp、Ki、Kd 组合，在仿真模型里表现很好，换一台设备就可能因为噪声、死区、延迟而振荡。PSO 的核心价值不在“搜出一个精确解”，而在于它不依赖梯度，能在有多个局部极小值的参数空间里做全局搜索，并且天然适合并行评估。这件事对手动试凑和 Ziegler-Nichols 整定来说是补充，不是替代。这篇文章会从 PSO 的算法骨架讲起，然后把一阶惯性加滞后对象的调参实验完整复现出来，最后落到 STM32 和 PLC 上真正能用的离散形式与防饱和处理。面向的读者是那些已经写过 PID 代码、却被参数调整折磨过，想用群体智能把“经验活”变成“半自动流程”的工程师。

2. 粒子位置装的就是 PID 参数：PSO 调参的骨架与目标函数

2.1 把 Kp、Ki、Kd 当成三维空间里的一个点

PSO 在做 PID 调参时，粒子群里的每个粒子都表示一组候选 PID 参数。比如粒子位置向量x = [x0, x1, x2]，直接映射为[Kp, Ki, Kd]。粒子的速度向量v表示这一组参数在每次迭代中的变化方向和幅度。算法初始化时，粒子散落在预设的参数边界内，比如Kp ∈ [0, 5]，Ki ∈ [0, 2]，Kd ∈ [0, 1]，这个范围直接决定了搜索空间的大小。

更新公式是 PSO 的灵魂：

v_i_new = w * v_i_old + c1 * r1 * (pbest_i - x_i) + c2 * r2 * (gbest - x_i) x_i_new = x_i_old + v_i_new

其中pbest_i是该粒子自己历史最优位置，gbest是整个群体共享的全局最优位置。w是惯性权重，控制上一代速度的保留程度；c1和c2是个体学习因子和社会学习因子；r1、r2是[0,1]之间的均匀随机数。每次迭代后，把x_i_new中超出边界的分量直接拉回边界，避免 Kp 变成负数这类物理上无意义的值。

这里的维度数量只有 3，但 PSO 完全可以直接扩展到更高维。比如你在做无人机串级 PID，外环位置环加上内环姿态环，需要同时调整的参数可能是 6 个甚至 9 个，那粒子的维度就对应增加。常见做法是每个粒子的前三维给内环，后三维给外环，适应度函数统一考察整个级联闭环的响应。这样写代码没有额外负担，唯一的代价是搜索空间变大，粒子数量和迭代次数要适当增加。

2.2 适应度函数决定“好”的标准

PSO 只知道找适应度的最小值或最大值，它根本不懂什么叫“超调 20% 太多了”。所有工程经验都必须折算成一个数值。常用的误差积分准则有这些：

我一般会优先选 ITAE，因为它在阶跃响应实验里能让调节时间更短，稳态误差收敛也更利落。做温控 PID 时，温度曲线尾部迟迟达不到设定点，ITAE 会因为这个尾部代价持续累积，迫使算法把 Ki 拉到合适大小。ISE 更偏向“快拉回来”，但容易带来大幅超调。单纯用 IAE 又会忽略长时间存在的低频扰动。

注意，适应度不是只有一个积分项。真实工程里一定要加入惩罚项，比如超调量超过 5% 时，在积分值上追加一个较大常数；或者执行器输出变化率过大时也追加代价。后面第三章会给出具体代码。这里先记住一个原则：适应度函数里没有的东西，PSO 优化出来的参数也不会替你保证。

2.3 最小可运行的 PSO 骨架代码

用 Python 写一个基础 PSO 类非常简单。下面这段代码不依赖任何第三方优化库，只用了numpy：

import numpy as np class PSO: def __init__(self, fitness_func, dim=3, n_particles=20, max_iter=30, w=0.6, c1=1.5, c2=1.5, lb=None, ub=None): self.fitness = fitness_func self.dim = dim self.n = n_particles self.max_iter = max_iter self.w = w self.c1 = c1 self.c2 = c2 self.lb = lb self.ub = ub self._init_swarm() def _init_swarm(self): lb = np.array(self.lb) ub = np.array(self.ub) self.x = lb + (ub - lb) * np.random.rand(self.n, self.dim) self.v = np.random.uniform(-1, 1, (self.n, self.dim)) self.pbest_x = self.x.copy() self.pbest_val = np.array([float("inf")] * self.n) self.gbest_x = np.zeros(self.dim) self.gbest_val = float("inf") def run(self): for it in range(self.max_iter): for i in range(self.n): val = self.fitness(self.x[i]) if val < self.pbest_val[i]: self.pbest_val[i] = val self.pbest_x[i] = self.x[i].copy() if val < self.gbest_val: self.gbest_val = val self.gbest_x = self.x[i].copy() r1 = np.random.rand(self.n, self.dim) r2 = np.random.rand(self.n, self.dim) self.v = (self.w * self.v + self.c1 * r1 * (self.pbest_x - self.x) + self.c2 * r2 * (self.gbest_x - self.x)) self.x = self.x + self.v self.x = np.clip(self.x, self.lb, self.ub) return self.gbest_x, self.gbest_val

这段逻辑很直白：fitness_func接收一个三维向量，返回误差准则值；_init_swarm在每个参数的上下界之间均匀随机撒点，速度在[-1, 1]范围开始；每次迭代先更新历史最优，再统一更新速度和位置，最后clip保证不越界。

实际使用中常见的参数设置如下：

2.4 多峰曲面：为什么 swarm 要限制搜索边界

PID 参数搜索的误差面通常是多峰的。比如一个非线性阀门，死区大小随压力变化，Kp 的某个区间内会产生极限环，适应度突然升高，形成很多“山峰沟壑”。PSO 不像梯度下降那样容易被某个局部极小值困住，但它也有早熟收敛的问题：所有粒子可能在迭代初期就被某个较好但并非最优的位置吸引，群体多样性迅速丢失。

限制搜索边界是最简单的对策，也是工程上最容易操作的。把 Kp 的上下限设得越窄，收敛越快，但风险是真正的最优解不在这个区间里。常见做法是先用 Ziegler-Nichols 或继电反馈测试拿到一组基准参数，然后把上下界设为基准值的 0.1 倍到 2 倍。这样既保留了全局搜索能力，又让粒子聚集在工程合理区域内。第 5 章还会再把这个技巧展开。

3. 复现 tunning-PID-by-PSO：一阶惯性加滞后对象上的完整实验

3.1 仿真对象与离散化方式

工业现场常见的对象可以近似成一阶惯性加纯滞后：

G(s) = K * exp(-L s) / (T s + 1)

这里K是稳态增益，T是惯性时间常数，L是纯滞后。比如加热棒通过可控硅驱动，功率到温度变化之间就符合这个形态。下面代码用一个典型的慢热对象做仿真：K=2，T=40秒，L=6秒，采样周期Ts=0.5秒。

import numpy as np K = 2.0 # 稳态增益 T = 40.0 # 时间常数 L = 6.0 # 纯滞后时间 Ts = 0.5 # 采样周期 sim_len = 4000 # 仿真步数 # 构建离散递推模型 A = np.exp(-Ts / T) B = K * (1 - A) def plant_step(u, y_prev, y_prev2, delay_buf): # 纯滞后通过固定长度缓冲实现 delay_val = delay_buf[-1] delay_buf.append(u) delay_buf.pop(0) y_new = A * y_prev + B * delay_val return y_new

上面的plant_step按一阶惯性差分方程递推：y[k] = A * y[k-1] + B * u[k - L/Ts]，其中u是控制器输出。纯滞后用循环数组存储，顺序是先取最旧的控制量再写入新控制量。注意L/Ts = 12步，所以初始化delay_buf为长度 13 的零数组比较合适。这样离散化比直接调用scipy.signal更直观，也便于以后移植到 C 代码。

3.2 把仿真器放到适应度函数里

适应度函数要把给定的一组[Kp, Ki, Kd]代入闭环系统，跑完整个仿真过程，返回 ITAE 值。这里一定要包含执行器饱和，否则 PSO 会给出一个需要无穷大能量的参数组合。

def pid_simulate(param): Kp, Ki, Kd = param e_prev = 0.0 e_prev2 = 0.0 integral = 0.0 y = 0.0 u = 0.0 setpoint = 1.0 u_min, u_max = 0.0, 1.0 delay_buf = [0.0] * 13 itae = 0.0 control_penalty = 0.0 for k in range(sim_len): e = setpoint - y integral += e * Ts derivative = (e - e_prev) / Ts u_raw = Kp * e + Ki * integral + Kd * derivative # 限制执行器输出 u = max(u_min, min(u_max, u_raw)) # 超限惩罚，防止 PSO 把参数搜到饱和区太深 control_penalty += 100.0 * abs(u_raw - u) * Ts y = plant_step(u, y, 0.0, delay_buf) itae += k * Ts * abs(e) * Ts e_prev2 = e_prev e_prev = e # 超调惩罚 overshoot = max(0.0, (y - setpoint) / setpoint) return itae + overshoot * 1000.0 + control_penalty

代码里的itae += k * Ts * abs(e) * Ts用的是简化矩形积分。实际运行时，k * Ts是当前时刻，积分权重随时间线性增长，这会让 PID 在后面迟迟不收敛时付出极高代价。control_penalty是执行器饱和惩罚，只要控制器输出被限幅，就用一个很大的系数去累积惩罚，所以 PSO 最终会避开那些一上来就把输出打到上限的参数。

3.3 跑一次 PSO 并解释输出

主程序代码很短：

pso = PSO( fitness_func=pid_simulate, dim=3, n_particles=25, max_iter=50, w=0.6, c1=1.5, c2=1.5, lb=[0.0, 0.0, 0.0], ub=[3.0, 1.0, 1.0] ) best_x, best_val = pso.run() print("best PID:", best_x) print("best ITAE:", best_val)

我给一组典型输出作为参考：

best PID: [0.812, 0.094, 0.237] best ITAE: 42.6

这个结果的解读要分两层。第一层，Ki=0.094不大，原因是被控对象时间常数很大，积分作用太强会让响应震荡很久。第二层，Kd=0.237提供了必要的阻尼，但它没有像手动调参时那样被调到很大，因为一阶惯性加滞后对象对微分噪声极度敏感，PSO 自动学会了一个折中方案。如果你把适应度函数里的超调惩罚系数调大，新结果里Kd会明显上升，同时Kp会小幅下降，这就是不同性能指标权衡的直接体现。

3.4 把惩罚项做进去，否则真机上百分之百发散

很多人第一次跑通 PSO 调参后，直接把仿真参数抄进真实设备，结果发现电机嗡嗡作响或者温控曲线大幅震荡。问题往往不在 PSO，而在适应度函数缺少两个工程细节。

第一个是执行器饱和积分饱和。仿真里控制器输出被限幅后，积分项还在继续累积，等到误差反向时，积分项需要很久才能退下来，这就是典型的风饱和。修正方法有两种：一是仿真时在饱和状态下停止积分累加，也就是抗积分饱和逻辑；二是像我上面的代码那样，给饱和程度加一个大的惩罚系数，让 PSO 知难而退。常见做法是两者同时做：仿真内部用抗饱和，外部再用惩罚帮助搜索快速避开饱和区。

第二个是噪声鲁棒性。适应度函数只用阶跃响应评估还不够，建议在误差信号里注入小幅白噪声，或者单独跑一次正弦扰动测试，把扰动最大偏差也折算成惩罚。增量式编码器测速在低速段噪声尤其大，不处理这个，PSO 优化的结果在低速段会出现高频抖振。你不需要把噪声模型做得很精确，只要噪声幅值接近真实传感器的十分之一，就能强迫算法不选择过大的Kd。

4. 从仿真到 STM32、PLC：PSO 结果的离散化与工程落地

4.1 离线调参与在线自整定的分工

真机环境里，PSO 通常不在设备内部跑。STM32 主频再高，跑上百次闭环仿真也不现实；三菱 PLC 更不适合做浮点密集型群体搜索。常见做法是分两步：第一步在 PC 上用 Python、MATLAB 或 Simulink 做离线调参，把被控对象的阶跃响应辨识成一阶模型；第二步把调好的参数写到控制器固件或 PLC 参数寄存器里。这个过程不需要修改设备程序，只需要一个通信接口把更新后的 Kp、Ki、Kd 写入。

对于参数变化较快的设备，比如说飞行器在悬停和高速飞行之间空气动力学参数差别很大，离线一批参数就不够用。这时可以在控制器里做一个“参数表切换”机制：上位机按不同工况分别跑 PSO，生成多组 PID 参数，设备根据飞行状态切换到对应的一组。每一组参数的得出过程仍然是离线算。pid pso的组合在自动驾驶和无人机场景里，多数是沿着这个思路落地的，而不是指望无人机在飞行中自己去跑粒子群。

4.2 离散增量式 PID 的移植要点

把 PSO 输出的连续域 PID 参数转换成离散增量式算法时，最容易错的是微分项系数和采样周期。增量式公式是：

Δu(k) = Kp * (e(k) - e(k-1)) + Ki * Ts * e(k) + Kd * (e(k) - 2*e(k-1) + e(k-2)) / Ts

注意：这里的Ki如果你在 PSO 仿真里已经乘过了Ts，那移植时就不能再乘；否则积分增益会被额外放大。我在仿真代码里用的是integral += e * Ts，也就是说Ki是连续的积分增益，移植到离散增量式就必须在公式里保留Ki * Ts。很多调参工具会把离散 PID 表示成K_i * e(k)这种形式，直接把Ki当离散增益，这时代码里就不能再乘采样周期。这是一个约定问题，统一要在工程文档开头写清楚。

下面是一段适合 STM32 / ARM DSP 平台的整数型增量式 PID C 代码：

typedef struct { float Kp; float Ki; // 连续积分增益 float Kd; // 连续微分增益 float Ts; // 采样周期 float e1; float e2; float out; float out_limit; } PID_INC_T; float pid_inc_calc(PID_INC_T *pid, float err) { float du, u_raw; du = pid->Kp * (err - pid->e1) + pid->Ki * pid->Ts * err + pid->Kd * (err - 2.0f * pid->e1 + pid->e2) / pid->Ts; u_raw = pid->out + du; if (u_raw > pid->out_limit) { u_raw = pid->out_limit; } else if (u_raw < -pid->out_limit) { u_raw = -pid->out_limit; } pid->e2 = pid->e1; pid->e1 = err; pid->out = u_raw; return u_raw; }

这段代码里没有处理积分饱和，因为增量式 PID 本身就包含一定的抗饱和特性：执行器限幅后，限制的是out，下一次计算时du会基于当前被截断的out继续累加，不会出现独立积分项的饱和堆积。但这个特性只对输出限幅有效，如果你用的是位置式 PID，必须显式加抗积分饱和逻辑。

调试时建议把pid->out_limit先设成实际执行器范围的 80%，观察响应再逐步放宽。PSO 在仿真里用的是[0, 1]的标幺输出，到了真实系统如果 PWM 占空比是[0%, 100%]，可以直接按比例放大，但要注意放大后Kp的物理单位也变了，不能只看数字大小。

4.3 串级 PID 与 PLC 场景：时间常数匹配和参数映射

无人机串级 PID 是另一个高频搜索词。外环位置环输出当作内环速度环的设定值，两个 PID 环串在一起，时间尺度必须拉开。经验法则是内环周期是外环周期的十分之一到二十分之一。如果你用 PSO 同时调内外环参数，适应度函数里要把两个环的周期差体现出来，否则内环频率过低，外环会看到明显的相位滞后，整个系统趋于振荡。

在三菱 PLC 上做自整定参数时，很多工程师习惯用 PLC 自带的 PID 自整定指令，让控制器自己去跑一个继电反馈过程，得到临界增益和临界周期。这个结果可以作为 PSO 的初值，也可以用 PSO 的输出直接覆盖到 PLC 的 PID 参数寄存器中。但需要注意 PLC 里 PID 指令的输入参数往往不是 Kp、Ki、Kd，而是增益、积分时间、微分时间。积分时间Ti = Kp / Ki，微分时间Td = Kd / Kp。换算之后，再确认 PLC 的采样时间是不是和 PSO 仿真里的Ts一致。很多三菱 FX 系列程序里 PID 执行放在定时中断或定时扫描上，扫描周期不固定，直接把仿真参数抄进去之前，至少要用示波器或 D 寄存器画一条实际 PV 曲线验证。

另外，PLC 的压力调节和在加热器上的参数经验值不能互相套用。气压系统响应快，惯性时间常数可能只有几秒；加热系统惯性可能有几十秒甚至几分钟。PSO 的优势恰恰是不怕对象差异大，你只需要修改仿真模型里的增益和时间常数，然后重新跑一遍，而不需要积累大量“经验试凑”的表格。这也是为什么当被控对象时常换线时，离线 PSO 工具比老师傅凭经验记忆更稳定。

5. 先用 Ziegler-Nichols 定范围，再让 PSO 细搜：一个提高收敛质量的小技巧

很多 PID 调参工具书都会花很大篇幅讲 Ziegler-Nichols 整定，但你实际用它调出来的一组参数往往超调偏大。一个很实用的混合方案是：用 Z-N 结果作为 PSO 搜索空间中心，再用 PSO 在这个局部空间里做细搜。这样既能利用 Z-N 的参数估计能力，又能借助 swarm 解决 Z-N 过度振荡的毛病。整个过程如下。

第一步，对被控对象做一个开环阶跃响应测试。记录从输入变化到输出开始变化的时间，即纯滞后L，再记录输出走到 63.2% 稳态值的时间，得到惯性时间常数T。用这个过程模型的 Z-N 公式求初始 PID 参数时，一般用的是反应曲线整定公式：

Kp = 1.2 * T / (K * L) Ti = 2.0 * L Td = 0.5 * L

这里的T、L、K就是从前面实验读出来的参数。代入之后，把 Kp、Ki = Kp/Ti、Kd = Kp * Td 作为 PSO 中心。

第二步，修改 PSO 的上下界。举个例子，如果 Z-N 算出的 Kp 是 1.6，那就把lb=0.8, ub=2.4，让搜索空间缩小到中心值的 50% 到 150%。这个范围足够覆盖模型不确定度，又远比直接在[0, 10]里撒点高效。实测中，这种初始化方式能让算法在不到 20 次迭代内稳定收敛，而完全随机初始化往往要 40 次以上，偶尔还会掉进无意义的超调解。

第三步，在适应度函数中加入超调约束。这里可以用一个简单技巧：把超调量大于阈值后的适应度值直接放大到原来的十倍，迫使粒子逃离振荡区域。因为 Z-N 参数往往超调偏大，PSO 的第一步优化基本就是在保留 Z-N 响应速度的同时压低超调，方向感非常明确。

最后，把得到的结果和初始 Z-N 参数放在同一个仿真器里对比阶跃响应，记录调节时间、超调量和 ITAE 这些可量化指标。如果 PSO 的结果没有显著提升，可能是边界范围设得太紧，把最优解排除在外，这时把上下界扩大到 0.5 倍到 2.5 倍重新跑一遍。这个“Z-N 中心 + PSO 细搜”的做法本质上是把随机搜索和工程先验结合起来，既不放弃 swarm 的全局性，也不让算法在明显不合理的参数区域浪费时间。对刚接触 PSO 调 PID 的人来说，这个技巧比单纯增加粒子数量更有效。

本文还有配套的精品资源，点击获取