遗传算法自整定 PID 参数（Simulink 实现）

PID参数的整定在实际工程里大多靠手动调节，费时费力，效果也依赖经验。今天整理一下用遗传算法（GA）来自动整定PID参数的思路和实现过程，主要在MATLAB/Simulink平台上完成。

需要先说明的是：遗传算法整定PID在实际工程中并不算实用——计算量较大，MCU往往扛不住——但作为一种研究和学习方法还是有一定价值的。

搭建Simulink仿真平台

首先搭建PID仿真的Simulink模型。

和普通的PID仿真模型相比，差别主要在于引入了GA函数来优化PID参数，因此需要一个输出口将仿真结果传出，同时需要一个输入口将GA函数得到的种群传进来。输入输出均链接到workspace，分别用simout和simin模块实现；simin的参数还需要通过总线单元拆分成三路，分别接到P、I、D端口。模型里的常数模块（值为100）用作滤波系数。

优化目标函数的设计

参考论文资料，选取如下优化目标函数：

其中e(t)为系统误差，u(t)为控制器输出（阶跃仿真中即为单位阶跃1(t)），tr为上升时间，w1、w2、w3为权值。为了抑制超调，引入惩罚机制：一旦产生超调，就把超调量作为目标函数的一项单独计入。

目标函数对应的MATLAB函数如下，仿真取参数w1=1.2，w2=0.001，w3=4.0，w4=100：

function aim = fcn(u)

time=(0:0.001:3);  % 时间序列：采样时间0.001秒，仿真时间3秒

output=u;

w1=1.2;
w2=0.001;
w3=4.0;
w4=100;

% 获取上升时间
i1=1;
while(output(i1)<0.1)
    i1=i1+1;
end
t1=time(i1);

i2=1;
while(output(i2)<0.9)
    i2=i2+1;
end
t2=time(i2);

tr=t2-t1;

% 获取积分量
int1=ones(1,3001);
int2=ones(1,3001);
int3=ones(1,3001);

for i=(1:3001)
    if output(i)<1
        int1(i)=0.001*abs(output(i)-1);
        int2(i)=0.001;
        int3(i)=0;
    else
        int1(i)=0;
        int2(i)=0.001;
        int3(i)=0.001*abs(output(i)-1);
    end
end

% 计算目标函数值
aim = w1*sum(int1)+w2*sum(int2)+w4*sum(int3)+w3*tr;

上升时间tr的计算采用10%到90%的标准定义。积分部分按是否超调分两段处理：输出未超过目标值时只计误差绝对值积分，一旦超调则将超调量单独计入并乘以较大权值w4=100，起到惩罚作用。

GA与Simulink的接口函数

有了目标函数之后，还需要一个整体的GA接口函数，输入为kp、ki、kd的种群，输出为每个个体对应的aim值。该函数本质上是一个隐式函数，内部调用Simulink模型完成仿真，再将仿真结果传入fcn()换算成aim值：

function z=ga(kp,ki,kd)

t=[0:0.001:3]';

% 定义simin结构体并声明为全局变量
field1='time';value1=ones(3001,1);
field3='values';value3=ones(3001,3);
field4='dimensions';value4=3;
field2='signals';value2=struct(field3,value3,field4,value4);

global simin;
simin=struct(field1,value1,field2,value2);

% 对40个个体逐一仿真，计算对应的aim值
for i=1:40
    gain=[kp(i),ki(i),kd(i)];
    gain_cal = [kp(i)*ones(3001,1),ki(i)*ones(3001,1),kd(i)*ones(3001,1)];

    simin.time = t;
    simin.signals.values = gain_cal;
    simin.signals.dimensions = 3;

    sim pid_test;

    z(i,1)=fcn(simout);
end

每个个体都要单独跑一次Simulink仿真，40个个体就是40次，这也是实际工程中难以在线使用的主要原因。

主函数：遗传算法迭代

最后是主函数，负责完成种群初始化、进化迭代（选择、交叉、变异、重插入）以及结果可视化：

global simin;

global NIND;  % 定义种群数量为NIND，初始种群40个个体
NIND=40;

NVAR=3;      % 变量维数为3
MAXGEN=100;  % 最大进化代数为100代
PRECI=30;    % 基因编码的二进制位数
GGAP=1;      % 代沟为1，种群个体总数不变

trace=zeros(2,MAXGEN);  % 算法性能追踪矩阵的初始化

Chrom=crtbp(NIND,PRECI*NVAR);  % 种群的初始化

FieldD=[PRECI PRECI PRECI;10 0 0;100 10 10;1 1 1;0 0 0;1 1 1;1 1 1];
% 种群进化域，其中10 0 0;100 10 10表示三个变量的上下限分别为[10,100],[0,10],[0,10]
% 后面四行分别为：1-标准二进制编码，0-标准刻度，1-下限闭区间，1-上限闭区间

gen=0;  % 初始代数为0

variable=bs2rv(Chrom,FieldD);  % 二进制种群转换为十进制自变量
ObjV=ga(variable(:,1),variable(:,2),variable(:,3));  % 计算初始种群的aim值

while gen<MAXGEN

    FitnV=ranking(ObjV);                    % 适应度排序，取最小值直接用ObjV
    SelCh=select('sus',Chrom,FitnV,GGAP);  % 选择运算
    SelCh=recombin('xovsp',SelCh,0.7);     % 交叉运算
    SelCh=mut(SelCh);                       % 变异运算

    variable=bs2rv(SelCh,FieldD);           % 二进制转换成十进制自变量
    ObjVSel=ga(variable(:,1),variable(:,2),variable(:,3));  % 计算进化后的aim值

    [Chrom, ObjV]=reins(Chrom,SelCh,1,1,ObjV,ObjVSel);    % 重插入运算

    gen=gen+1;

    [Y,I]=min(ObjV);hold on;  % 找出进化后种群里面的最小值

    trace(1,gen)=min(ObjV);               % 将最小值放入算法性能追踪矩阵
    trace(2,gen)=sum(ObjV)/length(ObjV);  % 种群平均aim值放入追踪矩阵

end

variable=bs2rv(Chrom,FieldD);

figure(1);
plot(variable(:,1),'o');grid

figure(2);
plot(trace(1,:));
hold on;
plot(trace(2,:),'-.');grid
legend('解的变化','种群均值的变化')  % 将aim值的变化和种群均值变化可视化

FieldD里设定的搜索范围是kp∈[10,100]、ki∈[0,10]、kd∈[0,10]，可以根据具体被控对象调整。整个流程跑完之后，trace矩阵记录了每一代的最优值和种群均值，从图上可以直观看出算法收敛情况。

整体来说，这套实现思路比较清晰：Simulink负责仿真，MATLAB端的GA负责寻优，两者通过全局变量simin/simout交换数据。主要的局限就是每次评估适应度都要调用一次完整仿真，速度上比较慢，后续如果有需要可以考虑用解析模型替代仿真来加快计算。