基于遗传算法的PID参数整定的MATLAB程序代码.docx
《基于遗传算法的PID参数整定的MATLAB程序代码.docx》由会员分享，可在线阅读，更多相关《基于遗传算法的PID参数整定的MATLAB程序代码.docx（19页珍藏版）》请在冰点文库上搜索。
基于遗传算法的PID参数整定的MATLAB程序代码
基于遗传算法的PID参数整定
1引言
传统的比例、积分、微分控制，即PID控制具有算法简单、鲁棒性好和可靠性高等优点，已经被广泛用于工业生产过程。
但工程实际中，PID控制器的比例、积分和微分调节参数往往采用实验加试凑的方法由人工整定。
这不仅需要熟练的技巧，往往还相当费时。
更为重要的是，当被控对象特性发生变化，需要控制器参数作相应调整时，PID控制器没有自适应能力，只能依靠人工重新整定参数，由于经验缺乏，整定结果往往达不到最优值，难以满足实际控制的要求。
考虑生产过程的连续性以及参数整定费事费力，这种整定实际很难进行。
所以，人们从工业生产实际需要出发，基于常规PID控制器的基本原理，对其进行了各种各样的改进。
近年来许多学者提出了基于各种智能算法的PID整定策略，比如模糊PID、神经元网络PID等等。
然而，这些先进算法都要求对被控对象有很多的先验知识，在实际应用中往往难于做到。
随着计算技术的发展，遗传算法有了很大的发展。
将遗传算法用于控制器参数整定，已成为遗传算法的重要应用之一。
本文介绍基于遗传算法的PID参数整定设计方法。
这是一种寻求全局最优的控制器优化方法，且无需对目标函数微分，可提高参数优化效果，简化计算过程。
仿真实例表明该方法与其他传统寻优方法相比，在优化效果上具有一定的优势。
2遗传算法简介
2.1遗传算法的基本原理
遗传算法是JohnH．Holland根据生物进化的模型提出的一种优化算法。
自然选择学说是进化论的中心内容。
根据进化论，生物的发展进化主要有三个原因：
即遗传、变异和选择。
遗传算法基于自然选择和基因遗传学原理的搜索方法，将“优胜劣汰，适者生存”的生物进化原理引入待优化参数形成的编码串群体中，按照一定的适应度函数及一系列遗传操作对各个体进行筛选，从而使适应度高的个体被保留下来，组成新的群体；新群体包含上一代的大量信息，并且引入了新的优于上一代的个体。
这样周而复始，群体中各个体适应度不断提高，直至满足一定的极限条件。
此时，群体中适应度最高的个体即为待优化问题的最优解。
遗传算法通过对参数空间编码并用随机选择作为工具来引导搜索过程朝着更高效的方向发展。
正是由于遗传算法独特的工作原理，使它能够在复杂空间进行全局优化搜索，具有较强的鲁棒性。
另外，遗传算法对于搜索空问，基本上不需要什么限制性的假设（如连续、可微及单峰等）。
而其它优化算法，如解析法，往往只能得到局部最优解而非全局最优解，且需要目标函数连续光滑及可微；枚举法虽然克服了这些缺点，但计算效率太低，对于一个实际问题常由于搜索空间太大而不能将所有情况都搜索到；即使很著名的动态规划法，也遇到“指数爆炸”问题，对于中等规模和适度复杂性的问题常常无能为力。
2.2遗传算法的特点
同常规优化算法相比，遗传算法有以下特点：
1）遗传算法是对参数的编码进行操作，而不对参数本身。
首先基于一个有限的字母表，把最优化问题的自然参数集编码为有限长度的字符串。
2）遗传算法是从许多点开始并行操作的，而不局限于一点，有效防止搜索过程收敛于局部最优解。
3）遗传算法通过目标函数来计算适应度，不需要其他推导和附加信息，对问题的依赖性较小。
4）遗传算法的寻优规则是由概率决定的，而非确定性的。
5）遗传算法在解空间进行高效启发式搜索，而非盲目地穷举或完全随机搜索。
6）遗传算法对于待寻优的函数基本无限制，它不要求函数连续和可微，可以是数学解析式所表达的显函数，又可以是映射矩阵甚至是神经网络等隐函数，因而应用范围较广。
7）遗传算法具有并行计算的特点，因而可以通过大规模并行计算来提高计算速度。
8）遗传算法计算简单，功能强，更适合大规模复杂问题的优化。
3遗传算法的基本操作
3.1复制
复制是从一个旧种群中选择生命力强的个体位串产生新种群的过程。
根据位串的适配值拷贝，也就是指具有高适配值的位串更有可能在下一代中产生一个或多个子孙。
它模仿了自然现象，应用了达尔文的适者生存理论。
复制操作可以通过随机方法来实现。
若用计算机程序来实现，可考虑首先产生0~1之间均匀分布的随机数，若某串的复制概率为40%，则产生的随机数在0~0.40之间时，该串被复制，否则被淘汰。
此外，还可以通过计算方法实现，其中较典型的几种方法为适应度比例法、期望值法、排位次序法等。
适应度比例法较常用。
3.2交叉
复制操作能从旧种群中选择出优秀者，但不能创造新的染色体。
而交叉模拟了生物进化过程中的繁殖现象，通过两个染色体的交换组合，来产生新的优良品种。
它的过程为：
在匹配池中任选两个染色体，随机选择一点或多点交换点位置；交换双亲染色体交换点右边的部分，即可得到两个新的染色体数字串。
交换体现了自然界中信息交换的思想。
交叉有一点交叉、多点交叉、还有一致交叉、顺序交叉和周期交叉。
一点交叉是最基本的方法，应用较广，它是指染色体切断点有一处，例如：
A:
1011001110->1011000101
B:
0010100101->0010101110
3.3变异
变异运算用来模拟生物在自然界的遗传环境中由于各种偶然因素的基因突变，它以很小的概率随机地改变遗传基因（表示染色体的符号串的某一位）的值。
在染色体以二进制编码的系统中，它随机地将染色体的某一基因由1变为0，或由0变为1。
若只有选择和交叉，而没有变异，则无法在初始基因组合以外的空间进行搜索，使进化过程在早期就陷入局部解而进入终止过程，从而影响解的质量。
为了在尽可能大的空间中获得质量较高的优化解，必须采用变异操作。
4基于遗传算法的PID控制设计
4.1遗传算法整定PID参数的优点
采用遗传算法进行PID三个系数的整定，具有以下优点：
1）与单纯形法相比，遗传算法同样具有良好的寻优特性，且它克服了单纯形法的初值敏感性。
在初始条件选择不当的情况下，遗传算法在不需要给出调节器初始参数的情况下，仍能寻找到合适的参数，使控制目标满足要求。
同时单纯形法难以解决多值函数问题以及在多参数寻优（如串级系统）中，容易造成寻优失败或时间过长，而遗传算法的特性决定了它能很好地克服以上问题。
2）与专家整定法相比，它具有操作方便，速度快的优点，不需要复杂的规则，只需要对字符串进行复制、交叉、变异运算，便可达到寻优。
避免了专家整定法中前期大量的知识库整理工作以及大量的仿真实验。
3）遗传算法是从许多点开始并行操作，在解空间进行高效启发式搜索，克服了从单点出发的弊端以及搜索的盲目性，从而提高运算速度，避免了过早陷入局部最优解（如图1所示）。
4）遗传算法不仅适用于单目标寻优，而且也适用于多目标寻优。
根据不同的控制系统，针对一个或多个目标，遗传算法均能在规定的范围内寻找到合适参数。
遗传算法作为一种全局优化算法，得到了越来越广泛的应用。
近年来，遗传算法在控制上的应用日益增多。
图1遗传算法的运算过程
4.2基于遗传算法整定PID参数的步骤
（1）确定决策变量和约束条件
我们知道，PID控制的三个系数为Kp、Ki、Kd，我们可以根据参数的物理意义和经验来确定它们的取值范围，即约束条件。
（2）建立优化模型
为获取满意的过渡过程动态特性，采用误差绝对值时间积分性能指标作为参数选择最小目标函数。
为防止控制量过大，在目标函数中加入控制输入的平方项。
最后选用以下公式作为参数选择的最优指标：
J=
式中，有e（t）为系统误差，u（t）为控制器输出，
为上升时间，
，
，
为加权值。
同时，为避免超调，采用遗传算法的处罚功能，即一旦超调，将超调量作为最优指标的一项，此时最优指标为：
如果e（t）<0则
J=
其中，
为加权值，且
»
。
（3）确定编码和解码方法。
遗传算法有二进制、格雷码、浮点数和实数编码方法。
用实数编码无需解码，但进行遗传操作时不太方便;用二进制编码方法，遗传操作方便，解码处理后即可得到最优解。
（4）确定个体评价方法，即确定个体适应度函数（FitnessFunction）。
评价个体适应度的一般过程是：
对个体编码串进行解码处理后，可得到个体的表现型;由表现型可计算出对应个体的目标函数值;根据最优问题的类型，有目标函数值按一定的转换规则可求出个体的适应度。
（5）确定遗传算法的运行参数
根据实际情况确定群体大小M、遗传代数kg、交叉概率Pc、变异概率Pm和权值
，
，
，
的大小。
本文取M=30，kg=100，Pc=0.9，Pm不固定，
=0.999
=0.001，
=2，
=100。
编写遗传算法程序，并进行仿真，从而得到优化后的PID参数。
5遗传算法整定PID参数的编程实现与仿真
5.1工业对象的数学建模
由于工业领域中的被控对象一般为一阶或二阶环节，因此，在本文里我们拟定受控对象的传递函数为如下：
G（s）=
其中采样时间为20s。
5.2遗传算法整定PID参数的实现与控制效果
5.2.1编写最优指标的M文件
利用遗传算法整定PID参数，首先，最核心的是要编写最优指标的m文件。
m文件里主要是一个函数，函数里包含了受控对象以及性能指标的程序实现。
它是联系遗传算法与PID整定环节的桥梁。
具体程序见附件。
根据4.2的步骤
（2）可以得出目标函数J如下：
J=
e（t）
0
J=
e（t）<0
MATLAB程序实现：
fori=1:
1:
P
Ji（i）=0.999*abs（error（i））+0.001*u（i）^2;
B=B+Ji（i）;
ifi>1
erry（i）=yout（i）-yout（i-1）;
iferry（i）<0
B=B+100*abs（erry（i））;
end
end
end
BsJ=B+2*tu;
在控制系统里要求目标函数的值越小越好。
然而，在遗传算法里是按照最大化方向进行搜索的。
因此，设计适应度函数为：
F=1/J.主程序的编程实现也比较简单，具体如下：
Ji=BsJi+1e-10;%避免分母为零
fi=1./Ji;
在MATLAB程序里用此句来调用函数：
[Kpidi,BsJ]=fitness（Kpidi,BsJ）;因此，在写完m文件时，应把文件命名为fitness.m。
另外一个也是非常重要的，那就是传递函数的生成与定义输出，实现如下：
由于传递函数为：
G（s）=
Z变换后可得G（z）=
程序实现：
ts=20;%采样时间20s
sys=tf（1,[60,1]，’inputdelay’,80）;%生成传递函数G（s）
dsys=c2d（sys,ts,'z'）;%Z变换G（z）
[num,den]=tfdata（dsys,'v'）;%返回z变换后的各项系数
其中：
num
（1）=0num
（2）=0.2835
den
（1）=1den
（2）=-0.7165
=G（z）=
，从而
U（z）*0.2835=Y（z）（z-0.7165）
即U（z）*num
（2）=Y（z）（z+den
（2））
写成差分方程为：
y（k+1）=-den
（2）y（k）+num
（2）u（k）
差分方程在程序可以如下实现：
yout（k）=-den
（2）*y_1+num
（2）*u_1;
u_1=u（k）;
y_1=yout（k）;
5.2.1编写主程序
主程序里是遗传算法的实现，其中最重要的是三个遗传操作：
复制、交叉、变异。
遗传算法的程序流程图如图2所示：
图2遗传算法的程序流程图
1）复制
根据第二章的介绍，编写复制的部分程序如下：
fi_sum=sum（fi）;
fi_Size=（Oderfi/fi_sum）*Size;
fi_S=floor（fi_Size）;
r=Size-sum（fi_S）;
Rest=fi_Size-fi_S;
[RestValue,Index]=sort（Rest）;
fori=Size:
-1:
Size-r+1
fi_S（Index（i））=fi_S（Index（i））+1;
end
k=1;
fori=Size:
-1:
1
forj=1:
1:
fi_S（i）
TempE（k,:
）=Kpid（Indexfi（i）,:
）;
k=k+1;
end
end
2）交叉
交叉概率在本文里设定为0.9，交叉的主要程序如下：
Pc=0.90;%交叉概率
fori=1:
2:
（Size-1）
temp=rand;
ifPc>temp
alfa=rand;
TempE（i,:
）=alfa*Kpid（i+1,:
）+（1-alfa）*Kpid（i,:
）;
TempE（i+1,:
）=alfa*Kpid（i,:
）+（1-alfa）*Kpid（i+1,:
）;
end
end
TempE（Size,:
）=BestS;
Kpid=TempE;
3）变异
变异部分程序的实现如下：
Pm=0.10-[1:
1:
Size]*（0.01）/Size;
Pm_rand=rand（Size,CodeL）;
Mean=（MaxX+MinX）/2;
Dif=（MaxX-MinX）;
fori=1:
1:
Size
forj=1:
1:
CodeL
ifPm（i）>Pm_rand（i,j）
TempE（i,j）=Mean（j）+Dif（j）*（rand-0.5）;
end
end
end
4）定义输出
本程序的定义了两个输出，一个是目标函数值J;另一个是PID参数整定后的具体值。
Bestfi
BestS
Best_J=BestJ（G）
figure
（1）;
plot（time,BestJ）;
xlabel（'Times'）;ylabel（'BestJ'）;
figure
（2）;
plot（timef,rin，'r',timef,yout,'b'）;
xlabel（'Time（s）'）;ylabel（'rin,yout1,yout'）;主程序的完整程序见附件。
5.2.2仿真结果与分析
根据遗传算法优化PID参数的步骤,运行程序，经过100代的进化，可获得优化参数Kp=1.0814、Ki=0.9093、Kd=0.0329，性能指标J=218.3443。
性能指标的优化过程和经过优化后的PID控制阶跃响应分别如图3和图4所示：
图3性能指标J的优化过程
图4整定后的PID阶跃响应
由上图可知，基于遗传算法整定PID的阶跃响应曲线能保证系统的稳定，调节时间较短，曲线能较快速地到达稳态，但是超调量较大。
6程序
（1）遗传算法程序
clearall;%清除了所有的变量，包括全局变量global
closeall;
globalrinyouttimef%定义全局变量
Size=30;%样本个数定为30
CodeL=3;
MinX
（1）=zeros
（1）;
MaxX
（1）=20*ones
（1）;%定义Kp的取值范围为0-20
MinX
（2）=zeros
（1）;
MaxX
（2）=1.0*ones
（1）;%定义Ki的取值范围为0-1
MinX（3）=zeros
（1）;
MaxX（3）=1.0*ones
（1）;%定义Kd的取值范围为0-1
Kpid（:
1）=MinX
（1）+（MaxX
（1）-MinX
（1））*rand（Size,1）;%取值为0-1的随机数，产生初始种群
Kpid（:
2）=MinX
（2）+（MaxX
（2）-MinX
（2））*rand（Size,1）;
Kpid（:
3）=MinX（3）+（MaxX（3）-MinX（3））*rand（Size,1）;
G=100;%100代进化
BsJ=0;%最优指标初始值
%***************开始运行***************
forkg=1:
1:
G
time（kg）=kg;
%******Step1:
估计最优指标******
fori=1:
1:
Size
Kpidi=Kpid（i,:
）;
[Kpidi,BsJ]=fitness（Kpidi,BsJ）;%调用m函数
BsJi（i）=BsJ;
end
[OderJi,IndexJi]=sort（BsJi）;
BestJ（kg）=OderJi
（1）;
BJ=BestJ（kg）;
Ji=BsJi+1e-10;%避免分母为零
fi=1./Ji;%转换为适应度函数
[Oderfi,Indexfi]=sort（fi）;%将适应度值从小到大排列
Bestfi=Oderfi（Size）;%Bestfi取得适应度最大值
BestS=Kpid（Indexfi（Size）,:
）;%LetBestS=E（m）,mistheIndexfibelongtomax（fi）
kg%输出
BJ
BestS
%******Step2:
选择与复制******
fi_sum=sum（fi）;
fi_Size=（Oderfi/fi_sum）*Size;
fi_S=floor（fi_Size）;%SelectingBiggerfivalue，floor（）朝负无穷方向舍入
r=Size-sum（fi_S）;
Rest=fi_Size-fi_S;
[RestValue,Index]=sort（Rest）;
fori=Size:
-1:
Size-r+1
fi_S（Index（i））=fi_S（Index（i））+1;
end
k=1;
fori=Size:
-1:
1
forj=1:
1:
fi_S（i）
TempE（k,:
）=Kpid（Indexfi（i）,:
）;%选择与复制
k=k+1;%k用于复制
end
end
%************Step3:
交叉************
Pc=0.90;%交叉概率
fori=1:
2:
（Size-1）
temp=rand;
ifPc>temp%交叉条件
alfa=rand;
TempE（i,:
）=alfa*Kpid（i+1,:
）+（1-alfa）*Kpid（i,:
）;
TempE（i+1,:
）=alfa*Kpid（i,:
）+（1-alfa）*Kpid（i+1,:
）;
end
end
TempE（Size,:
）=BestS;
Kpid=TempE;
%************Step4:
变异**************
Pm=0.10-[1:
1:
Size]*（0.01）/Size;%Pm应取很小值
Pm_rand=rand（Size,CodeL）;
Mean=（MaxX+MinX）/2;
Dif=（MaxX-MinX）;
fori=1:
1:
Size
forj=1:
1:
CodeL
ifPm（i）>Pm_rand（i,j）%变异条件
TempE（i,j）=Mean（j）+Dif（j）*（rand-0.5）;
end
end
end
TempE（Size,:
）=BestS;
Kpid=TempE;
end
Bestfi
BestS
Best_J=BestJ（G）
figure
（1）;
plot（time,BestJ）;
xlabel（'Times'）;ylabel（'BestJ'）;
figure
（2）;
plot（timef,rin,'r',timef,yout,'b'）;
xlabel（'Time（s）'）;ylabel（'rin,yout1,yout'）;
（2）适应度函数程序（fitness.m）：
function[Kpidi,BsJ]=pid_gaf（Kpidi,BsJ）
globalrinyouttimef
ts=20;%采样时间为20s
sys=tf（1,[60,1],’inputdelay’,80）;%传递函数公式
dsys=c2d（sys,ts,'z'）;
[num,den]=tfdata（dsys,'v'）;
rin=1.0;
u_1=0.0;u_2=0.0;
y_1=0.0;y_2=0.0;
x=[0,0,0]';
B=0;
error_1=0;
tu=1;
s=0;
P=100;
fork=1:
1:
P
timef（k）=k*ts;
r（k）=rin;
u（k）=Kpidi
（1）*x
（1）+Kpidi
（2）*x
（2）+Kpidi（3）*x（3）;
ifu（k）>=10%防止积分饱和
u（k）=10;
end
ifu（k）<=-10
u（k）=-10;
end
yout（k）=-den
（2）*y_1+num
（2）*u_1;
error（k）=r（k）-yout（k）;
%------------ReturnofPIDparameters-------------
u_1=u（k）;
y_1=yout（k）;
x
（1）=error（k）;%计算P
x
（2）=（error（k）-error_1）/ts;%计算D
x（3）=x（3）+error（k）*ts;%计算I
error_2=error_1;
error_1=error（k）;
ifs==0
ifyout（k）>0.95&yout（k）<1.05
tu=timef（k）;
s=1;
end
end
end
fori=1:
1:
P%目标函数J公式实现
Ji（i）=0.999*abs（error（i））+0.001*u（i）^2;
B=B+Ji（i）;
ifi>1
erry（i）=yout（i）-yout（i-1）;
iferry（i）<0
B=B+100*abs（erry（i））;
end
end
end
BsJ=B+2*tu;
升级会员 