pso.m

自己实现的粒子群算法

function [R,g,d]=PSO(x,point)
w=0.9;
c1=2;
c2=2;
Vmax=0.2  
gbest=[10 10 3 3 25];
p=x;
for i=1:20
    for j=1:5
v(i,j)=0;
    end
end
for i=1:500
   
        w=0.5;
  
    for j=1:20
        temp=f(x(j,:),point);
        dis(i,:)=temp;
        if( temp(3)<x(j,5))
            x(j,5)=temp(3);
            p(j,:)=x(j,:);
        end
    end
    gbest=getbest(x,gbest)
    gback(i)=gbest(5);
    for j=1:20
        for k=1:4
          v(j,k)=w*v(j,k)+c1*rand()*(p(j,k)-x(j,k))+c2*rand()*(gbest(k)-x(j,k));
          if (v(j,k)>Vmax)
              v(j,k)=Vmax;
          elseif(v(j,k)<-Vmax)
              v(j,k)=-Vmax;
          end
          
          x(j,k)=x(j,k)+v(j,k);
            
          if(k==1)
              if (x(j,k)>25)
                  x(j,k)=25;
              elseif (x(j,k)<-25)
                  x(j,k)=-25;
              end
            end
            if(k==2)
              if (x(j,k)>25)
                  x(j,k)=25;
              elseif (x(j,k)<-25)
                  x(j,k)=-25;
              end
            end
            if(k==3)
              if (x(j,k)>25)
                  x(j,k)=25;
              elseif (x(j,k)<-25)
                  x(j,k)=-25;
              end
            end
            if(k==4)
              if (x(j,k)>1)
                  x(j,k)=1;
              elseif (x(j,k)<-1)
                  x(j,k)=-1;
              end
            end
            if(k==5)
              if (x(j,k)>1)
                  x(j,k)=1;
              elseif (x(j,k)<-1)
                  x(j,k)=-1;
              end
            end
            if(k==6)
              if (x(j,k)>1)
                  x(j,k)=1;
              elseif (x(j,k)<-1)
                  x(j,k)=-1;
              end
            end
            if(k==7)
              if (x(j,k)>5)
                  x(j,k)=5;
              elseif (x(j,k)<-5)
                  x(j,k)=-5;
              end
            end
            if(k==8)
              if (x(j,k)>5)
                  x(j,k)=5;
              elseif (x(j,k)<-5)
                  x(j,k)=-5;
              end
            end
        end
    end
    if (i>50)
        if (abs(gback(i)-gback(i-50))<0.0000001)
            break
        end
    end
end
R=gbest;
g=gback;
d=dis;

function [g]=getbest(a,gb)
temp=gb;
for i=1:20
    if (temp(5)>a(i,5))
        temp=a(i,:);
    end
end
g=temp;


 function t=f(x,point)
dmax=0;
rmin=100;
n1(1)=x(1);
n1(2)=x(2);
n2(1)=x(3);
n2(2)=x(4);
for i=1:24
    nt(1)=(point(i,2)-n1(2))-point(i,3)*n2(2);
    nt(2)=(point(i,1)-n1(1))-point(i,3)*n2(1);
    nt(3)=(point(i,1)-n1(1))*n2(2)-(point(i,2)-n1(2))*n2(1);
    d=absM(nt);
    d=d/(absM(n2));
    if (d>dmax)      
        dmax=d;
    end
    if (d<rmin)
        rmin=d;
    end
end
t(1)=dmax;
t(2)=rmin;
t(3)=dmax-rmin;


function r=absM(a)
t=sqrt(a(1)*a(1)+a(2)*a(2)+1);
r=t;