26号粒子群

赧然的贼158

已于 2024-08-07 13:58:04 修改

阅读量202

点赞数 2

文章标签：人工智能算法

于 2024-07-30 13:51:50 首次发布

本文链接：https://blog.csdn.net/weixin_65744632/article/details/140794609

版权

选址定容


%%%%情况1   不同位置 相同容量
Pwind=50;
Qwind=Pwind*tan(acos(0.90));
% P_bloss=[];
Pdg=zeros(33,3);
Pdg([16 17 31 32],1) = [20,20,50,70];
Pdg([16 17 31 32],2) = [40,40,100,140];
Pdg([16 17 31 32],3) = [60,60,150,210];

Vot=zeros(33,3);
P_bloss=[];VV=[];
for i=1:3
mpc = case33bw;
mpc.bus(:,3:4)=mpc.bus(:,3:4)-0.001*[Pdg(:,i),-Pdg(:,i)*tan(acos(0.90))];
OPT = mpoption('verbose', 0,  'out.all', 0);
results=runpf(mpc,OPT);
if results.success==1
P_bloss(i,1)=sum(results.branch(:,14)+results.branch(:,16));
Vot(:,i)=results.bus(:,8);
VV=[VV,min(results.bus(:,8))];
Pline(:,i)=results.branch(:,14)+results.branch(:,16);
end
end
figure
plot(P_bloss,'-o')
title('3不同位置相同容量');
xlabel('支路编号')
ylabel('网损LOSS(MW)')
figure
plot(P_bloss,'-o')
title('4不同位置相同容量');
xlabel('支路编号')
ylabel('支路LOSS(MW)')

%%%情况3  相同位置，不同容量
Pwind=50;
Qwind=Pwind*tan(acos(0.90));
Pdg=zeros(33,11);
Vot=zeros(33,11);
P_bloss=[];VV=[];Pline=[];
for i=1:11
Pdg(17,i)=[i*200-200];
Pdg(32,i)=[i*300-300];
mpc = case33bw;
mpc.bus(:,3:4)=mpc.bus(:,3:4)-0.001*[Pdg(:,i),-Pdg(:,i)*tan(acos(-0.90))];
OPT = mpoption('verbose', 0,  'out.all', 0);
results=runpf(mpc,OPT);
if results.success==1
P_bloss(:,i)=sum(results.branch(:,14)+results.branch(:,16));
Pline(:,i)=results.branch(:,14)+results.branch(:,16);
Vot(:,i)=results.bus(:,8);
VV=[VV,min(results.bus(:,8))];

end
end
figure
plot(1:i,P_bloss,'-o')
title('节点17');
xlabel('接入容量*200')
ylabel('网损LOSS(kW)')
figure
plot(Pline,'-o')
title('线损');
xlabel('支路编号')
ylabel('线损LOSS(kW)')
legend('接入50KW','接入100KW')

粒子群

适应度函数

function y = Fac(x)
global flag;
% x 代表DG容量
  x1 = x(1,1);%%光伏        
  x2 = x(1,2);%%风电
  %%对应容量接入的位置
  location1=fix(x(3));
  location2=fix(x(4));
  
  Cinv=0.394+0.271;
  %%年化收益率
  r=0.2;
  m=20;%%寿命
  Cr=1/r-1/(r*(1+r)^m);
  %%01投资成本
 Cost1=(x1+x2)*Cinv/Cr;
 Pv_1=[0
0
0
0
0
0.6
2.4
10.5
30
69
69.9
95.4
129.9
111
120.9
99
71.4
39.9
12.9
0.9
0
0
0
0]/200;
 kk=[0.768517191206980,0.881154126188142,0.895711895637549,0.908639736589063,0.757725433040895,0.757919510382514,0.823958189376985,0.760534450529473,0.819170680011629,0.780249729625146,0.894024450939335,0.886520299947759,0.802227298535102,0.921023535768218,0.801396777905281,0.908975427583852,0.878915510970877,0.947742849873102,0.910218455419634,0.995899131838041,0.872054190496853,0.984074396278727,0.820703260897279,0.803732355323619];
Pw_1=[2.67
2.67
2.34
3.12
3.29
4.76
4.77
4.24
3.81
4.59
3.9
4.94
3.55
4.33
3.21
3.29
3.03
3.64
3.73
2.6
3.38
3.12
3.46
3.64]/5;
P_bloss=0;VV=0;
for i=1:24
mpc = case33bw;
mpc.bus(:,3:4)=mpc.bus(:,3:4)*kk(i);%%模拟负荷变化
mpc.bus(location1,3:4)=mpc.bus(location1,3:4)-0.001*[Pv_1(i)*x1,Pv_1(i)*x1*tan(acos(0.90))];%%%光伏接入
mpc.bus(location2,3:4)=mpc.bus(location2,3:4)-0.001*[Pw_1(i)*x2,Pw_1(i)*x2*tan(acos(0.95))];%%%光伏接入

OPT = mpoption('verbose', 0,  'out.all', 0);
results=runpf(mpc,OPT);
if results.success==1
P_bloss=P_bloss+sum(results.branch(:,14)+results.branch(:,16));
Vot=results.bus(:,8);
Vmax=1.05;Vmin=0.95;
VV=VV+sum(abs(Vot-1)/(Vmax-Vmin));
else
   P_bloss=50000; VV=5000;%%%如果容量配置不合理，则进行惩罚，加快迭代。
end
end

 %%%02网损
Cost2= P_bloss*0.45;
 %%03电压偏差
Cost3=VV*1.5;

%%%最终的目标函数
y=10000/(0.25*Cost1+0.5*Cost2*0+0.25*Cost3);%%%粒子群一般求最小值，取倒数

if flag==1
 disp(['01投资成本:',num2str(Cost1),'万元'])
   disp(['02网损:',num2str(Cost2),'kW']) 
   disp(['03电压偏差:',num2str(VV),'p.u.'])
   
%    V0=[];V1=[]; PL0=[];PL1=[];
%    for i=1:24
% mpc = case33bw;
% mpc.bus(:,3:4)=mpc.bus(:,3:4)*kk(i);%%模拟负荷变化
% OPT = mpoption('verbose', 0,  'out.all', 0);
% results0=runpf(mpc,OPT);
% V0(1,i)=min(results0.bus(:,8));
% PL0(1,i)=sum(results0.branch(:,14)+results0.branch(:,16));
% 
% mpc.bus(location1,3:4)=mpc.bus(location1,3:4)-0.001*[Pv_1(i)*x1,Pv_1(i)*x1*tan(acos(0.90))];%%%光伏接入
% mpc.bus(location2,3:4)=mpc.bus(location2,3:4)-0.001*[Pw_1(i)*x2,Pw_1(i)*x2*tan(acos(0.95))];%%%光伏接入
% 
% OPT = mpoption('verbose', 0,  'out.all', 0);
% results=runpf(mpc,OPT);
% if results.success==1
% P_bloss=P_bloss+sum(results.branch(:,14)+results.branch(:,16));
% PL1(1,i)=sum(results.branch(:,14)+results.branch(:,16));
% Vot=results.bus(:,8);
% V1(1,i)=min(Vot);
% 
% Vmax=1.05;Vmin=0.95;
% VV=VV+sum(abs(Vot-1)/(Vmax-Vmin));
% else
%    P_bloss=50000; VV=5000;%%%如果容量配置不合理，则进行惩罚，加快迭代。
% end
%    end
% figure
% plot(1:24,[V1;V1],'-d')
%  legend('原始','含DG');
%  title('电压对比')
%    figure
% plot(1:24,[PL0;PL1],'-d')
%  legend('原始','含DG');
%   title('网损对比')
%  
% end
end

主函数


d=4;
n=5;
smin=[0 0 2 2];
smax=[5000 5000 33 33];
wmin=0.2;
wmax=1.2;
w=1;
wdamp=0.99;
alpha=2;
beta=2;
max_iter=13;
global flag;%%%画图标识
flag=0;
% 函数表达式
cost_func=@(x) Fac(x);
% 
% figure
% [x0_1, x0_2,y0]=DrawGriewank;
% hold on;

part=struct();
part.pos=[];
part.cost=Inf;
part.pbest=-Inf;
part.pos_pbest=[];
part.vel=[];

swarm=repmat(part,n,1);
TrACE=[];

gbest=-Inf;
pos_gbest=[];

for i=1:n
%   swarm(i).pos=unifrnd(smin,smax,1,d);
swarm(i).pos=smin+rand(1,d).*(-smin+smax);
  swarm(i).cost=cost_func(swarm(i).pos);
  swarm(i).vel=zeros(1,d);
  swarm(i).pos_pbest=swarm(i).pos;
  swarm(i).pbest=cost_func(swarm(i).pos_pbest);
  if swarm(i).pbest>gbest
    gbest=swarm(i).pbest;
    pos_gbest=swarm(i).pos_pbest;
  end
end
% TrACE=[];
for t=1:max_iter
  
  for i=1:n
  
    [~, SortOrder]=sort([swarm.cost]);
    swarm=swarm(SortOrder);
    swarm=flipud(swarm);
  %calculation of new inertia weight factor  
    avg=sum([swarm.cost])/n;
%     if swarm(i).cost>avg
%       w=wmin+((wmax-wmin)*(swarm(i).cost-swarm(1).cost)/(avg-swarm(1).cost));
%     else
%       w=wmax;
%     end 
%    disp([num2str(swarm(i).pbest)]);
    swarm(i).vel=w*swarm(i).vel+alpha*rand(1,d).*(pos_gbest-swarm(i).pos)+beta*rand(1,d).*(swarm(i).pos_pbest-swarm(i).pos);  
    swarm(i).pos=swarm(i).pos+swarm(i).vel;
    swarm(i).pos=max(swarm(i).pos,smin);
    swarm(i).pos=min(swarm(i).pos,smax);
    swarm(i).cost=cost_func(swarm(i).pos);
%    disp([num2str(swarm(i).cost < swarm(i).pbest)]);
    if (swarm(i).cost > swarm(i).pbest)
      swarm(i).pos_pbest=swarm(i).pos;
      swarm(i).pbest=swarm(i).cost;
    end
    if (swarm(i).pbest > gbest)
      pos_gbest=swarm(i).pos_pbest;
      gbest=swarm(i).pbest;
    end
  
  end
  w=w*wdamp;
  disp(['Iteration : ' num2str(t) ' Minimum value = ' num2str(gbest) ' Best solution is : ' num2str(pos_gbest)]);
  TrACE(1,t)=gbest;
end
   %% 这段用来画出最后的结果
figure();
%在搜索最优个体的时候，最佳的位置存储在中
  plot(TrACE,'->')
  title('迭代过程');
%%%输出结果
pos_gbest(3:4)=fix(pos_gbest(3:4));
disp(['光伏最佳接入位置:',num2str(pos_gbest(3)),'节点','最佳接入容量',num2str((pos_gbest(1))),'kW'])
disp(['风电最佳接入位置:',num2str(pos_gbest(4)),'节点','最佳接入容量',num2str(fix(pos_gbest(2))),'kW'])
flag=1;
cost_func(pos_gbest);