m基于PSO粒子群优化的可靠性和费用优化问题求解MATLAB仿真

1.算法概述 PSO是粒子群优化算法（——Particle Swarm Optimization）的英文缩写，是一种基于种群的随机优化技术，由Eberhart和Kennedy于1995年提出。粒子群算法模仿昆虫、兽群、鸟群和鱼群等的群集行为，这些群体按照一种合作的方式寻找食物，群体中的每个成员通过学习它自身的经验和其他成员的经验来不断改变其搜索模式。。这里理论部分和论文《基于配电系统可靠性的开关优化配置研究》完全相同，项目就仿真代码的算法流程做简单的介绍。

传统的PSO如下：

本课题对PSO做了如下两点改进：

第一，引入惯性因子：

那么上面的式子可以表示为：

第二，引入变异因子：

    PSO初始化为一群随机粒子（随机解）。然后通过迭代找到最优解，在每一次迭代中，粒子通过跟踪两个“极值”来更新自己。第一个就是粒子本身所找到的最优解，这个解叫做个体极值。另一个极值是整个种群目前找到的最优解，这个极值是全局机制。另外也可以不用整个种群而只是用其中一部分作为粒子的邻居，那么在所有邻居中的极值就是局部极值。

PSO 过程: 第1 步: 种群随机初始化。 第2 步: 对种群内的每一个个体计算适应值(fitness value) , 适应值与最优解的距离直接有关。 第3 步: 种群根据适应值进行复制。 第4 步: 如果终止条件满足, 则停止; 否则转到第2 步。 PSO 算法中每个优化问题的解都有是搜索空间中的一只鸟, 称为粒子。与其他进化计算技术不同的是群体中的每个粒子可以记忆自己到过的最优位置, 并能感知邻近群体已达到的最优位置, 每个粒子能够根据自身到过的最优位置和邻近群体已到过的最优位置来更新自己, 然后粒子们不断地追随当前的最优粒子在解空间搜索。

2.部分程序 .................................................................. popsize = 10; num = 5; %注意，如下两组参数，在实际中，可有配电网的公式得到，这里重点是对PSO在可靠性满足的情况下进行费用优化，故直接通过参数进行输入 %定义价格初始参数 %F = min(alphaASUI + betaf); %F = [0.025,0.075,0.04,0.035,0.05];

%计算f = CS + CM + CL %CS开关设备投资费用 %CS开关设备投资费用 M = 5; %开关类型的总数 N = [2 10 3 10 24]; %每个类型的开关的安装台数,该参数的设置和M值有关 Cs= [0.0003 0.0004 0.0005 0.0006 0.0005];%每个类型的开关的单价,该参数的设置和M值有关 i = 0.6; %贴现率 P = [2 2 2 2 2]; %第j中开关的使用寿命,该参数的设置和M值有关 for j = 1:M CS_sub(j) = N(j) * Cs(j) * ( ((1+i)^P(j))i )/( ((1+i)^P(j)) - 1 ); end CS = 0; for j = 1:M CS = CS + CS_sub(j);%公式4.3 end %CM %CM Cso = 0.1;%投资费用 r = 0.3; %运行费用所在投资费用的比例 CM = Cso * r; %CL %CL %系统年停电损失费用 LP = num; T = 4; for i = 1:LP for j = 1:T EENS(i,j) = 0.05N(i)rand(); C(i,j) = 0.05rand(); end end for i = 1:LP for j = 1:T CL_sub(i,j) = EENS(i,j)*C(i,j);
end CL(i) = sum(CL_sub(i,:)); end for i = 1:LP f(i) = CS + CM + CL(i); end

%ASUI NT = 1000; for i = 1:LP U(i) = round(100*rand); ASUI(i) = U(i)N(i)100 / NT/8760; end alpha = 0.5; beta = 0.5; for i = 1:LP F(i) = alphaASUI(i) + betaf(i); end

%在不同价格指数下的可靠性指标 Reliability =[0.95,0.92,0.75,0.88,0.83]; F = F'; r0 = 0.9; Pop_sub = zeros(popsize,num+2);
New_pop_sub = zeros(popsize,num+2);
PSO_V = zeros(popsize,num); %在论文中，这两个参数定义为2.05 c1 = 2.05;%加速度常数 c2 = 2.05;%加速度常数 MaxIteration= 50;%最大迭代次数 %最大最小速度 V_max = 1; V_min = -1; %最大最小个体 P_max = 3; P_min = 0 ; %% 粒子参数初始化 %粒子参数初始化，且初始化得满足可靠性指标 [BEST,Pop_sub,New_pop_sub] = PSO_initial(Pop_sub,New_pop_sub,popsize,num,Reliability,F,r0);

%速度初始化 for i=1:popsize for j=1:num PSO_V(i,j)=rands(1); end end

%% 主程序，通过PSO进行费用优化 ws = zeros(MaxIteration,1);

for k=1:MaxIteration k ws(k) = 0.9-0.5*k/MaxIteration;%惯性权重,对应论文公式4.17

%最大值计算
for i=1:popsize                                                      
    if Pop_sub(i,num+2) < New_pop_sub(i,num+2)  
       New_pop_sub(i,:) = Pop_sub(i,:);  
    end
    
    if New_pop_sub(i,num+2) < BEST(num+2)  
       BEST = New_pop_sub(i,:);  
    end
end
    
%将计算状态进行保存
PSO_state(k,1) = k;

PSO_state(k,2) = BEST(num+2);

PSO_state(k,3) = min(Pop_sub(:,num+2));

PSO_state(k,4) = max(Pop_sub(:,num+2));

%更新位置与速度  
for i=1:popsize
    for j=1:num
        %更新速度,对应论文公式4.15
        PSO_V(i,j)     = ws(k)*PSO_V(i,j) + c1*rand*(New_pop_sub(i,j)-Pop_sub(i,j)) + c2*rand*(BEST(j)-Pop_sub(i,j)); 
        
        %对速度进行限速
        if PSO_V(i,j) > V_max
           PSO_V(i,j) = V_max;%小于最大速度
        end
        if PSO_V(i,j) < V_min
           PSO_V(i,j) = V_min;%大于最小速度    
        end            
   
        
        %更新位置,对应论文公式4.16
        Pop_sub(i,j) = Pop_sub(i,j)+PSO_V(i,j);                           
        
        %粒子变异
        tmp{k} = Pop_sub(i,j);
        if k > 1
           %随机选择两个点作为父母点
           p = rand(1,1);            
           Pop_sub(i,j) = p*tmp{k-1} + (1-p)*tmp{k-1};
        end
        
        %由父母例子产生子例子，并带有5%的变异概率
        if rand>0.95
           %变异
           Pop_sub(i,j) = 0.3*rand;
        end           
        

        Pop_sub(i,j)=(Pop_sub(i,j));
        

        if Pop_sub(i,j) > V_max
           Pop_sub(i,j) = V_max;
        end
        if Pop_sub(i,j) < V_min
           Pop_sub(i,j) = V_min;  
        end           
        
    end
end

 %保证可靠性的要求，同时计算价格
 for i=1:popsize,
     
     Pop_sub(i,num+1)=1;
     for j=1:num
         Pop_sub(i,num+1)=Pop_sub(i,num+1)*(1-(1-Reliability (j))^Pop_sub(i,j));
     end
     %保证可靠性要求
     if Pop_sub(i,num+1)<r0  
         Pop_sub(i,:)=func_Reliability(Pop_sub(i,:),F,Reliability ,r0,num);  
     end   
     
     Pop_sub(i,num+2)=Pop_sub(i,1:5)*F;
 end

end

......................................................... 02-005m

3.算法部分仿真结果图