matlab遗传算法pid大作业

遗传算法在调节控制系统参数中的应用【摘要】自动化控制系统多采用PID 控制器来调节系统稳定性和动态性，PID 的Kp,Ki,Kd 参数需要合理选择方能达到目标。

遗传算法是一种模拟生物进化寻求最优解的有效算法，本文通过利用GAbx 工具箱实现对控制电机的PID 进行参数优化，利用matlab 的仿真功能可以观察控制效果。

1. 直流伺服电机模型 1.1物理模型图1 直流伺服电机的物理模型αu ---电枢输入电压（V ） a R ---电枢电阻（Ω） S L ---电枢电感（H ） q u ---感应电动势（V ） g T ---电机电磁转矩（N m ⋅） J---转动惯量（2m kg ⋅）B---粘性阻尼系数（s m N ⋅⋅） g i ---流过电枢的电流（A ） θ---电机输出的转角（rad ）1.2传递函数利用基尔霍夫定律和牛顿第二定律得出电机基本方程并进行拉布拉斯变换)()()()()()()()()()()(2s s K s U K s I s T s Bs s Js s T s I s L R s I s U s U e q t a g g a a a a q a θθθ⋅=⋅=⋅+⋅=⋅+⋅=-式中：t K 为电机的转动常数（m N ⋅）A ；e K 为感应电动势常数（s V ⋅）radaa R s L +1S1 BJs +1i KC K)(s U a)(s U q)(s I a )(s T g)(s Ω)(s θ图2 直流伺服电机模型方框图消去中间变量得系统的开环传递函数：sK K B Js R s L K s U s s G C t a d ta ]))([()()()(+++==θ系统参数如下：s m uN B m mg J ⋅⋅=⋅=51.3,23.32A m N K K uH L R e t a a )(03.0,75.2,4⋅===Ω=2. PID 校正图3 PID 校正s K sK K s G d ip c ++=)( Kp,Ki,Kd 为比例，积分，微分系数 令Kp=15、Ki=0.8 、Kd=0.6M 文件：J=3.23E-6;B=3.51E-6; Ra=4;La=2.75E-6; Kt=0.03; num= Kt;den=[(J*La) ((J*Ra)+(La*B)) ((B*Ra)+Kt*Kt) 0]; t=0:0.001:0.2; step(num,den,t); Kp=15; Ki=0.8; Kd=0.6;numcf=[Kd Kp Ki]; dencf=[1 0];numf=conv(numcf,num); denf=conv(dencf,den);[numc,denc]=cloop(numf,denf); t=0:0.001:0.04; step(numc,denc,t);matlab 进行仿真，我们可以看出不恰当的PID 参数并不能使系统达到控制系统的要求，因此需要对PID参数进行优化。

图4 系统阶跃响应（Kp=15，Ki=0.8，Kd=0.6）3.遗传算法3.1 遗传算法和工具箱简介遗传算法（GA)是基于自然选择和基因遗传学原理的优化搜索方法。

它借鉴了达尔文的进化论和孟德尔的遗传学说。

其本质上是一种高效、并行、全局搜索的方法，它能在搜索过程中自动获取和积累有关搜索空间的知识，并自适应地控制搜索过程以求得最优解。

遗传算法操作使用适者生存的原则，在潜在的解决方案种群中逐次产生一个近似最优的方案。

在遗传算法的每一代中，根据个体在问题域中的适应度值和从自然遗传学中借鉴来的再造方法进行个体选择，产生一个新的近似解。

在这个过程导致种群中个体的进化，得到的新个体比原个体更能适应环境，就像自然界中的改造一样。

表3.1遗传学和遗传算法中基本用语对照表遗传学遗传算法染色体（Chromosome) 解的编码（算法的操作对象）基因（Gene）解中每一分量等位基因（Allele）特性值基因座（Locus）二进制串中位置基因型（Genptype）结构表现型（Phenotype）参数集、候选解个体（Individual）解适者生存在算法停止时，最优目标值的解有最大可能被留住适应性（Fitness）适应度函数值群体（Population）选定的一组解复制（Reproduction）根据适应度函数值选取的一组解交配（Crossover）通过交配产生一组新解的过程变异（Mutation）编码的某一个分量发生变化的过程英国谢菲尔德大学开发的遗传算法工具箱把参数，选择，交叉，变异等过程封装成函数进行操作，其基本搜索过程不变。

表3.2 遗传算法工具箱常用函数创建种群crtbase 创建基向量crtbp 创建任意离散随即种群crtrp 创建实值初始种群（bs2rv）适应度计算ranking 常用的基于秩的适应度计算scaling 比率适应度计算选择函数reins一致随机和基于适应度的重插入rws 轮盘选择select高级选择例程sus 随机遍历采样变异算子mut 离散变异mutate高级变异函数mutbga 实值变异交叉算子recdis 离散重组recint 中间重组reclin 线性重组recmut 具有变异特征的线性重组recombine高级重组算子xovdp 两点交叉算子xovdprs 减少代理的两点交叉xovmp 通常多点交叉xovsh 洗牌交叉xovshrs 减少代理的洗牌交叉xovsp 单点交叉xovsprs 减少代理的单点交叉子种群的支持migrate 在子种群间交换个体实用函数bs2rv 二进制串到实值的转换rep 矩阵的复制3.2 利用遗传算法优化过程1）根据遗传算法优缺点这里选择遗传代数为100，种群大小为30，变量维数为3其中10≤Kp≤20,0≤Ki≤1,0≤Kd≤1，要求精度0.0001所以二进制串编码长度为17故Kp 精度为(20-10)/(-1) Ki,Kd精度为（1-0）/(-1),代沟为0.9，交叉概率为0.6，变异概率为0.01。

2）适应度函数的设计PID 优化设计的目的是使系统某些性能指标最优，然而，单纯的误差性能指标很难同时满足系统对快速性、稳定性和鲁棒性的要求，因此在适应度函数中引入超调量、上升时间和累计绝对误差指标项。

设σ=(|y p -y ref |/y ref )×100%为系统的超调量，y p 和y ref 分别为输出峰值和输入参考值;t *r 为上升时间，将其定义为输出从0第1次达到0.95yref 的时间;e k =y k -y ref 为采样时刻k 的输出误差;w j (j=1，2，3)为权重系数。

则多目标适应度函数为:∑=++=nk k r e w t w w f 13*21σ (5)通过对权重系数的调整，可以改变系统对快速性和稳定性的要求。

如系统要求较小的超调，可以适当增大w 1;若系统要求快速的动态响应，则可以适当增大w 2。

计算机控制是一种采样控制，它只能根据采样时刻的偏差值计算控制量。

因此连续PID 控制算法不能直接使用，需要采用离散化方法。

在计算机PID 控制中，使用的是数字PID 控制器。

图3.2 增量式PID 控制系统增量式PID 控制是数字控制器的输出只是控制器的增量△u(k)。

当执行机构需要的是控制量的增量时，应采用增量式PID 控制，采样时刻点kT （T 为采样周期）代表连续时间t 。

做如下近似变换来离散化：()()00()()(1)k tj t kT e t Te j de t e k e k dtT =⎧⎪⎪⎪⎪⎪⎨⎪⎪⎪⎪⎪⎩=≈--≈∑⎰ 可得离散的PID 表达式：1(1)(2)(1)(1)()k p i dj e k e k u k k e k k e j T k T==----=-++∑增量式PID 控制算法：()2(1)(2)()(()(1))()e k e k e k u k kp e k e k kie k T kdT--+-∆=--++增量式PID 控制算法不需多次累加，控制增量△u(k)仅与最近3次得采样有关，所以误动作时的影响相对较小图3.3 基于GA 算法参数调整的PID 控制系统采样时间取1ms ，为获取满意的过渡过程动态特性，采用误差绝对值时间积分性能指标作为参数选择的目标函数。

为防止控制量过大，在目标函数中加入控制输入的平方项。

选用下式作为参数选取的最优指标20(0.999|()|0.001())2uJ e t u t dt t ∞=++⎰()e t 为系统偏差，()u t 为PID 控制器输出，u t为上升时间。

为避免超调了，采用惩罚功能，一旦产生超调，将超调量作为最优指标的一项，此时最优指标为20()0(0.999|()|0.001()100|()|)2uif e t J e t u t e t dt t ∞<=+++⎰J 值越小，对应的PID 参数越优良。

但是，由于GA 算法用于求解最大值问题，所以需要对J 做一个简单变形，从而得到算法的适应度函数：1f J =3）依据遗传算法的步骤编写程序并仿真图3.4 BestKpKiKd优化值时的阶跃响应图3.5 J随进化代数的变化曲线图3.6 优化得到BestKpKiKd及此时的LeastJ4.程序代码1.GA.m文件%GA(Generic Algorithm) Program to optimize Parameters of PIDclcclear all;close all; %清除命令，变量，关闭窗口global rin yout timef %输入，输出，离散时间%******Initialization******%MAXGEN=100; %遗传代数NIND=30; %种群大小NVAR=3; %变量维数PRECI=17; %变量精度GGAP=0.9; %代沟(Generation gap)Pc=0.6; %交叉概率Pm=0.01; %变异概率FieldD=[PRECI,PRECI,PRECI;10,0,0;20,1,1; 1,1,1;0,0,0;1,1,1;1,1,1]; %区域描述器Chrom=crtbp(NIND,NVAR*PRECI); %初始种群for i=1:1:MAXGENgen(i)=i; %遗传代数%****** Step 1 : Evaluate LeastJ ******KpKiKdlsj=bs2rv(Chrom,FieldD); %各代种群由二进制串转化为实值for j=1:1:NINDJ=0;KpKiKd=KpKiKdlsj(j,:);[KpKiKd,J]=pidf(KpKiKd,J); %计算本代种群中各个个体的目标函数值Lsj(j,1)=J;end[OderLsj,IndexLsj]=sort(Lsj); %本代种群中各个个体的目标函数值升序排列genlstj(i)=OderLsj(1); %升序排列中第一个最小，为本代最优BestKpKiKd=KpKiKdlsj(IndexLsj(1),:);%BestKpKiKdLsj=Lsj+1e-10; %Avoiding deviding zero 1*10^(-10)Msj=1./Lsj; %设置适应度函数[OderMsj,IndexMsj]=sort(Msj); %Arranging fi small to bigger适应度值升序排列genmost=OderMsj(NIND) ; % Let Bestfi=max(fi)升序中最后一个最大，为本代最优moststring=Chrom(IndexMsj(NIND),:); %本代中对应最优的PID参数%****** Step 2 : Select and Reproduct Operation******FitnV = ranking(Lsj); %分配适应度值(Assign fitness values)SelCh=select('sus',Chrom,FitnV,GGAP); %选择%************ Step 3 : Crossover Operation ************SelCh=recombin('xovsp', SelCh,Pc); %交叉%************ Step 4: Mutation Operation **************SelCh=mut(SelCh,Pm); %变异SelCh(NIND,:)=moststring; %保证交叉后的新种群Temp2包含原种群中的最优个体%*************Step 5：Reinsert Operation***************ObjV=FitnV;KpKiKdlsj=bs2rv(SelCh,FieldD); %计算ObjVSelfor j=1:1:NINDJ=0;KpKiKd=KpKiKdlsj(j,:);[KpKiKd,J]=pidf(KpKiKd,J);Lsj(j,1)=J;endLsj=Lsj+1e-10;Msj=1./Lsj;ObjVSel=Msj;[Chrom ObjV] = reins(Chrom, SelCh, 1, 1, ObjV,ObjVSel); %重插入endLeastJ=genlstj(MAXGEN) %最优目标函数值（性能指标）BestKpKiKd %最优PID参数（实值）figure(1);plot(gen,genlstj);xlabel('gen(s)');ylabel('genlstj');figure(2);plot(timef,rin,'r',timef,yout,'b');xlabel('Time(s)');ylabel('rin,yout');2.pidf.m函数文件function [KpKiKd,J]=pidf(KpKiKd,J)global rin yout timeftz=0.001;sys=tf(0.03,[8.8825e-12,1.2920e-05,9.1404e-04,0]);zsys=c2d(sys,tz,'z');[num,den]=tfdata(zsys,'v');rin=1.0;u_1=0.0;u_2=0.0;u_3=0.0;y_1=0.0;y_2=0.0;y_3=0.0;K=[0,0,0]';error_1=0;Tup=1;m=0;final=100;for t=1:1:finaltimef(t)=t*tz;r(t)=rin;u(t)=KpKiKd(1)*K(1)+KpKiKd(2)*K(2)+KpKiKd(3)*K(3);if u(t)>=10u(t)=10;endif u(t)<=-10u(t)=-10;endyout(t)=-den(2)*y_1-den(3)*y_2-den(4)*y_3+num(2)*u_1+num(3)*u_2+num(4)*u_3; error(t)=r(t)-yout(t);%------------ Return of PID parameters -------------u_3=u_2;u_2=u_1;u_1=u(t);y_3=y_2;y_2=y_1;y_1=yout(t);K(1)=error(t); % Calculating PK(2)=(error(t)-error_1)/tz; % Calculating DK(3)=K(3)+error(t)*tz; % Calculating Ierror_1=error(t);if m==0if yout(t)>0.95&yout(t)<1.05Tup=timef(t);m=1;endendendfor t=1:1:finalJ=J+0.999*abs(error(t))+0.001*u(t)^2; if t>1erry(t)=yout(t)-yout(t-1);if erry(t)<0J=J+100*abs(erry(t));endendendJ=J+2*Tup;。