这是我个人的程序库中的程序,是本人在大二时为求解电网崩溃节点概率所写。根据每个节点的崩溃概率,求解最大崩溃路径,据此在设计方案时尽可能避免路径的产生。
随即结果与最优结果对比,同时给出相应崩溃概率
MATLAB源程序:
- %% 该函数演示多目标perota优化问题
- %清空环境
- %close all
- %clear,clc
- load qzjz %权值矩阵
- %% 初始参数
- %初始化程序
- Dim=20; %粒子维数
- PopSize=100; %种群个数
- MaxIt=500; %迭代次数
- c1=0.5; %算法参数
- c2=0.7; %算法参数
- wmax=1.2; %惯性因子
- wmin=0.1; %惯性因子
- NumToNext=1; %记录最优适应值个数
- NumToQuit=1; %记录最优适应度重复出现的个数
- w=0;
- for i=1:PopSize %粒子初始化
- Group(i,:)=randperm(20);
- end
- v=zeros(1,Dim); %速度初始化
- %定义位置
- x=[10 10 5 20 30 20 10 20 25 5 20 5 10 10 0 30 15 15 0 5 10];
- y=[2.5 3 2 2 0 1.5 1 0.5 1 0.5 0 0 0.5 0 1.5 0.5 1.5 2 2 1.5 2.5];
- % 粒子适应度值
- IndivdualFitness=ones(1,PopSize); %单个粒子适应度
- %计算初始目标向量
- for i=1:PopSize
- for j=1:Dim-1
- IndivdualFitness(i) = IndivdualFitness(i)*qzjz(Group(i,j),Group(i,j+1)); %粒子概率
- end
- end
- % 最优值初始化
- IndivdualBestFitness=IndivdualFitness; %粒子最大适应度
- IndivdualBestPath=Group; %个体最佳值
- GlobalBestPath=Group(1,:); %粒子群最佳位置
- %记录上一次各粒子的适应值
- OldIndivdualFitness=IndivdualFitness;
- %记录初始最优适应值
- OldGlobalBestFitness=max(IndivdualBestFitness);
- %随机解
- sign={'1','2','3','4','5','6','7','8','9','10','11','12','13','14','15','16','17','18','19','20'};
- figure;
- subplot(2,2,1);
- plot(x,y,'r*');
- for i=1:Dim
- text(x(i)+0.1,y(i)+0.1,sign{i}) %利用点的坐标添加对应标注
- end
- line(x,y);
- title('随机解')
- grid on;
- %% 循环迭代
- for iter=1:MaxIt
-
- % 权值更新
- w=wmax-(wmax-wmin)*iter/MaxIt;
- %% 群体更新
- for i=1:PopSize
- %速度更新
- v=w*v+c1*rand(1,1)*(IndivdualBestPath(i,:)-Group(i,:))+c2*rand(1,1)*(GlobalBestPath-Group(i,:));
-
- %位置更新
- Group(i,:)=Group(i,:)+v;
- Group(i,:)=ceil(Group(i,:));
- index0=[];
- for j=1:Dim
- flag=0; %没有不符合条件的值
- if ((Group(i,j)<=0)||(Group(i,j)>Dim))
- flag=1;
- end
- for k=j+1:Dim
- if Group(i,j)==Group(i,k)
- flag=1;
- break;
- end
- end
- if flag==1
- index0=[index0,j];
- end
- end
-
- t=1;
- m=randperm(Dim);
- NumOfSize=size(index0,2);
- for j=1:Dim
- flag=0; %新产生的随机数与原矩阵中不相同
- for k=1:Dim
- flag0=0;
- for q=1:NumOfSize
- if k==index0(q)
- flag0=1;
- end
- end
- if flag0==1
- continue;
- end
- if m(1,j)==Group(i,k)
- flag=1;
- break;
- end
- end
-
- if flag==0
- Group(i,index0(t))=m(1,j);
- t=t+1;
- end
- if t>length(index0)
- break;
- end
- end
- end
- %% 计算更新后个体适应度
- IndivdualFitness(:)=1;
- for i=1:PopSize
- for j=1:Dim-1 %控制类别
- IndivdualFitness(i) = IndivdualFitness(i)*qzjz(Group(i,j),Group(i,j+1)); %计算粒子i 概率
- end
- end
-
- IsChange=IndivdualFitness>=OldIndivdualFitness; %判断新粒子路径的适应值是否大于原来的
- IndivdualBestPath(find(IsChange),:) = Group(find(IsChange),:);%更新个体最佳路径
- OldIndivdualFitness = OldIndivdualFitness.*( ~IsChange) + IndivdualFitness.*IsChange;%更新个体最佳路径距离
-
- [GlobalBestFitness,Index]=max(IndivdualFitness);
- GlobalBestPath=Group(Index,:);
- if GlobalBestFitness>OldGlobalBestFitness
- NumToNext=NumToNext+1;
- OldGlobalBestFitness(1,NumToNext)=GlobalBestFitness;
- else
- NumToQuit=NumToQuit+1;
- end
- if NumToQuit>=20 %重复20次种群最优适应度为改变就退出迭代
- break;
- end
- end
- %最优解
- for i=1:Dim
- x1(i)=x(GlobalBestPath(i));
- y1(i)=y(GlobalBestPath(i));
- end
- subplot(2,2,2);
- plot(x,y,'r*');
- for i=1:Dim
- text(x(i)+0.1,y(i)+0.1,sign{i}) %利用点的坐标添加对应标注
- end
- line(x1,y1);
- title('最优解');
- grid on;
- %适应值曲线
- subplot(2,2,3:4);
- s=size(OldGlobalBestFitness,2);
- plot(1:s,OldGlobalBestFitness,'p');
- line(1:s,OldGlobalBestFitness);
- %%给出计算结果
- disp('最优路径为:');
- for i=1:Dim
- fprintf('%d',GlobalBestPath(i));
- if i~=Dim
- fprintf('%s','-');
- end
- end
- fprintf('\n');
复制代码
下载:
PSO求解最大崩溃路径.zip
(2.25 KB, 下载次数: 8)
|