Optimal-economy.rar_economic_matlab 最优化_optimization_power optim

clear; tic; % 算例为 IEEE33 %NodeData(i,3)=NodeData(i,3)-Q(i,2); SS=[1 0 1 0.0922 0.047 100 60 2 0 22 2 1 2 0.493 0.2511 90 40 2 0 15 3 2 3 0.366 0.1864 120 80 1 0 20 4 3 4 0.3811 0.1941 60 30 1 0 15 5 4 5 0.8190 0.707 60 20 1 0 36 6 5 6 0.1872 0.6188 200 100 7 0 38 7 6 7 0.7114 0.2351 200 100 1 0 27 8 7 8 1.03 0.74 60 20 8 0 15 9 8 9 1.044 0.74 60 20 2 0 22 10 9 10 0.1966 0.065 45 30 1 0 14 11 10 11 0.3744 0.1238 60 35 4 0 22 12 11 12 1.468 1.155 60 35 3 0 20 13 12 13 0.5416 0.7129 120 80 9 0 23 14 13 14 0.5910 0.526 60 10 3 0 20 15 14 15 0.7463 0.5450 60 20 2 0 36 16 15 16 1.289 1.7210 60 20 14 0 35 17 16 17 0.7320 0.574 90 40 2 0 29 18 1 18 0.164 0.1565 90 40 6 0 28 19 18 19 1.5042 1.3554 90 40 2 0 33 20 19 20 0.4095 0.4784 90 40 4 0 27 21 20 21 0.7089 0.9373 90 40 6 0 25 22 2 22 0.4512 0.3083 90 50 5 0 22 23 22 23 0.8980 0.7091 420 200 5 0 15 24 23 24 0.8960 0.7011 420 200 3 0 25 25 5 25 0.2030 0.1034 60 25 7 0 45 26 25 26 0.2842 0.1447 60 25 2 0 62 27 26 27 1.059 0.9337 60 20 7 0 47 28 27 28 0.8042 0.7006 120 70 1 0 51 29 28 29 0.5075 0.2585 200 600 4 0 47 30 30 29 0.9744 0.9630 150 70 3 0 44 31 30 31 0.3105 0.3619 210 100 4 0 74 32 31 32 0.3410 0.5302 60 40 2 0 24 33 20 7 2.0 2.0 90 40 4 0 27 34 8 14 2.0 20.0 60 10 5 0 20 35 21 11 2.0 2.0 90 40 4 0 25 36 17 32 0.5 0.5 60 40 3 0 24 37 24 28 0.5 0.5 120 70 4 0 51]; a=0.1;%主干线故障率； b=3;%主干线平均修复时间； c=0.25;%分支线的故障率； d=1;%分支线的平均修复时间； MM=1; DIE=1; T=0; k=1; hua=1; huatu(1,1)=202.7 hua=hua+1; %可断开关集 D1=[33 6 19 18 7 20]; D2=[9 10 11 21 35]; D3=[12 13 14 34]; D4=[22 23 24 37]; D5=[15 16 17 25 26 27 29 30 31 32 36]; [x1,ch1]=size(D1); [x2,ch2]=size(D2); [x3,ch3]=size(D3); [x4,ch4]=size(D4); [x5,ch5]=size(D5); for xunhuan=1:1000 if T==1 break; end sj1=randint(1,1,[1 ch1]); d1=D1(1,sj1); sj2=randint(1,1,[1 ch2]); d2=D2(1,sj2); sj3=randint(1,1,[1 ch3]); d3=D3(1,sj3); sj4=randint(1,1,[1 ch4]); d4=D4(1,sj4); sj5=randint(1,1,[1 ch5]); d5=D5(1,sj5); FF=SS; FF(FF(:,1)==d1,:)=[]; FF(FF(:,1)==d2,:)=[]; FF(FF(:,1)==d3,:)=[]; FF(FF(:,1)==d4,:)=[]; FF(FF(:,1)==d5,:)=[]; QZ=FF; %%%找出T节点并存储到TT矩阵中 UB = 12.66; % 电压基准 kV SB = 10; % 功率基准 MVA ZB = UB^2/SB; % 阻抗基准 ohm QZ(:,[4,5]) = QZ(:,[4,5]) / ZB; % 阻抗标幺化 QZ(:,[6,7]) = QZ(:,[6,7]) / SB / 1000;% 功率标幺化 S = [QZ(:,6) + j*QZ(:,7)];% 形成 S %disp(S) ZL = [QZ(:,4) + j*QZ(:,5)]; % 形成 ZL %disp('ZL'); %disp(ZL); b=32;%节点数 k=0; V=ones(b+1,1); t=0; %%%找出T节点并存储到TT矩阵中 bb=1;%判断是否为T节点 TT=zeros(1,10);%存储T节点 for r=1:b t=0; for p=1:b if (QZ(p,2)==r)|(QZ(p,3)==r) t=t+1; end end if t>2 TT(bb,1)=r; bb=bb+1; t=0; end end for i=1:bb-1 x=2; for jj=1:b if QZ(jj,2)==TT(i,1) TT(i,x)=QZ(jj,3); x=x+1; end end for jj=1:b if QZ(jj,3)==TT(i,1) TT(i,x)=QZ(jj,2); x=x+1; end end end %disp('T节点'); %disp(TT); bbb=1;%判断是否为末梢节点 MS=zeros(1,10);%存储末梢节点 for r=1:b t=0; for p=1:b if (QZ(p,2)==r)|(QZ(p,3)==r) t=t+1; end end if t==1 MS(1,bbb)=r; bbb=bbb+1; t=0; end end %disp('末梢节点'); %disp(MS); %disp(bbb); A=0;%首节点 B=0;%末节点 C=0; W=0;%交换 %调整网络结构矩阵 m=b; G=zeros(1,m);%节点是否包含在末节点中 G(1,1)=1; H=0;%基点的上一个节点 A=1;%基点 %进行纵向搜索 for jj=1:m H=0;%基点的上一个节点 A=1;%基点 for g=2:m if (QZ(g,2)==A)&(QZ(g,3)~=H)&(G(1,QZ(g,3))==0) H=A; A=QZ(g,3); G(1,QZ(g,3))=1; end for g=2:m if (QZ(g,3)==A)&(G(1,QZ(g,2))==0) A=QZ(g,2);H=QZ(g,3); W=QZ(g,2);QZ(g,2)=QZ(g,3);QZ(g,3)=W; G(1,QZ(g,3))=1; end end end end %disp('********************G'); %disp(QZ); %搜索T节点进行横向搜索 for jg=1:2 for jt=1:bb-1 A=TT(jt,1); for gt=1:m for it=2:m if (QZ(it,2)==A)&(G(1,QZ(it,3))==0) A=QZ(it,3);H=QZ(it,2);G(1,QZ(it,3))=1; end end for it=2:m if (QZ(it,3)==A)&(G(1,QZ(it,2))==0) A=QZ(it,2);H=QZ(it,3); W=QZ(it,2);QZ(it,2)=QZ(it,3);QZ(it,3)=W; G(1,QZ(it,3))=1; end end end end end %%%%%%%%%%%%%%%%%% for jj=1:bbb-1 A=MS(1,jj); for i=1:m if (QZ(i,2)==MS(1,jj))&(G(1,QZ(i,3))==0) W=QZ(i,2); QZ(i,2)=QZ(i,3); QZ(i,3)=W; A=QZ(i,2); B=QZ(i,3); G(1,QZ(B,3))=1; end end for i=1:m if (QZ(i,2)==A&QZ(i,3)~=B) W=QZ(i,2); QZ(i,2)=QZ(i,3); QZ(i,3)=W; A=QZ(i,2); B=QZ(i,3); end end end %%矩阵排序 Huan=zeros(1,10); for i=1:m-1 for jj=i:m if QZ(i,3)>QZ(jj,3) Huan=QZ(i,:); QZ(i,:)=QZ(jj,:); QZ(jj,:)=Huan; end end end %迭代开始处 Delta = 1; % 收敛判据赋初值 TempV = V; % 赋初值，用于记忆上次迭代结果 kk=0; while (Delta > 1e-8)&(kk<40) %算节点注入电流 for l=1:b jj=QZ(l,3); ua=V(jj+1,1); I(jj,1)=conj(S(jj,1)/ua); end TJ=zeros(b,1);%标出T节点 for i=1:bb-1 TJ(TT(i,1),1)=1; end %回推算支路电流 J=zeros(b,1); for y=1:bbb-1 A=MS(1,y); B=b; for yy=1:b for yyy=1:b if (QZ(yyy,3)==A)&(B~=QZ(yyy,2))&(TJ(QZ(yyy,3),1)~=1) A=QZ(yyy,2); J(yyy,1)=J(B,1)+I(yyy,1); B=QZ(yyy,3); end end end end for ju=1:bb-1 for i=1:bb-1 if (J(TT(i,2),1)~=0)&(J(TT(i,3),1)~=0)&(J(TT(i,1),1)==0) J(TT(i,1),1)=J(TT(i,2),1)+J(TT(i,3),1)+I(TT(i,1),1); A=QZ(TT(i,1),2); B=TT(i,1); for iii=1:b for ii=1:b if (QZ(ii,3)==A)&(B~=QZ(ii,2))&(TJ(QZ(ii,3),1)~=1) A=QZ(ii,2); J(ii,1)=J(B,1)+I(ii,1); B=QZ(ii,3); end end end end end end %disp(J); %前推算节点电压 for l=1:b jj=QZ(l,3)+1; i=QZ(l,2)+1; V(jj,1)=V(i,1)-ZL(l,1)*J(l,1); end Delta = max(abs(V-TempV)); % 更新收敛判据 TempV = V; % 记忆迭代结果 kk=kk+1; end %求支路电流最小值 minI=1;%支路最小电流 for i=1:m if minI>real(J(i,1)) minI=real(J(i,1)); end end minV=1;%节点最小电压 VV=abs(V); for i=1:m if minV>VV(i,1) minV=VV(i,1); end end tiaojian=0; if (minI>0)&(minV>0.9) tiaojian=1; end %计算支路功率损耗 plo=2; if tiaojian==1 for l=1:b jj=QZ(l,3)+1; i=QZ(l,2)+1; dV(l,1)=V(i,1)-V(jj,1); end %ploss=zeros(1,32);%支路功率损耗 %pl=(dV.*conj(J)).*10; ploss=sum((dV.*conj(J)).*10); plo=real(ploss); end if (MM>plo)&(tiaojian==1) MM=plo; duan1=d1; duan2=d2; duan3=d3; duan4=d4; duan5=d5; DIE=1; huatu(hua,1)=MM*1000; hua=hua+1; elseif tiaojian==1 DIE=DIE-1; end T=0; disp('*****************程序运行中********************') if DIE==0 disp('输出网络初始损耗'); disp(202.7); disp('输出断开的支路编号'); disp(duan1); disp(duan2); disp(duan3); disp(duan4); disp(duan5); disp('输出节点电压'); d