pid_lianxi_descrete.m

这是个利用粒子群算法(PSO)优化PID控制器参数的matlab程序

clear;
x0=[2.813,1.719,1.151];%x=[Kp,Ki,Kd]
c1=1.495;c2=1.495;n=3;group=50;Dmax=100;
% Vmm=[0.4 1.2;0.35 1.05;0.05 0.15];
Vmm=[-2.5 2.5;-2.5 2.5;-2.5 2.5];%可以为负吗

Xmm=[0 5;0 5;0 5];
X=zeros(n,group);V=zeros(n,group);
% omega=zeros(1,group);
% e=0.001;

X(1,:)=x0(1)+5*rand(1,group);V(1,:)=Vmm(1,1)+5*rand(1,group);
X(2,:)=x0(2)+5*rand(1,group);V(2,:)=Vmm(2,1)+5*rand(1,group);
X(3,:)=x0(3)+5*rand(1,group);V(3,:)=Vmm(3,1)+5*rand(1,group);
P_p=X;globe=zeros(n,Dmax);

%%%评价每个粒子适应值，寻找出 P_globle
for j=1:group
    xx=X(:,j)';
%     fz(j)=J_ITAE(xx);
    fz(j)=J_ITAE_discrete(xx);
end
[P_g,I]=min(fz);%P_g  1*1 ?
globe(:,1)=X(:,I);

for itrtn=1:Dmax
    omega=0.739;
    r1=rand(1);r2=rand(1);
    for i=1:group
        V(:,i)=omega*V(:,i)+c1*r1*(P_p(:,i)-X(:,i))+c2*r2*(globe(itrtn)-X(:,i));
    end
    X=X+V;%先虑出位置，如果先是速度虑出，则位置虑出的受速度影响。
    for i=1:group
    %compare each dimension of V
        for row=1:n
            if X(row,i)>Xmm(row,2)
                X(row,i)=Xmm(row,2);
            elseif X(row,i)<Xmm(row,1)
                X(row,i)=Xmm(row,1);
            else
            end
        end
    end
    for i=1:group
        %compare each dimension of V
        for row=1:n
            if V(row,i)>Vmm(row,2)
                V(row,i)=Vmm(row,2);
            elseif V(row,i)<Vmm(row,1)
                V(row,i)=Vmm(row,1);
            else
            end
        end
    end
    
    for k=1:group
        xx=X(:,j)';
%         fz1(k)=J_ITAE(xx);
        fz1(j)=J_ITAE_discrete(xx);
        if fz1(k)<fz(k)
            P_p(:,k)=X(:,k);
            fz(k)=fz1(k);
        end
        if fz(k)<P_g
            P_g=fz(k);
            
        end
    end
    [P_g1 I]=min(fz);
    Ess(itrtn)=P_g1;
%     globe=P_p(:,I);
    globe(:,itrtn+1)=X(:,I)
%     if P_g1<e    
%         break;
%     end
end
kp=globe(1,100);ti=globe(2,100)/globe(1,100);td=globe(3,100)/globe(2,100);T=0.5;
   globe_0=[1.719,2.813,1.151];
kp_0=globe_0(1);ti_0=globe_0(2)/globe_0(1);td_0=globe_0(3)/globe_0(2);

numpid=[td*ti,ti,1];denpid=[0,ti/kp,0];G_pid=tf(numpid,denpid);
  numpid_0=[td_0*ti_0,ti_0,1];denpid_0=[0,ti_0/kp_0,0];G_pid_0=tf(numpid_0,denpid_0);

[numz,denz]=pade(T,4);G_exp=tf(numz,denz);

numd=([0,0,1]);dend=([1,2,1]);G_d=tf(numd,dend);


G_c=feedback(G_d*G_pid,G_exp);step(G_c,'b');hold on
    G_c_0=feedback(G_d*G_pid_0,G_exp);step(G_c_0,'r') ;