parent
115c8b540b
commit
6285adc0d0
23
AssignXX.m
23
AssignXX.m
|
|
@ -1,23 +0,0 @@
|
||||||
function [deltZ,deltL,deltW,deltU,deltX,deltY]=AssignXX1(XX,ContrlCount,RestraintCount,Busnum)
|
|
||||||
% deltX=XX(1:14);
|
|
||||||
% deltY=XX(15:24);
|
|
||||||
% deltZ=XX(25:38);
|
|
||||||
% deltW=XX(39:52);
|
|
||||||
% deltL=XX(53:66);
|
|
||||||
% deltU=XX(67:80);
|
|
||||||
deltX=XX(1:ContrlCount);
|
|
||||||
k1=ContrlCount+2*Busnum;
|
|
||||||
deltY=XX(ContrlCount+1:k1);
|
|
||||||
k2=k1;
|
|
||||||
k1=k2+RestraintCount;
|
|
||||||
deltZ=XX(k2+1:k1);
|
|
||||||
k2=k1;
|
|
||||||
k1=k2+RestraintCount;
|
|
||||||
deltW=XX(k2+1:k1);
|
|
||||||
k2=k1;
|
|
||||||
k1=k2+RestraintCount;
|
|
||||||
deltL=XX(k2+1:k1);
|
|
||||||
k2=k1;
|
|
||||||
k1=k2+RestraintCount;
|
|
||||||
deltU=XX(k2+1:k1);
|
|
||||||
end
|
|
||||||
13
CalCost.m
13
CalCost.m
|
|
@ -1,13 +0,0 @@
|
||||||
function CalCost(GenC,PG,PGi)
|
|
||||||
cost=GenC(:,2).*PG(PGi).^2+GenC(:,3).*PG(PGi)+GenC(:,4);
|
|
||||||
% Org_PG=[5;
|
|
||||||
% 2.5794];
|
|
||||||
% book_PG=[5.5056;
|
|
||||||
% 2.1568];
|
|
||||||
% cost2=GenC(:,2).*Org_PG(1:2).^2+GenC(:,3).*Org_PG(1:2)+GenC(:,4);
|
|
||||||
% cost3=GenC(:,2).*book_PG(1:2).^2+GenC(:,3).*book_PG(1:2)+GenC(:,4);
|
|
||||||
fprintf('总花费为%f\n',sum(cost,1));
|
|
||||||
% fprintf('PF总花费为%f\n',sum(cost2,1));
|
|
||||||
% fprintf('书上OPF总花费为%f\n',sum(cost3,1));
|
|
||||||
% fprintf('较书上减少费用为为%f\n',sum(cost3,1)-sum(cost,1));
|
|
||||||
end
|
|
||||||
|
|
@ -1,6 +0,0 @@
|
||||||
function DrawGap(plotGap)
|
|
||||||
x=find(plotGap);
|
|
||||||
ts=size(x,2);
|
|
||||||
|
|
||||||
|
|
||||||
end
|
|
||||||
|
|
@ -1,6 +0,0 @@
|
||||||
function DrawGap(plotGap)
|
|
||||||
x=find(plotGap);
|
|
||||||
ts=size(x,2);
|
|
||||||
plot(1:ts,plotGap(1:ts));
|
|
||||||
|
|
||||||
end
|
|
||||||
12
FormAA.asv
12
FormAA.asv
|
|
@ -1,12 +0,0 @@
|
||||||
function AA=FormAA(L_1Z,deltG,U_1W,Hcoma,deltH)
|
|
||||||
tOnes=eye(14);
|
|
||||||
tZeros=zeros(14);
|
|
||||||
AA=[
|
|
||||||
tOnes,L_1Z,tZeros,tZeros,tZeros,zeros(14,10);
|
|
||||||
tZeros,tOnes,tZeros,tZeros,-deltG',zeros(14,10);
|
|
||||||
tZeros,tZeros,tOnes,U_1W,tZeros,zeros(14,10);
|
|
||||||
tZeros,tZeros,tZeros,tOnes,deltG,zeros(14,10);
|
|
||||||
tZeros,tZeros,tZeros,tZeros,Hcoma,deltH;
|
|
||||||
zeros(10,14),zeros(10,14),zeros(10,14),zeros(10,14),deltH',zeros(10,10);
|
|
||||||
];
|
|
||||||
end
|
|
||||||
|
|
@ -1,9 +0,0 @@
|
||||||
function AlphaD=FormAlphaD(Init_Z,deltZ,Init_W,deltW)
|
|
||||||
tdeltZinx=find(deltZ<0);
|
|
||||||
tdeltWinx=find(deltW>0);
|
|
||||||
t1=-Init_Z(tdeltZinx)./deltZ(tdeltZinx)';
|
|
||||||
t2=-Init_W(tdeltWinx)./deltW(tdeltWinx)';
|
|
||||||
t3=[t1,t2];
|
|
||||||
t4=min(t3);
|
|
||||||
AlphaD=0.9995*min([t4 1]);
|
|
||||||
end
|
|
||||||
|
|
@ -1,7 +0,0 @@
|
||||||
function AlphaP=FormAlphaP(Init_L,deltL,Init_U,deltU)
|
|
||||||
ti=deltL(delt)
|
|
||||||
t1=-Init_L./deltL';
|
|
||||||
t2=-Init_U./deltU';
|
|
||||||
t3=[t1,t2];
|
|
||||||
t4=t3()
|
|
||||||
end
|
|
||||||
|
|
@ -1,9 +0,0 @@
|
||||||
function AlphaP=FormAlphaP(Init_L,deltL,Init_U,deltU)
|
|
||||||
tdeltLinx=find(deltL<0);
|
|
||||||
tdeltUinx=find(deltU<0);
|
|
||||||
t1=-Init_L(tdeltLinx)./deltL(tdeltLinx)';
|
|
||||||
t2=-Init_U(tdeltUinx)./deltU(tdeltUinx)';
|
|
||||||
t3=[t1,t2];
|
|
||||||
t4=min(t3);
|
|
||||||
AlphaP=0.9995*min([t4 1]);
|
|
||||||
end
|
|
||||||
26
FormG.asv
26
FormG.asv
|
|
@ -1,26 +0,0 @@
|
||||||
function Mat_G=FormG(Volt,PVi,PGi,PG,QG,PD,QD,Loadi)
|
|
||||||
PD_=PD(Loadi);
|
|
||||||
QD_=QD(Loadi);
|
|
||||||
Gsize=length(PGi)+length(PVi)+2*length(Loadi);
|
|
||||||
Mat_G=sdpvar(Gsize,1);
|
|
||||||
add=1;
|
|
||||||
for I=1:length(PGi)
|
|
||||||
Mat_G(add)=PG(I).value
|
|
||||||
add=add+1;
|
|
||||||
end
|
|
||||||
for I=1:length(PVi)
|
|
||||||
Mat_G(add)=QG(I).value
|
|
||||||
add=add+1;
|
|
||||||
end
|
|
||||||
for I=1:length(Loadi)
|
|
||||||
Mat_G(add)=QG(I).value
|
|
||||||
add=add+1;
|
|
||||||
end
|
|
||||||
% Mat_G=[
|
|
||||||
% [PG(PGi).value]';
|
|
||||||
% [QG(PVi).value]';
|
|
||||||
% [PD(Loadi).value]';
|
|
||||||
% [QD(Loadi).value]';
|
|
||||||
% Volt;
|
|
||||||
% ];
|
|
||||||
end
|
|
||||||
34
FormG.m
34
FormG.m
|
|
@ -1,34 +0,0 @@
|
||||||
function Mat_G=FormG(Volt,PVi,PGi,PG,QG,PD,QD,Loadi)
|
|
||||||
PD_=PD(Loadi);
|
|
||||||
QD_=QD(Loadi);
|
|
||||||
Gsize=length(PGi)+length(PVi)+2*length(Loadi)+length(Volt);
|
|
||||||
Mat_G=sdpvar(Gsize,1);
|
|
||||||
add=1;
|
|
||||||
for I=1:length(PGi)
|
|
||||||
Mat_G(add)=PG(I).value;
|
|
||||||
add=add+1;
|
|
||||||
end
|
|
||||||
for I=1:length(PVi)
|
|
||||||
Mat_G(add)=QG(I).value;
|
|
||||||
add=add+1;
|
|
||||||
end
|
|
||||||
for I=1:length(Loadi)
|
|
||||||
Mat_G(add)=PD(I).value;
|
|
||||||
add=add+1;
|
|
||||||
end
|
|
||||||
for I=1:length(Loadi)
|
|
||||||
Mat_G(add)=QD(I).value;
|
|
||||||
add=add+1;
|
|
||||||
end
|
|
||||||
for I=1:length(Volt)
|
|
||||||
Mat_G(add)=Volt(I);
|
|
||||||
add=add+1;
|
|
||||||
end
|
|
||||||
% Mat_G=[
|
|
||||||
% [PG(PGi).value]';
|
|
||||||
% [QG(PVi).value]';
|
|
||||||
% [PD(Loadi).value]';
|
|
||||||
% [QD(Loadi).value]';
|
|
||||||
% Volt;
|
|
||||||
% ];
|
|
||||||
end
|
|
||||||
31
FormH.asv
31
FormH.asv
|
|
@ -1,31 +0,0 @@
|
||||||
function Mat_H=FormH(Busnum,Volt,PG,PD,QG,QD,Y,UAngel,r,c,Angle,Loadi,notLoadi,PGi,PVi,PD0,QD0)
|
|
||||||
%%
|
|
||||||
%QDcos=textread('D:\Project\×îС»¯³±Á÷\×îС³±Á÷ËãÀý\Ïɺ£919PDQDglys.txt');
|
|
||||||
%QD(QD~=0)=PD(QD~=0)./tan(QDcos);
|
|
||||||
%QD(QD_NON_ZERO_IND)=QD_NON_ZERO;
|
|
||||||
%%
|
|
||||||
%PD(Loadi)=QD(Loadi)./tan(acos(0.98));
|
|
||||||
AngleIJ=sparse(r,c,UAngel(r)-UAngel(c)-Angle,Busnum,Busnum);
|
|
||||||
Loadi=find(QD0~=0 | PD0~=0);
|
|
||||||
PDi=find(PD0~=0);
|
|
||||||
QDi=find(QD0~=0);
|
|
||||||
notLoadi=setdiff(1:Busnum,Loadi);
|
|
||||||
PGOnly=setdiff(PGi,PDi);
|
|
||||||
QGOnly=setdiff(PVi,QDi);
|
|
||||||
PDOnly=setdiff(PDi,PGi);
|
|
||||||
QDOnly=setdiff(QDi,PVi);
|
|
||||||
PGPD=intersect(PDi,PGi);
|
|
||||||
QGQD=intersect(QDi,PVi);
|
|
||||||
noPGPDQGQD=setdiff(1:Busnum,union(union(PDi,PGi),union(QDi,PVi)));
|
|
||||||
YCos=Y.*cos(AngleIJ);
|
|
||||||
YSin=Y.*sin(AngleIJ);
|
|
||||||
dP1=PG(ismember(PGi,PGOnly))-diag(Volt(PGOnly))*YCos(PGOnly,:)*Volt(PGOnly);
|
|
||||||
dQ1=QG(ismember(PVi,QGOnly))-diag(Volt(QGOnly))*YSin(QGOnly,:)*Volt(QGOnly);
|
|
||||||
dP2=PD(ismember(PDi,PDOnly))-diag(Volt(PDOnly))*YCos(PDOnly,:)*Volt(PDOnly);
|
|
||||||
dQ2=QD(ismember(QDi,QDOnly))-diag(Volt(QOnly))*YSin(QDOnly,:)*Volt(QDOnly);
|
|
||||||
dP3=PG(ismember(PGi,PGOnly))-diag(Volt(PGOnly))*YCos(PGOnly,:)*Volt(PGOnly);
|
|
||||||
dQ=QG-QD-diag(Volt)*Y.*sin(AngleIJ)*Volt;
|
|
||||||
|
|
||||||
Mat_H=[dP;dQ;];
|
|
||||||
|
|
||||||
end
|
|
||||||
40
FormH.m
40
FormH.m
|
|
@ -1,40 +0,0 @@
|
||||||
function Mat_H=FormH(Busnum,Volt,PG,PD,QG,QD,Y,UAngel,r,c,Angle,Loadi,notLoadi,PGi,PVi,PD0,QD0,dP,dQ)
|
|
||||||
%%
|
|
||||||
%QDcos=textread('D:\Project\×îС»¯³±Á÷\×îС³±Á÷ËãÀý\Ïɺ£919PDQDglys.txt');
|
|
||||||
%QD(QD~=0)=PD(QD~=0)./tan(QDcos);
|
|
||||||
%QD(QD_NON_ZERO_IND)=QD_NON_ZERO;
|
|
||||||
%%
|
|
||||||
%PD(Loadi)=QD(Loadi)./tan(acos(0.98));
|
|
||||||
AngleIJ=sparse(r,c,UAngel(r)-UAngel(c)-Angle,Busnum,Busnum);
|
|
||||||
Loadi=find(QD0~=0 | PD0~=0);
|
|
||||||
PDi=find(PD0~=0);
|
|
||||||
QDi=find(QD0~=0);
|
|
||||||
notLoadi=setdiff(1:Busnum,Loadi);
|
|
||||||
PGOnly=setdiff(PGi,PDi);
|
|
||||||
QGOnly=setdiff(PVi,QDi);
|
|
||||||
PDOnly=setdiff(PDi,PGi);
|
|
||||||
QDOnly=setdiff(QDi,PVi);
|
|
||||||
PGPD=intersect(PDi,PGi);
|
|
||||||
QGQD=intersect(QDi,PVi);
|
|
||||||
noPGPDQGQD=setdiff(1:Busnum,union(union(PDi,PGi),union(QDi,PVi)));
|
|
||||||
%%
|
|
||||||
YCos=Y.*cos(AngleIJ);
|
|
||||||
YSin=Y.*sin(AngleIJ);
|
|
||||||
% dP1=PG(ismember(PGi,PGOnly))-diag(Volt(PGOnly))*YCos(PGOnly,:)*Volt(PGOnly);
|
|
||||||
% dQ1=QG(ismember(PVi,QGOnly))-diag(Volt(QGOnly))*YSin(QGOnly,:)*Volt(QGOnly);
|
|
||||||
% dP2=PD(ismember(PDi,PDOnly))-diag(Volt(PDOnly))*YCos(PDOnly,:)*Volt(PDOnly);
|
|
||||||
% dQ2=QD(ismember(QDi,QDOnly))-diag(Volt(QOnly))*YSin(QDOnly,:)*Volt(QDOnly);
|
|
||||||
% dP3=PG(ismember(PGi,PGOnly))-diag(Volt(PGOnly))*YCos(PGOnly,:)*Volt(PGOnly);
|
|
||||||
% dQ=QG-QD-diag(Volt)*Y.*sin(AngleIJ)*Volt;
|
|
||||||
%%
|
|
||||||
|
|
||||||
for I=1:Busnum
|
|
||||||
dP(I)=PG(I).value-PD(I).value-Volt(I)*YCos(I,:)*Volt;
|
|
||||||
dQ(I)=QG(I).value-QD(I).value-Volt(I)*YSin(I,:)*Volt;
|
|
||||||
end
|
|
||||||
%%
|
|
||||||
% dP=[PG(:).value]'-[PD(:).value]'-diag(Volt)*YCos*Volt;
|
|
||||||
% dQ=[QG(:).value]'-[QD(:).value]'-diag(Volt)*YSin*Volt;
|
|
||||||
Mat_H=[dP;dQ;];
|
|
||||||
|
|
||||||
end
|
|
||||||
16
FormLw.asv
16
FormLw.asv
|
|
@ -1,16 +0,0 @@
|
||||||
function Lw=FormLw(Mat_G,Init_U,GenU,Busnum,PVQU,PD0,Loadi,KK)
|
|
||||||
%KK=999;
|
|
||||||
PU=1*GenU(:,2);%发电机有功上界
|
|
||||||
QU=PVQU(:,1);%发电机无功上界
|
|
||||||
%VoltU=(1.5+1/exp(KK))*ones(1,Busnum);
|
|
||||||
VoltU=1.5*ones(1,Busnum);
|
|
||||||
PDU=PD0(Loadi);
|
|
||||||
% PDU(PDU>0)=1200*PDU(PDU>0);
|
|
||||||
% PDU(PDU<0)=-800*PDU(PDU<0);
|
|
||||||
% PDU(PDU==0)=200;
|
|
||||||
PDU=1000000*ones(length(Loadi),1);
|
|
||||||
t1=([PU',QU',PDU',VoltU])';
|
|
||||||
t2=Mat_G+Init_U'-t1;
|
|
||||||
Lw=t2;
|
|
||||||
|
|
||||||
end
|
|
||||||
25
FormLw.m
25
FormLw.m
|
|
@ -1,25 +0,0 @@
|
||||||
function Lw=FormLw(Mat_G,GenU,Busnum,PVQU,PD0,QD0,Loadi,PF)
|
|
||||||
%PU=GenU(:,2);%发电机有功上界
|
|
||||||
PU=8*ones(length(GenU(:,2)),1);
|
|
||||||
%QU=PVQU(:,1);%发电机无功上界
|
|
||||||
QU=8*ones(length(PVQU(:,1)),1);
|
|
||||||
VoltU=(1.2)*ones(1,Busnum);
|
|
||||||
%VoltU=10*ones(1,Busnum);
|
|
||||||
PDU=PD0(Loadi);
|
|
||||||
PDU(PDU>0)=1.200*PDU(PDU>0);
|
|
||||||
PDU(PDU<0)=0.800*PDU(PDU<0);
|
|
||||||
PDU(PDU==0)=0.400;
|
|
||||||
%PDU=10*ones(length(Loadi),1);
|
|
||||||
QDU=QD0(Loadi);
|
|
||||||
QDU(QDU>0)=1.200*QDU(QDU>0);
|
|
||||||
QDU(QDU<0)=0.800*QDU(QDU<0);
|
|
||||||
QDU(QDU==0)=0.200;
|
|
||||||
%PF_=1.1*PF;
|
|
||||||
%PF_(PF_>1)=1;
|
|
||||||
%powerFactorU=2*log(1).*ones(length(Loadi),1);
|
|
||||||
%t1=([PU',QU',PDU',QDU',VoltU,powerFactorU'])';
|
|
||||||
t1=([PU',QU',PDU',QDU',VoltU])';
|
|
||||||
t2=Mat_G-t1;
|
|
||||||
Lw=t2;
|
|
||||||
|
|
||||||
end
|
|
||||||
4
FormLx.m
4
FormLx.m
|
|
@ -1,4 +0,0 @@
|
||||||
function Lx=FormLx(deltF,deltH,Init_Y,deltG,Init_Z,Init_W)
|
|
||||||
t1=deltF-deltH*Init_Y'-deltG*(Init_Z'+Init_W');
|
|
||||||
Lx=t1;
|
|
||||||
end
|
|
||||||
11
FormLxComa.m
11
FormLxComa.m
|
|
@ -1,11 +0,0 @@
|
||||||
function LxComa=FormLxComa(deltF,deltG,deltH,Init_L,Luu,Lul,Init_Z,Init_Y,Lz,Init_U,Init_W,Lw,Lx)
|
|
||||||
%t1=deltG*(Init_Z'+Init_W');%%
|
|
||||||
t2=Lul+diag(Init_Z)*Lz;
|
|
||||||
t3=inv(diag(Init_L));
|
|
||||||
t4=t3*t2;%
|
|
||||||
t5=Luu-diag(Init_W)*Lw;
|
|
||||||
t6=inv(diag(Init_U));
|
|
||||||
t7=t6*t5;%
|
|
||||||
t8=deltG*(t4+t7);%%
|
|
||||||
LxComa=Lx+t8;
|
|
||||||
end
|
|
||||||
22
FormLz.asv
22
FormLz.asv
|
|
@ -1,22 +0,0 @@
|
||||||
function Lz=FormLz(Mat_G,Init_L,GenL,Busnum,PVQL,PD0,QD0,Loadi,KK)
|
|
||||||
KK=999;
|
|
||||||
PL=GenL(:,2);%发电机有功下界
|
|
||||||
QL=PVQL(:,1);%发电机无功下界
|
|
||||||
VoltL=(0.8)*ones(1,Busnum);
|
|
||||||
%VoltL=-10*ones(1,Busnum);
|
|
||||||
PDL=PD0(Loadi);
|
|
||||||
PDL(PDL>0)=0.800*PDL(PDL>0);
|
|
||||||
PDL(PDL<0)=1.200*PDL(PDL<0);
|
|
||||||
PDL(PDL==0)=-0.400;
|
|
||||||
%PDL=-10*ones(length(Loadi),1);
|
|
||||||
% QDL=QD0(Loadi);
|
|
||||||
% QDL(QDL>0)=0.800*QDL(QDL>0);
|
|
||||||
% QDL(QDL<0)=1.200*QDL(QDL<0);
|
|
||||||
% QDL(QDL==0)=-0.200;
|
|
||||||
powerFactorL=2*log(cos(0.90))*ones(length(Loadi),1);
|
|
||||||
%t1=([PL',QL',PDL',QDL',VoltL,powerFactorL'])';
|
|
||||||
t1=([PL',QL',PDL',QDL,VoltL,powerFactorL'])';
|
|
||||||
t2=Mat_G-Init_L'-t1;
|
|
||||||
Lz=t2;
|
|
||||||
|
|
||||||
end
|
|
||||||
23
FormLz.m
23
FormLz.m
|
|
@ -1,23 +0,0 @@
|
||||||
function Lz=FormLz(Mat_G,GenL,Busnum,PVQL,PD0,QD0,Loadi,PF)
|
|
||||||
%PL=GenL(:,2);%发电机有功下界
|
|
||||||
PL=-60*ones(length(GenL(:,2)),1);
|
|
||||||
%QL=PVQL(:,1);%发电机无功下界
|
|
||||||
QL=-60*ones(length(PVQL(:,1)),1);
|
|
||||||
VoltL=(0.8)*ones(1,Busnum);
|
|
||||||
%VoltL=-10*ones(1,Busnum);
|
|
||||||
PDL=PD0(Loadi);
|
|
||||||
PDL(PDL>0)=0.800*PDL(PDL>0);
|
|
||||||
PDL(PDL<0)=1.200*PDL(PDL<0);
|
|
||||||
PDL(PDL==0)=-0.400;
|
|
||||||
%PDL=-10*ones(length(Loadi),1);
|
|
||||||
QDL=QD0(Loadi);
|
|
||||||
QDL(QDL>0)=0.800*QDL(QDL>0);
|
|
||||||
QDL(QDL<0)=1.200*QDL(QDL<0);
|
|
||||||
QDL(QDL==0)=-0.200;
|
|
||||||
%powerFactorL=2*log(0.75).*ones(length(Loadi),1);
|
|
||||||
%t1=([PL',QL',PDL',QDL',VoltL,powerFactorL'])';
|
|
||||||
t1=([PL',QL',PDL',QDL',VoltL])';
|
|
||||||
t2=Mat_G-t1;
|
|
||||||
Lz=t2;
|
|
||||||
|
|
||||||
end
|
|
||||||
|
|
@ -1,5 +0,0 @@
|
||||||
function FormYY(Init_L,Lul.Lz,Init_U,Luu,Lw,Lz,LxComa)
|
|
||||||
t=[
|
|
||||||
-
|
|
||||||
]
|
|
||||||
end
|
|
||||||
10
FormYY.m
10
FormYY.m
|
|
@ -1,10 +0,0 @@
|
||||||
function YY=FormYY(Lul,Lz,Ly,Luu,Lw,Lx)
|
|
||||||
YY=[
|
|
||||||
Lx;
|
|
||||||
-Ly;
|
|
||||||
-Lz;
|
|
||||||
-Lw;
|
|
||||||
-Lul;
|
|
||||||
-Luu;
|
|
||||||
];
|
|
||||||
end
|
|
||||||
17
FormddgZW.m
17
FormddgZW.m
|
|
@ -1,17 +0,0 @@
|
||||||
function [ ddgZW ] = FormddgZW(PGi,PVi,Busnum,ContrlCount,Loadi,PD,QD,Init_Z,Init_W)
|
|
||||||
|
|
||||||
ddg1=sparse(size(PVi,1)+size(PGi,1),ContrlCount);
|
|
||||||
%¼ÆËãWºÍZÓ¦¸Ã³ËµÄλÖÃ
|
|
||||||
rePos1=length(PGi)+length(PVi)+length(Loadi)*2+Busnum+1;
|
|
||||||
rePos2=rePos1+length(Loadi)-1;
|
|
||||||
t1=sparse(length(Loadi),size(PVi,1)+size(PGi,1));
|
|
||||||
t2=diag(-2./(PD(Loadi).^2)-2*(QD(Loadi).^2-PD(Loadi).^2)./(PD(Loadi).^2+QD(Loadi).^2))*diag(Init_Z(rePos1:rePos2)+Init_W(rePos1:rePos2));
|
|
||||||
t3=4*diag(PD(Loadi).*QD(Loadi)./((PD(Loadi).^2+QD(Loadi).^2).^2))*diag(Init_Z(rePos1:rePos2)+Init_W(rePos1:rePos2));
|
|
||||||
t4=sparse(length(Loadi),2*Busnum);
|
|
||||||
ddg2=[t1,t2,t3,t4];
|
|
||||||
t2=diag(-2*(PD(Loadi).^2-QD(Loadi).^2)./(PD(Loadi).^2+QD(Loadi).^2))*diag(Init_Z(rePos1:rePos2)+Init_W(rePos1:rePos2));
|
|
||||||
ddg3=[t1,t3,t2,t4];
|
|
||||||
ddg4=sparse(Busnum*2,ContrlCount);
|
|
||||||
ddgZW=[ddg1;ddg2;ddg3;ddg4];
|
|
||||||
end
|
|
||||||
|
|
||||||
28
Lineloss.asv
28
Lineloss.asv
|
|
@ -1,28 +0,0 @@
|
||||||
%% 计算线损
|
|
||||||
function Lineloss(Linei,Linej,Liner,Linex,Lineb2,Transi,Transj,Transr,Transx,k0,Volt,Angle)
|
|
||||||
%format long
|
|
||||||
fprintf('功率为有名值\n');
|
|
||||||
fprintf('节点号\t节点号\t有功')
|
|
||||||
cmpVolt=Volt.*cos(Angle)+1i*Volt.*sin(Angle);
|
|
||||||
cmpVolt=cmpVolt';
|
|
||||||
y0=1i*Lineb2;
|
|
||||||
%yj0=1i*standardinput(:,7);
|
|
||||||
yij=1./(Liner+1i*Linex);
|
|
||||||
%% 线路损耗
|
|
||||||
Sij=cmpVolt(Linei).*conj( cmpVolt(Linei) .* y0 + ( cmpVolt(Linei)- cmpVolt(Linej) ).*yij );
|
|
||||||
Sji=cmpVolt(Linej).*conj( cmpVolt(Linej) .*y0 + ( cmpVolt(Linej)- cmpVolt(Linei) ).*yij );
|
|
||||||
%Sij1==Sij2
|
|
||||||
deltLineS=Sij+Sji;
|
|
||||||
dispLineloss=[Linei Linej real(deltLineS)*100 imag(deltLineS)*100];
|
|
||||||
%full(dispLineloss)
|
|
||||||
dispLineloss=sortrows(dispLineloss,-3);
|
|
||||||
full(dispLineloss)
|
|
||||||
%% 以下是变压器损耗
|
|
||||||
yij=1./(Transr+1i*Transx);
|
|
||||||
Sij=cmpVolt(Transi)./k0.*conj( ( cmpVolt(Transi)./k0- cmpVolt(Transj) ).*yij );
|
|
||||||
Sji=cmpVolt(Transj).*conj( ( cmpVolt(Transj)- cmpVolt(Transi)./k0 ).*yij );
|
|
||||||
deltTransS=Sij+Sji;
|
|
||||||
dispTransloss=[Transi Transj real(deltTransS)*100 imag(deltTransS)*100];
|
|
||||||
dispTransloss=sortrows(dispTransloss,-3);
|
|
||||||
full(dispTransloss)
|
|
||||||
end
|
|
||||||
|
|
@ -1,28 +0,0 @@
|
||||||
function [Init_Z,Init_L,Init_W,Init_U,Init_Y,PG,QG,Volt,UAngel,PD,QD]=Modification(Init_Z,Init_L,Init_W,Init_U,Init_Y,deltZ,deltL,deltW,deltU,deltX,deltY,PG,QG,Volt,UAngel,PVi,ContrlCount,Balance,Busnum,PGi,PD,QD,Loadi)
|
|
||||||
AlphaP=FormAlphaP(Init_L,deltL,Init_U,deltU);
|
|
||||||
%fprintf('AlphaP %f\n',full(AlphaP));
|
|
||||||
AlphaD=FormAlphaD(Init_Z,deltZ,Init_W,deltW);
|
|
||||||
%fprintf('AlphaD %f\n',full(AlphaD));
|
|
||||||
|
|
||||||
Init_Z=Init_Z+AlphaD*deltZ';
|
|
||||||
Init_L=Init_L+AlphaP*deltL';
|
|
||||||
Init_W=Init_W+AlphaD*deltW';
|
|
||||||
Init_U=Init_U+AlphaP*deltU';
|
|
||||||
Init_Y=Init_Y+AlphaD*deltY';
|
|
||||||
%PG(PGi)=PG(PGi)+deltX(1:size(PGi,1));
|
|
||||||
PG(PGi)=PG(PGi)+AlphaP*deltX(1:size(PGi,1));
|
|
||||||
%QG(PVi)=QG(PVi)+deltX(size(PGi,1)+1:size(PVi,1)+size(PGi,1) );
|
|
||||||
QG(PVi)=QG(PVi)+AlphaP*deltX(size(PGi,1)+1:size(PVi,1)+size(PGi,1) );
|
|
||||||
t=deltX(size(PVi,1)+size(PGi,1)+1:size(PVi,1)+size(PGi,1)+size(Loadi,1)*2);
|
|
||||||
PD(Loadi)=PD(Loadi)+AlphaP*t(1:length(Loadi));
|
|
||||||
%PD(PD<0)=-PD(PD<0);
|
|
||||||
QD(Loadi)=QD(Loadi)+AlphaP*t(length(Loadi)+1:length(Loadi)*2);
|
|
||||||
t=deltX(size(PVi,1)+size(PGi,1)+size(Loadi,1)*2+1:ContrlCount)';
|
|
||||||
t(Busnum+Balance)=0;
|
|
||||||
%Volt=Volt+AlphaP*t(2:2:2*Busnum);ÔÝʱ¸ÄÒ»ÏÂ20111227
|
|
||||||
%UAngel=UAngel+AlphaP*t(1:2:2*Busnum);ÔÝʱ¸ÄÒ»ÏÂ20111227
|
|
||||||
balVolt=Volt(Balance);
|
|
||||||
Volt=Volt+AlphaP*t(1:Busnum);
|
|
||||||
Volt(Balance)=balVolt;
|
|
||||||
UAngel=UAngel+AlphaP*t(Busnum+1:2*Busnum);
|
|
||||||
end
|
|
||||||
101
OPF.asv
101
OPF.asv
|
|
@ -1,101 +0,0 @@
|
||||||
tic
|
|
||||||
clc
|
|
||||||
clear
|
|
||||||
[kmax,Precision,UAngel,Volt,Busnum,PVi,PVu,Balance,Y,Angle,P0,Q0,r,c,GB,Linei,Linej,Transfori,Transforj,GenU,GenL,GenC,PG,QG,PD,QD,CenterA,PGi,PVQU,PVQL,Liner,Linex,Lineb,Transforr,Transforx,Transfork0]= ...
|
|
||||||
pf('ieee4.dat');
|
|
||||||
%pf('D:\Project\青秀降损项目\最小化潮流\最小潮流算例\原始\津头站津视922(3-1)_0.5_120%.txt');
|
|
||||||
%pf('D:\Project\最小化潮流\最小潮流算例\仙海919.txt');
|
|
||||||
%pf('c:/file31.txt');
|
|
||||||
|
|
||||||
|
|
||||||
%% 计算功率因数
|
|
||||||
%atan(PD(QD~=0 | PD~=0)./QD(QD~=0 | PD~=0));
|
|
||||||
Volt;
|
|
||||||
UAngel*180/3.1415926;
|
|
||||||
%% 通过潮流计算PG
|
|
||||||
AngleIJ=sparse(r,c,UAngel(r)-UAngel(c)-Angle',Busnum,Busnum);
|
|
||||||
PGBal=diag(Volt)*Y.*cos(AngleIJ)*Volt';
|
|
||||||
QGBal=diag(Volt)*Y.*sin(AngleIJ)*Volt';
|
|
||||||
%% 初值-即测量值
|
|
||||||
PG0=PG;
|
|
||||||
QG0=QG;
|
|
||||||
PD0=PD;
|
|
||||||
QD0=QD;
|
|
||||||
PDReal=PD;%真值
|
|
||||||
%PD0(12)=PD0(12)+0.001;
|
|
||||||
%%
|
|
||||||
PG0(Balance)=PGBal(Balance);
|
|
||||||
QG0(Balance)=QGBal(Balance);
|
|
||||||
%%
|
|
||||||
[Volt,UAngel,Init_Z,Init_W,Init_L,Init_U,Init_Y,PG,QG,RestraintCount,wPG,wQG,wD,PD,PD0,QD,randPDind,Loadi]=OPF_Init(Busnum,Balance,PG,QG,Volt,GenU,GenL,PVi,PGi,PVQU,PVQL,PD0,QD0,QD,PD);
|
|
||||||
Gap=(Init_L*Init_Z'-Init_U*Init_W');
|
|
||||||
KK=0;
|
|
||||||
plotGap=zeros(1,50);
|
|
||||||
ContrlCount=size(PVi,1)+size(PGi,1)+size(Loadi,1)*2+Busnum*2;
|
|
||||||
kmax=60;
|
|
||||||
%% 20120523 临时
|
|
||||||
QD_NON_ZERO=QD(PD==0 & QD~=0);
|
|
||||||
QD_NON_ZERO_IND=find(PD==0 & QD~=0);
|
|
||||||
%%
|
|
||||||
while(abs(Gap)>Precision)
|
|
||||||
if KK>kmax
|
|
||||||
break;
|
|
||||||
end
|
|
||||||
plotGap(KK+1)=Gap;
|
|
||||||
Init_u=Gap/2/RestraintCount*CenterA;
|
|
||||||
AngleIJMat=0;
|
|
||||||
%% 开始计算OPF
|
|
||||||
%% 形成等式约束的雅克比
|
|
||||||
deltH=func_deltH(Busnum,Volt,PVi,Y,PGi,UAngel,r,c,Angle,Loadi);
|
|
||||||
%% 形成不等式约束的雅克比
|
|
||||||
deltG=func_deltG(Busnum,PVi,PGi,Loadi);
|
|
||||||
%%
|
|
||||||
L_1Z=diag(Init_Z./Init_L);
|
|
||||||
U_1W=diag(Init_W./Init_U);
|
|
||||||
%% 形成海森阵
|
|
||||||
deltdeltF=func_deltdeltF(PVi,wPG,wQG,wD,ContrlCount);
|
|
||||||
%% 形成ddHy
|
|
||||||
ddh=func_ddh(Volt,Init_Y,Busnum,PVi,PGi,Y,UAngel,r,c,Angle,Loadi,ContrlCount);
|
|
||||||
%% 开始构建ddg
|
|
||||||
ddg=func_ddg(PGi,PVi,Busnum,RestraintCount,Loadi);
|
|
||||||
%% 开始构建deltF
|
|
||||||
deltF=func_deltF(PG,QG,PVi,PGi,wPG,wQG,wD,PG0,QG0,PD0,PD,Busnum,Loadi);
|
|
||||||
|
|
||||||
%% 形成方程矩阵
|
|
||||||
Luu=Init_U'.*Init_W'+Init_u*ones(RestraintCount,1);
|
|
||||||
Lul=Init_L'.*Init_Z'-Init_u*ones(RestraintCount,1);
|
|
||||||
Mat_G=FormG(Volt,PVi,PGi,PG,QG,PD,QD,Loadi);
|
|
||||||
Mat_H=FormH(Busnum,Volt,PG,PD,QG,QD,Y,UAngel,r,c,Angle,QD_NON_ZERO,QD_NON_ZERO_IND);
|
|
||||||
Ly=Mat_H;
|
|
||||||
Lz=FormLz(Mat_G,Init_L,GenL,Busnum,PVQL,PD0,QD0,Loadi,KK);
|
|
||||||
Lw=FormLw(Mat_G,Init_U,GenU,Busnum,PVQU,PD0,QD0,Loadi,KK);
|
|
||||||
Lx=FormLx(deltF,deltH,Init_Y,deltG,Init_Z,Init_W);
|
|
||||||
YY=FormYY(Lul,Lz,Ly,Luu,Lw,Lx);
|
|
||||||
%% 开始解方程
|
|
||||||
fprintf(')
|
|
||||||
XX=SolveIt(deltF,deltG,Init_L,Init_Z,Init_U,Init_W,deltdeltF,ddh,ddg,deltH,Init_Y,Ly,Lz,ContrlCount,Lw,Lul,Luu,RestraintCount,Lx,Balance,PVi,PGi,Busnum,Loadi);
|
|
||||||
%%取各分量
|
|
||||||
[deltZ,deltL,deltW,deltU,deltX,deltY]=AssignXX(XX,ContrlCount,RestraintCount,Busnum);
|
|
||||||
[Init_Z,Init_L,Init_W,Init_U,Init_Y,PG,QG,Volt,UAngel,PD,QD]=Modification(Init_Z,Init_L,Init_W,Init_U,Init_Y,deltZ,deltL,deltW,deltU,deltX,deltY,PG,QG,Volt,UAngel,PVi,ContrlCount,Balance,Busnum,PGi,PD,QD,Loadi);
|
|
||||||
Gap=(Init_L*Init_Z'-Init_U*Init_W');
|
|
||||||
KK=KK+1;
|
|
||||||
end
|
|
||||||
fprintf('迭代次数%d\n',KK);
|
|
||||||
%fprintf('目标值%f\n',full(ObjectiveFun(PG,PG0,PGi,PD,PD0,wG,wD,Loadi)));
|
|
||||||
DrawGap(plotGap);
|
|
||||||
%%
|
|
||||||
%Volt=full(Volt');
|
|
||||||
%PD=full(PD);
|
|
||||||
%% 统计PD误差
|
|
||||||
% absPDLoad=abs( (PD(Loadi)-PD0(Loadi))./PD0(Loadi) );
|
|
||||||
absPDLoad=abs( (PD(Loadi)-PDReal(Loadi))./PDReal(Loadi) );
|
|
||||||
maxPDError=max(absPDLoad);
|
|
||||||
disp('index');
|
|
||||||
Loadi(absPDLoad==maxPDError);
|
|
||||||
%% 计算总线损
|
|
||||||
totalLoss=(sum(PG)-sum(PD(Loadi)))*100;
|
|
||||||
fprintf('总的损耗为%f(MW 有名值)\n',full(totalLoss));
|
|
||||||
%% 计算各线损
|
|
||||||
Lineloss(Linei,Linej,Liner,Linex,Lineb,Transfori,Transforj,Transforr,Transforx,Transfork0,Volt,UAngel);
|
|
||||||
toc
|
|
||||||
|
|
||||||
121
OPF.m
121
OPF.m
|
|
@ -1,121 +0,0 @@
|
||||||
tic
|
|
||||||
clc
|
|
||||||
clear
|
|
||||||
[kmax,Precision,UAngel,Volt,Busnum,PVi,PVu,Balance,Y,Angle,P0,Q0,r,c,GB,Linei,Linej,Transfori,Transforj,GenU,GenL,GenC,PG,QG,PD,QD,CenterA,PGi,PVQU,PVQL,Liner,Linex,Lineb,Transforr,Transforx,Transfork0]= ...
|
|
||||||
pf('ieee14.dat');
|
|
||||||
%pf('D:\Project\青秀降损项目\最小化潮流\最小潮流算例\原始\津头站津视922(3-1)_0.5_120%.txt');
|
|
||||||
%pf('D:\Project\最小化潮流\最小潮流算例\仙海919.txt');
|
|
||||||
%pf('c:/file31.txt');
|
|
||||||
|
|
||||||
Loadi=QD~=0 | PD~=0;
|
|
||||||
%% 计算功率因数
|
|
||||||
PF=sqrt(PD(Loadi).^2./(QD(Loadi).^2+PD(Loadi).^2))
|
|
||||||
Volt;
|
|
||||||
UAngel*180/3.1415926;
|
|
||||||
%% 通过潮流计算PG
|
|
||||||
AngleIJ=sparse(r,c,UAngel(r)-UAngel(c)-Angle',Busnum,Busnum);
|
|
||||||
PGBal=PD+diag(Volt)*Y.*cos(AngleIJ)*Volt';
|
|
||||||
QGBal=QD+diag(Volt)*Y.*sin(AngleIJ)*Volt';
|
|
||||||
%% 初值-即测量值
|
|
||||||
PG0=PG;
|
|
||||||
QG0=QG;
|
|
||||||
PD0=PD;
|
|
||||||
QD0=QD;
|
|
||||||
PDReal=PD;%真值
|
|
||||||
QDReal=QD;%真值
|
|
||||||
%PD0(12)=PD0(12)+0.001;
|
|
||||||
%%
|
|
||||||
PG0(Balance)=PGBal(Balance);
|
|
||||||
QG0(Balance)=QGBal(Balance);
|
|
||||||
QG0(PVi)=QGBal(PVi);
|
|
||||||
PG(Balance)=PGBal(Balance);
|
|
||||||
QG(PVi)=QGBal(PVi);
|
|
||||||
%%
|
|
||||||
AngleIJ=sparse(r,c,UAngel(r)-UAngel(c)-Angle',Busnum,Busnum);
|
|
||||||
dP=PG-PD-diag(Volt)*Y.*cos(AngleIJ)*Volt';
|
|
||||||
dPD=abs(dP./PD);
|
|
||||||
dQ=QG-QD-diag(Volt)*Y.*sin(AngleIJ)*Volt';
|
|
||||||
dQD=abs(dQ./QD);
|
|
||||||
maxdPQ=max([dPD(dPD<10);dQD(dQD<10)])
|
|
||||||
|
|
||||||
%%
|
|
||||||
[Volt,UAngel,Init_Z,Init_W,Init_L,Init_U,Init_Y,PG,QG,RestraintCount,wPG,wQG,wPD,wQD,PD,PD0,QD,randPDind,Loadi]=OPF_Init(Busnum,Balance,PG,QG,Volt,GenU,GenL,PVi,PGi,PVQU,PVQL,PD0,QD0,QD,PD);
|
|
||||||
Gap=(Init_L*Init_Z'-Init_U*Init_W');
|
|
||||||
KK=0;
|
|
||||||
plotGap=zeros(1,60);
|
|
||||||
ContrlCount=size(PVi,1)+size(PGi,1)+size(Loadi,1)*2+Busnum*2;
|
|
||||||
kmax=150;
|
|
||||||
%% 20120523 临时
|
|
||||||
QD_NON_ZERO=QD(PD==0 & QD~=0);
|
|
||||||
QD_NON_ZERO_IND=find(PD==0 & QD~=0);
|
|
||||||
%%
|
|
||||||
while(abs(Gap)>Precision)
|
|
||||||
if KK>kmax
|
|
||||||
break;
|
|
||||||
end
|
|
||||||
plotGap(KK+1)=Gap;
|
|
||||||
Init_u=Gap/2/RestraintCount*CenterA;
|
|
||||||
AngleIJMat=0;
|
|
||||||
%% 开始计算OPF
|
|
||||||
%% 形成等式约束的雅克比
|
|
||||||
deltH=func_deltH(Busnum,Volt,PVi,Y,PGi,UAngel,r,c,Angle,Loadi);
|
|
||||||
%% 形成不等式约束的雅克比
|
|
||||||
deltG=func_deltG(Busnum,PVi,PGi,Loadi,PD,QD);
|
|
||||||
%%
|
|
||||||
L_1Z=diag(Init_Z./Init_L);
|
|
||||||
U_1W=diag(Init_W./Init_U);
|
|
||||||
%% 形成海森阵
|
|
||||||
deltdeltF=func_deltdeltF(PVi,wPG,wQG,wPD,wQD,ContrlCount);
|
|
||||||
%% 形成ddHy
|
|
||||||
ddh=func_ddh(Volt,Init_Y,Busnum,PVi,PGi,Y,UAngel,r,c,Angle,Loadi,ContrlCount);
|
|
||||||
%% 开始构建ddg
|
|
||||||
ddg=func_ddg(PGi,PVi,Busnum,ContrlCount,RestraintCount,Loadi,PD,QD,Init_Z,Init_W);
|
|
||||||
%% 开始构建deltF
|
|
||||||
deltF=func_deltF(PG,QG,PVi,PGi,wPG,wQG,wPD,wQD,PG0,QG0,PD0,PD,QD0,QD,Busnum,Loadi);
|
|
||||||
|
|
||||||
%% 形成方程矩阵
|
|
||||||
Luu=Init_U'.*Init_W'+Init_u*ones(RestraintCount,1);
|
|
||||||
Lul=Init_L'.*Init_Z'-Init_u*ones(RestraintCount,1);
|
|
||||||
Mat_G=FormG(Volt,PVi,PGi,PG,QG,PD,QD,Loadi);
|
|
||||||
Mat_H=FormH(Busnum,Volt,PG,PD,QG,QD,Y,UAngel,r,c,Angle,QD_NON_ZERO,QD_NON_ZERO_IND,Loadi);
|
|
||||||
Ly=Mat_H;
|
|
||||||
Lz=FormLz(Mat_G,Init_L,GenL,Busnum,PVQL,PD0,QD0,Loadi,KK,PF);
|
|
||||||
Lw=FormLw(Mat_G,Init_U,GenU,Busnum,PVQU,PD0,QD0,Loadi,KK,PF);
|
|
||||||
Lx=FormLx(deltF,deltH,Init_Y,deltG,Init_Z,Init_W);
|
|
||||||
YY=FormYY(Lul,Lz,Ly,Luu,Lw,Lx);
|
|
||||||
%% 开始解方程
|
|
||||||
fprintf('迭代次数 %d Gap %f\n',KK+1,plotGap(KK+1));
|
|
||||||
XX=SolveIt(deltF,deltG,Init_L,Init_Z,Init_U,Init_W,deltdeltF,ddh,ddg,deltH,Init_Y,Ly,Lz,ContrlCount,Lw,Lul,Luu,RestraintCount,Lx,Balance,PVi,PGi,Busnum,Loadi);
|
|
||||||
%%取各分量
|
|
||||||
[deltZ,deltL,deltW,deltU,deltX,deltY]=AssignXX(XX,ContrlCount,RestraintCount,Busnum);
|
|
||||||
[Init_Z,Init_L,Init_W,Init_U,Init_Y,PG,QG,Volt,UAngel,PD,QD]=Modification(Init_Z,Init_L,Init_W,Init_U,Init_Y,deltZ,deltL,deltW,deltU,deltX,deltY,PG,QG,Volt,UAngel,PVi,ContrlCount,Balance,Busnum,PGi,PD,QD,Loadi);
|
|
||||||
Gap=(Init_L*Init_Z'-Init_U*Init_W');
|
|
||||||
KK=KK+1;
|
|
||||||
end
|
|
||||||
fprintf('迭代次数%d\n',KK);
|
|
||||||
fprintf('目标值%f\n',full(ObjectiveFun(PG,PG0,PGi,QG,QG0,PVi,PD,PD0,QD,QD0,wPG,wQG,wPD,wQD,Loadi)));
|
|
||||||
DrawGap(plotGap);
|
|
||||||
%%
|
|
||||||
%Volt=full(Volt');
|
|
||||||
%PD=full(PD);
|
|
||||||
%% 统计PD误差
|
|
||||||
% absPDLoad=abs( (PD(Loadi)-PD0(Loadi))./PD0(Loadi) );
|
|
||||||
absPDLoad=abs( (PD(Loadi)-PDReal(Loadi))./PDReal(Loadi) );
|
|
||||||
absQDLoad=abs( (QD(Loadi)-QDReal(Loadi))./QDReal(Loadi) );
|
|
||||||
maxPDError=max([absPDLoad])
|
|
||||||
maxQDError=max([absQDLoad])
|
|
||||||
disp('index');
|
|
||||||
Loadi(absPDLoad==maxPDError);
|
|
||||||
%% 计算总线损
|
|
||||||
totalLoss=(sum(PG)-sum(PD(Loadi)))*100;
|
|
||||||
fprintf('总的损耗为%f(MW 有名值)\n',full(totalLoss));
|
|
||||||
%% 计算各线损
|
|
||||||
Lineloss(Linei,Linej,Liner,Linex,Lineb,Transfori,Transforj,Transforr,Transforx,Transfork0,Volt,UAngel);
|
|
||||||
toc
|
|
||||||
AngleIJ=sparse(r,c,UAngel(r)-UAngel(c)-Angle',Busnum,Busnum);
|
|
||||||
dP=PG-PD-diag(Volt)*Y.*cos(AngleIJ)*Volt';
|
|
||||||
dPD=abs(dP./PD);
|
|
||||||
dQ=QG-QD-diag(Volt)*Y.*sin(AngleIJ)*Volt';
|
|
||||||
dQD=abs(dQ./QD);
|
|
||||||
maxdPQ=max([dPD(dPD<10);dQD(dQD<10)])
|
|
||||||
PF=sqrt(PD(Loadi).^2./(QD(Loadi).^2+PD(Loadi).^2))
|
|
||||||
|
|
@ -1,106 +0,0 @@
|
||||||
clc
|
|
||||||
|
|
||||||
clear
|
|
||||||
for II=1:16
|
|
||||||
tic
|
|
||||||
[kmax,Precision,UAngel,Volt,Busnum,PVi,PVu,Balance,Y,Angle,P0,Q0,r,c,GB,Linei,Linej,Transfori,Transforj,GenU,GenL,GenC,PG,QG,PD,QD,CenterA,PGi,PVQU,PVQL]= ...
|
|
||||||
pf('D:\Project\最小化潮流\最小潮流算例\金湖924(2-1)_0.5_85%.txt');
|
|
||||||
%pf('c:/file31.txt');
|
|
||||||
%pf('ieee118PG.dat');
|
|
||||||
|
|
||||||
%% 计算功率因数
|
|
||||||
atan(PD(QD~=0 | PD~=0)./QD(QD~=0 | PD~=0));
|
|
||||||
Volt;
|
|
||||||
UAngel*180/3.1415926;
|
|
||||||
%% 通过潮流计算PG
|
|
||||||
AngleIJ=sparse(r,c,UAngel(r)-UAngel(c)-Angle',Busnum,Busnum);
|
|
||||||
PGBal=diag(Volt)*Y.*cos(AngleIJ)*Volt';
|
|
||||||
|
|
||||||
%% 初值-即测量值
|
|
||||||
PG0=PG;
|
|
||||||
PD0=PD;
|
|
||||||
|
|
||||||
PDReal=PD;%真值
|
|
||||||
|
|
||||||
%%
|
|
||||||
PG0(Balance)=PGBal(Balance);
|
|
||||||
%%
|
|
||||||
[Volt,UAngel,Init_Z,Init_W,Init_L,Init_U,Init_Y,PG,QG,RestraintCount,wG,wD,PD,PD0,randPDind,Loadi]=OPF_Init(Busnum,Balance,PG,QG,Volt,GenU,GenL,PVi,PGi,PVQU,PVQL,PD0,QD,PD);
|
|
||||||
PD0(Loadi(II))=PD0(Loadi(II))*(1+0.018);
|
|
||||||
Gap=(Init_L*Init_Z'-Init_U*Init_W');
|
|
||||||
KK=0;
|
|
||||||
plotGap=zeros(1,50);
|
|
||||||
ContrlCount=size(PVi,1)+size(PGi,1)+size(Loadi,1)+Busnum*2;
|
|
||||||
kmax=60;
|
|
||||||
%% 20120523 临时
|
|
||||||
QD_NON_ZERO=QD(PD==0 & QD~=0);
|
|
||||||
QD_NON_ZERO_IND=find(PD==0 & QD~=0);
|
|
||||||
%%
|
|
||||||
while(abs(Gap)>Precision)
|
|
||||||
if KK>kmax
|
|
||||||
break;
|
|
||||||
end
|
|
||||||
plotGap(KK+1)=Gap;
|
|
||||||
Init_u=Gap/2/RestraintCount*CenterA;
|
|
||||||
AngleIJMat=0;
|
|
||||||
%% 开始计算OPF
|
|
||||||
%% 形成等式约束的雅克比
|
|
||||||
deltH=func_deltH(Busnum,Volt,PVi,Y,PGi,UAngel,r,c,Angle,Loadi);
|
|
||||||
%% 形成不等式约束的雅克比
|
|
||||||
deltG=func_deltG(Busnum,PVi,PGi,Loadi);
|
|
||||||
%%
|
|
||||||
L_1Z=diag(Init_Z./Init_L);
|
|
||||||
U_1W=diag(Init_W./Init_U);
|
|
||||||
%% 形成海森阵
|
|
||||||
deltdeltF=func_deltdeltF(PVi,wG,wD,ContrlCount);
|
|
||||||
%% 形成ddHy
|
|
||||||
ddh=func_ddh(Volt,Init_Y,Busnum,PVi,PGi,Y,UAngel,r,c,Angle,Loadi,ContrlCount);
|
|
||||||
%% 开始构建ddg
|
|
||||||
ddg=func_ddg(PGi,PVi,Busnum,RestraintCount,Loadi);
|
|
||||||
%% 开始构建deltF
|
|
||||||
deltF=func_deltF(PG,PVi,PGi,wG,wD,PG0,PD0,PD,Busnum,Loadi);
|
|
||||||
|
|
||||||
%% 形成方程矩阵
|
|
||||||
Luu=Init_U'.*Init_W'+Init_u*ones(RestraintCount,1);
|
|
||||||
Lul=Init_L'.*Init_Z'-Init_u*ones(RestraintCount,1);
|
|
||||||
Mat_G=FormG(Volt,PVi,PGi,PG,QG,PD,Loadi);
|
|
||||||
Mat_H=FormH(Busnum,Volt,PG,PD,QG,QD,Y,UAngel,r,c,Angle,QD_NON_ZERO,QD_NON_ZERO_IND);
|
|
||||||
Ly=Mat_H;
|
|
||||||
Lz=FormLz(Mat_G,Init_L,GenL,Busnum,PVQL,PD0,Loadi,KK);
|
|
||||||
Lw=FormLw(Mat_G,Init_U,GenU,Busnum,PVQU,PD0,Loadi,KK);
|
|
||||||
Lx=FormLx(deltF,deltH,Init_Y,deltG,Init_Z,Init_W);
|
|
||||||
YY=FormYY(Lul,Lz,Ly,Luu,Lw,Lx);
|
|
||||||
%% 开始解方程
|
|
||||||
XX=SolveIt(deltF,deltG,Init_L,Init_Z,Init_U,Init_W,deltdeltF,ddh,ddg,deltH,Init_Y,Ly,Lz,ContrlCount,Lw,Lul,Luu,RestraintCount,Lx,Balance,PVi,PGi,Busnum,Loadi);
|
|
||||||
%%取各分量
|
|
||||||
[deltZ,deltL,deltW,deltU,deltX,deltY]=AssignXX(XX,ContrlCount,RestraintCount,Busnum);
|
|
||||||
[Init_Z,Init_L,Init_W,Init_U,Init_Y,PG,QG,Volt,UAngel,PD]=Modification(Init_Z,Init_L,Init_W,Init_U,Init_Y,deltZ,deltL,deltW,deltU,deltX,deltY,PG,QG,Volt,UAngel,PVi,ContrlCount,Balance,Busnum,PGi,PD,Loadi);
|
|
||||||
Gap=(Init_L*Init_Z'-Init_U*Init_W');
|
|
||||||
KK=KK+1;
|
|
||||||
end
|
|
||||||
fprintf('迭代次数%d\n',KK);
|
|
||||||
ObjectiveFun(PG,PG0,PGi,PD,PD0,wG,wD,Loadi);
|
|
||||||
%DrawGap(plotGap);
|
|
||||||
%%
|
|
||||||
%Volt=full(Volt');
|
|
||||||
%PD=full(PD);
|
|
||||||
%% 统计PD误差
|
|
||||||
% absPDLoad=abs( (PD(Loadi)-PD0(Loadi))./PD0(Loadi) );
|
|
||||||
absPDLoad=abs( (PD(Loadi)-PDReal(Loadi))./PDReal(Loadi) );
|
|
||||||
maxPDError=max(absPDLoad);
|
|
||||||
%disp('index')
|
|
||||||
LoadiArray=Loadi(absPDLoad==maxPDError);
|
|
||||||
if length(LoadiArray)>1
|
|
||||||
disp('没找出')
|
|
||||||
Loadi(II)
|
|
||||||
end
|
|
||||||
if length(LoadiArray)==1
|
|
||||||
if LoadiArray~=Loadi(II)
|
|
||||||
disp('没找出')
|
|
||||||
Loadi(II)
|
|
||||||
end
|
|
||||||
end
|
|
||||||
toc;
|
|
||||||
|
|
||||||
end
|
|
||||||
|
|
||||||
|
|
@ -1,106 +0,0 @@
|
||||||
clc
|
|
||||||
|
|
||||||
clear
|
|
||||||
for II=1:53
|
|
||||||
tic
|
|
||||||
[kmax,Precision,UAngel,Volt,Busnum,PVi,PVu,Balance,Y,Angle,P0,Q0,r,c,GB,Linei,Linej,Transfori,Transforj,GenU,GenL,GenC,PG,QG,PD,QD,CenterA,PGi,PVQU,PVQL]= ...
|
|
||||||
pf('D:\Project\最小化潮流\最小潮流算例\仙海919.txt');
|
|
||||||
%pf('c:/file31.txt');
|
|
||||||
%pf('ieee118PG.dat');
|
|
||||||
|
|
||||||
%% 计算功率因数
|
|
||||||
atan(PD(QD~=0 | PD~=0)./QD(QD~=0 | PD~=0));
|
|
||||||
Volt;
|
|
||||||
UAngel*180/3.1415926;
|
|
||||||
%% 通过潮流计算PG
|
|
||||||
AngleIJ=sparse(r,c,UAngel(r)-UAngel(c)-Angle',Busnum,Busnum);
|
|
||||||
PGBal=diag(Volt)*Y.*cos(AngleIJ)*Volt';
|
|
||||||
|
|
||||||
%% 初值-即测量值
|
|
||||||
PG0=PG;
|
|
||||||
PD0=PD;
|
|
||||||
|
|
||||||
PDReal=PD;%真值
|
|
||||||
|
|
||||||
%%
|
|
||||||
PG0(Balance)=PGBal(Balance);
|
|
||||||
%%
|
|
||||||
[Volt,UAngel,Init_Z,Init_W,Init_L,Init_U,Init_Y,PG,QG,RestraintCount,wG,wD,PD,PD0,randPDind,Loadi]=OPF_Init(Busnum,Balance,PG,QG,Volt,GenU,GenL,PVi,PGi,PVQU,PVQL,PD0,QD,PD);
|
|
||||||
PD0(Loadi(II))=PD0(Loadi(II))*(1+0.086);
|
|
||||||
Gap=(Init_L*Init_Z'-Init_U*Init_W');
|
|
||||||
KK=0;
|
|
||||||
plotGap=zeros(1,50);
|
|
||||||
ContrlCount=size(PVi,1)+size(PGi,1)+size(Loadi,1)+Busnum*2;
|
|
||||||
kmax=60;
|
|
||||||
%% 20120523 临时
|
|
||||||
QD_NON_ZERO=QD(PD==0 & QD~=0);
|
|
||||||
QD_NON_ZERO_IND=find(PD==0 & QD~=0);
|
|
||||||
%%
|
|
||||||
while(abs(Gap)>Precision)
|
|
||||||
if KK>kmax
|
|
||||||
break;
|
|
||||||
end
|
|
||||||
plotGap(KK+1)=Gap;
|
|
||||||
Init_u=Gap/2/RestraintCount*CenterA;
|
|
||||||
AngleIJMat=0;
|
|
||||||
%% 开始计算OPF
|
|
||||||
%% 形成等式约束的雅克比
|
|
||||||
deltH=func_deltH(Busnum,Volt,PVi,Y,PGi,UAngel,r,c,Angle,Loadi);
|
|
||||||
%% 形成不等式约束的雅克比
|
|
||||||
deltG=func_deltG(Busnum,PVi,PGi,Loadi);
|
|
||||||
%%
|
|
||||||
L_1Z=diag(Init_Z./Init_L);
|
|
||||||
U_1W=diag(Init_W./Init_U);
|
|
||||||
%% 形成海森阵
|
|
||||||
deltdeltF=func_deltdeltF(PVi,wG,wD,ContrlCount);
|
|
||||||
%% 形成ddHy
|
|
||||||
ddh=func_ddh(Volt,Init_Y,Busnum,PVi,PGi,Y,UAngel,r,c,Angle,Loadi,ContrlCount);
|
|
||||||
%% 开始构建ddg
|
|
||||||
ddg=func_ddg(PGi,PVi,Busnum,RestraintCount,Loadi);
|
|
||||||
%% 开始构建deltF
|
|
||||||
deltF=func_deltF(PG,PVi,PGi,wG,wD,PG0,PD0,PD,Busnum,Loadi);
|
|
||||||
|
|
||||||
%% 形成方程矩阵
|
|
||||||
Luu=Init_U'.*Init_W'+Init_u*ones(RestraintCount,1);
|
|
||||||
Lul=Init_L'.*Init_Z'-Init_u*ones(RestraintCount,1);
|
|
||||||
Mat_G=FormG(Volt,PVi,PGi,PG,QG,PD,Loadi);
|
|
||||||
Mat_H=FormH(Busnum,Volt,PG,PD,QG,QD,Y,UAngel,r,c,Angle,QD_NON_ZERO,QD_NON_ZERO_IND);
|
|
||||||
Ly=Mat_H;
|
|
||||||
Lz=FormLz(Mat_G,Init_L,GenL,Busnum,PVQL,PD0,Loadi,KK);
|
|
||||||
Lw=FormLw(Mat_G,Init_U,GenU,Busnum,PVQU,PD0,Loadi,KK);
|
|
||||||
Lx=FormLx(deltF,deltH,Init_Y,deltG,Init_Z,Init_W);
|
|
||||||
YY=FormYY(Lul,Lz,Ly,Luu,Lw,Lx);
|
|
||||||
%% 开始解方程
|
|
||||||
XX=SolveIt(deltF,deltG,Init_L,Init_Z,Init_U,Init_W,deltdeltF,ddh,ddg,deltH,Init_Y,Ly,Lz,ContrlCount,Lw,Lul,Luu,RestraintCount,Lx,Balance,PVi,PGi,Busnum,Loadi);
|
|
||||||
%%取各分量
|
|
||||||
[deltZ,deltL,deltW,deltU,deltX,deltY]=AssignXX(XX,ContrlCount,RestraintCount,Busnum);
|
|
||||||
[Init_Z,Init_L,Init_W,Init_U,Init_Y,PG,QG,Volt,UAngel,PD]=Modification(Init_Z,Init_L,Init_W,Init_U,Init_Y,deltZ,deltL,deltW,deltU,deltX,deltY,PG,QG,Volt,UAngel,PVi,ContrlCount,Balance,Busnum,PGi,PD,Loadi);
|
|
||||||
Gap=(Init_L*Init_Z'-Init_U*Init_W');
|
|
||||||
KK=KK+1;
|
|
||||||
end
|
|
||||||
fprintf('迭代次数%d\n',KK);
|
|
||||||
ObjectiveFun(PG,PG0,PGi,PD,PD0,wG,wD,Loadi);
|
|
||||||
%DrawGap(plotGap);
|
|
||||||
%%
|
|
||||||
%Volt=full(Volt');
|
|
||||||
%PD=full(PD);
|
|
||||||
%% 统计PD误差
|
|
||||||
% absPDLoad=abs( (PD(Loadi)-PD0(Loadi))./PD0(Loadi) );
|
|
||||||
absPDLoad=abs( (PD(Loadi)-PDReal(Loadi))./PDReal(Loadi) );
|
|
||||||
maxPDError=max(absPDLoad);
|
|
||||||
%disp('index')
|
|
||||||
LoadiArray=Loadi(absPDLoad==maxPDError);
|
|
||||||
if length(LoadiArray)>1
|
|
||||||
disp('没找出')
|
|
||||||
Loadi(II)
|
|
||||||
end
|
|
||||||
if length(LoadiArray)==1
|
|
||||||
if LoadiArray~=Loadi(II)
|
|
||||||
disp('没找出')
|
|
||||||
Loadi(II)
|
|
||||||
end
|
|
||||||
end
|
|
||||||
toc;
|
|
||||||
|
|
||||||
end
|
|
||||||
|
|
||||||
46
OPF_Init.asv
46
OPF_Init.asv
|
|
@ -1,46 +0,0 @@
|
||||||
function [Volt,UAngel,Init_Z,Init_W,Init_L,Init_U,Init_Y,PG,QG,RestraintCount,wPG,wQG,wPD,wQD,PD,PD0,QD,randPDind,Loadi,notLoadi]=OPF_Init(Busnum,Balance,PG,QG,Volt,GenU,GenL,PVi,PGi,PVQU,PVQL,PD0,QD0,QD,PD)
|
|
||||||
Loadi=find(QD~=0 | PD~=0);
|
|
||||||
PDi=find(PD~=0);
|
|
||||||
QDi=find(QD~=0);
|
|
||||||
notLoadi=setxor(1:Busnum,Loadi);
|
|
||||||
PGOnly=setdiff(PGi,PDi);
|
|
||||||
QGOnly=setdiff(PVi,QDi);
|
|
||||||
PDOnly=setdiff(PDi,PGi);
|
|
||||||
QDOnly=setdiff(QDi,PVi);
|
|
||||||
PGPD=intersec(PDi,PGi);
|
|
||||||
QGQD=intersec(QDi,PVi);
|
|
||||||
noPGPDQGQD=setdiff(1:B)
|
|
||||||
%Loadi=[1:Busnum]';
|
|
||||||
RestraintCount=size(PVi,1)+size(PGi,1)+size(Loadi,1)*1+Busnum*1; %约束条件数,放开所有QD,不单独约束QD上下限
|
|
||||||
RestraintCount=RestraintCount+length(Loadi); %加上功率因数,用功率因数约束QD
|
|
||||||
t_Bal_volt=Volt(Balance);
|
|
||||||
Volt=sparse(1*ones(1,Busnum));
|
|
||||||
Volt(Balance)=t_Bal_volt;
|
|
||||||
UAngel=sparse(1,Busnum);
|
|
||||||
Init_Z=sparse(ones(1,RestraintCount));
|
|
||||||
Init_W=sparse(-1*ones(1,RestraintCount));
|
|
||||||
Init_L=1*sparse(ones(1,RestraintCount));
|
|
||||||
Init_U=1*sparse(ones(1,RestraintCount));
|
|
||||||
Init_Y=sparse(1,2*Busnum);%与学姐一致
|
|
||||||
tPU=sparse(GenU(:,2));% 发电机有功上限
|
|
||||||
tQU=sparse(PVQU(:,1));% 无功上限
|
|
||||||
tPL=sparse(GenL(:,2));% 发电机有功下限
|
|
||||||
tQL=sparse(PVQL(:,1));% 无功下限
|
|
||||||
PG(PGi)=(tPU+tPL)/2;
|
|
||||||
QG(PVi)=(tQU+tQL)/2;
|
|
||||||
wPG=1*ones(size(PGi,1),1);
|
|
||||||
wQG=1*ones(size(PVi,1),1);
|
|
||||||
%randInt=randperm(size(Loadi,1));
|
|
||||||
%randPDind=randInt(1:10);
|
|
||||||
randPDind=0;
|
|
||||||
wPD=1*ones(size(Loadi,1),1);
|
|
||||||
wQD=1.5*ones(size(Loadi,1),1);
|
|
||||||
%wD(randPDind)=0;%一些负荷不约束
|
|
||||||
%wD(7)=0;
|
|
||||||
% wD(11)=0;
|
|
||||||
PD=1*PD0;
|
|
||||||
%powerFacter=0.98;
|
|
||||||
%QD=PD*sqrt((1-powerFacter^2)/powerFacter^2);
|
|
||||||
QD=QD0;
|
|
||||||
|
|
||||||
end
|
|
||||||
|
|
@ -1,19 +0,0 @@
|
||||||
function [out_arg]=ObjectiveFun(PG,PG0,PGi,QG,QG0,PVi,PD,PD0,QD,QD0,wPG,wQG,wPD,wQD,Loadi)
|
|
||||||
% t1=[PG(PGi).value]'-PG0(PGi);
|
|
||||||
% t2=wPG.*t1.*t1;
|
|
||||||
t2=sdpvar(1,1);
|
|
||||||
for I=1:length(PGi)
|
|
||||||
t2=t2+wPG(I)*(PG(I).value-PG0(I));
|
|
||||||
end
|
|
||||||
|
|
||||||
% t1=[QG(PVi).value]'-QG0(PVi);
|
|
||||||
% t3=wQG.*t1.*t1;
|
|
||||||
t3=sdpvar(1,1);
|
|
||||||
for I=1:length(PVi)
|
|
||||||
t3=t3+wPG(I)*(QG(I).value-QG0(I));
|
|
||||||
end
|
|
||||||
% t4=wPD.*(([PD(Loadi).value]'-PD0(Loadi)).*([PD(Loadi).value]'-PD0(Loadi)));
|
|
||||||
for I=1:length(L)
|
|
||||||
% t5=wQD.*(([QD(Loadi).value]'-QD0(Loadi)).*([QD(Loadi).value]'-QD0(Loadi)));
|
|
||||||
out_arg= sum(t2)+sum(t3)+sum(t4)+sum(t5);
|
|
||||||
end
|
|
||||||
|
|
@ -1 +0,0 @@
|
||||||
对照学姐给的公式
|
|
||||||
|
|
@ -5,7 +5,7 @@ tic
|
||||||
[kmax,Precision,UAngel,Volt,Busnum,PVi,PVu,Balance,Y,Angle,P0,Q0,r,c,GB, ...
|
[kmax,Precision,UAngel,Volt,Busnum,PVi,PVu,Balance,Y,Angle,P0,Q0,r,c,GB, ...
|
||||||
Linei,Linej,Transfori,Transforj,GenU,GenL,GenC,PG,QG,PD,QD,CenterA,PGi,PVQU,PVQL, ...
|
Linei,Linej,Transfori,Transforj,GenU,GenL,GenC,PG,QG,PD,QD,CenterA,PGi,PVQU,PVQL, ...
|
||||||
Liner,Linex,Lineb,Transforr,Transforx,Transfork0]= ...
|
Liner,Linex,Lineb,Transforr,Transforx,Transfork0]= ...
|
||||||
pf('ieee4.dat');
|
pf('ieee30.dat');
|
||||||
|
|
||||||
%% 潮流等式
|
%% 潮流等式
|
||||||
AngleIJ=sparse(r,c,UAngel(r)-UAngel(c)-Angle',Busnum,Busnum);
|
AngleIJ=sparse(r,c,UAngel(r)-UAngel(c)-Angle',Busnum,Busnum);
|
||||||
|
|
@ -37,40 +37,52 @@ maxdPQ=max([dPD(dPD<10);dQD(dQD<10)]);
|
||||||
%% 定义变量
|
%% 定义变量
|
||||||
Volt=sdpvar(Busnum,1);
|
Volt=sdpvar(Busnum,1);
|
||||||
UAngel=sdpvar(Busnum,1);
|
UAngel=sdpvar(Busnum,1);
|
||||||
PG=sdpvar(Busnum,1);
|
% PG=sdpvar(Busnum,1);
|
||||||
QG=sdpvar(Busnum,1);
|
% QG=sdpvar(Busnum,1);
|
||||||
PD=sdpvar(Busnum,1);
|
PD=sdpvar(Busnum,1);
|
||||||
QD=sdpvar(Busnum,1);
|
QD=sdpvar(Busnum,1);
|
||||||
%AngleIJ=sdpvar(Busnum,Busnum,'full');
|
AngleIJ=sdpvar(Busnum,Busnum,'full');
|
||||||
%% 目标函数
|
%% 目标函数
|
||||||
Objective=ObjectiveFun(PG,PG0,PGi,QG,QG0,PVi,PD,PD0,QD,QD0,wPG,wQG,wPD,wQD,Loadi);
|
Objective=ObjectiveFun(PG,PG0,PGi,QG,QG0,PVi,PD,PD0,QD,QD0,wPG,wQG,wPD,wQD,Loadi);
|
||||||
AngleIJ=sparse(r,c,UAngel(r)-UAngel(c)-Angle,Busnum,Busnum);
|
%AngleIJ=sparse(r,c,UAngel(r)-UAngel(c)-Angle,Busnum,Busnum);
|
||||||
|
%% 赋初值,可以加快求解速度。
|
||||||
|
assign(Volt(:),1);
|
||||||
|
assign(UAngel(:),0);
|
||||||
|
assign(PD(:),PD0(:));
|
||||||
|
assign(QD(:),QD0(:));
|
||||||
|
% assign(PG(:),PG0(:));
|
||||||
|
% assign(QG(:),QG0(:));
|
||||||
%% YALMIP部分
|
%% YALMIP部分
|
||||||
Constraints = [%AngleIJ==sparse(r,c,UAngel(r)-UAngel(c)-Angle,Busnum,Busnum), ...
|
dP=PG0-PD-diag(Volt)*Y.*cos( sparse(r,c,UAngel(r)-UAngel(c)-Angle,Busnum,Busnum) )*Volt;
|
||||||
AngleIJ==0
|
dQ=QG0-QD-diag(Volt)*Y.*sin( sparse(r,c,UAngel(r)-UAngel(c)-Angle,Busnum,Busnum) )*Volt;
|
||||||
PG-PD-diag(Volt)*Y.*cos(AngleIJ)*Volt==0, ...
|
Loadi=PD0~=0 | QD0~=0;
|
||||||
QG-QD-diag(Volt)*Y.*sin(AngleIJ)*Volt==0, ...
|
Constraints = [%AngleIJ-sparse(r,c,UAngel(r)-UAngel(c)-Angle,Busnum,Busnum)==0, ...
|
||||||
|
dP==0, ...
|
||||||
|
dQ==0, ...
|
||||||
PG(setxor(1:Busnum, PGi) )==0, ...
|
PG(setxor(1:Busnum, PGi) )==0, ...
|
||||||
QG(setxor(1:Busnum, PVi) )==0, ...
|
QG(setxor(1:Busnum, PVi) )==0, ...
|
||||||
PD(setdiff(1:Busnum,PD0~=0))==0, ...
|
PD(PD0==0)==0, ...
|
||||||
QD(setdiff(1:Busnum,QD0~=0))==0, ...
|
QD(QD0==0)==0, ...
|
||||||
0.9*ones(Busnum,1)<=Volt<=1.1*ones(Busnum,1), ...
|
0.9*ones(Busnum,1)<=Volt<=1.1*ones(Busnum,1), ...
|
||||||
GenL(:,2)<=PG(PGi)<=GenU(:,2), ...
|
-60*ones(Busnum,1)<=PG<=60*ones(Busnum,1), ...
|
||||||
PVQL<=QG(PVi)<=PVQU, ...
|
-60*ones(Busnum,1)<=QG<=60*ones(Busnum,1)
|
||||||
];
|
];
|
||||||
options = sdpsettings('verbose',2,'showprogress',1,'debug',0,'solver','ipopt');
|
options = sdpsettings('verbose',2,'showprogress',1,'debug',0,'solver','ipopt','usex0','1');
|
||||||
sol = solvesdp(Constraints,Objective,options);
|
sol = solvesdp(Constraints,Objective,options);
|
||||||
if sol.problem == 0
|
if sol.problem == 0
|
||||||
fprintf('Volt\n');
|
fprintf('Volt\n');
|
||||||
dvolt=double(Volt);
|
dvolt=double(Volt)
|
||||||
fprintf('VoltAngle\n');
|
fprintf('VoltAngle\n');
|
||||||
dVangle=double(UAngel);
|
dVangle=double(UAngel)
|
||||||
fprintf('ojb\n');
|
fprintf('ojb\n');
|
||||||
optimalObj=double(Objective)
|
optimalObj=double(Objective)
|
||||||
|
double(PG)-PG0;
|
||||||
sol
|
sol
|
||||||
else
|
else
|
||||||
display('Hmm, something went wrong!');
|
display('Hmm, something went wrong!');
|
||||||
sol.info
|
sol.info
|
||||||
|
sol.solveroutput
|
||||||
yalmiperror(sol.problem)
|
yalmiperror(sol.problem)
|
||||||
end
|
end
|
||||||
|
|
||||||
toc
|
toc
|
||||||
|
|
|
||||||
36
Run_YALMIP.m
36
Run_YALMIP.m
|
|
@ -37,42 +37,52 @@ maxdPQ=max([dPD(dPD<10);dQD(dQD<10)]);
|
||||||
%% ¶¨Òå±äÁ¿
|
%% ¶¨Òå±äÁ¿
|
||||||
Volt=sdpvar(Busnum,1);
|
Volt=sdpvar(Busnum,1);
|
||||||
UAngel=sdpvar(Busnum,1);
|
UAngel=sdpvar(Busnum,1);
|
||||||
PG=sdpvar(Busnum,1);
|
% PG=sdpvar(Busnum,1);
|
||||||
QG=sdpvar(Busnum,1);
|
% QG=sdpvar(Busnum,1);
|
||||||
PD=sdpvar(Busnum,1);
|
PD=sdpvar(Busnum,1);
|
||||||
QD=sdpvar(Busnum,1);
|
QD=sdpvar(Busnum,1);
|
||||||
AngleIJ=sdpvar(Busnum,Busnum,'full');
|
AngleIJ=sdpvar(Busnum,Busnum,'full');
|
||||||
%% Ä¿±êº¯Êý
|
%% Ä¿±êº¯Êý
|
||||||
Objective=ObjectiveFun(PG,PG0,PGi,QG,QG0,PVi,PD,PD0,QD,QD0,wPG,wQG,wPD,wQD,Loadi);
|
Objective=ObjectiveFun(PG,PG0,PGi,QG,QG0,PVi,PD,PD0,QD,QD0,wPG,wQG,wPD,wQD,Loadi);
|
||||||
%AngleIJ=sparse(r,c,UAngel(r)-UAngel(c)-Angle,Busnum,Busnum);
|
%AngleIJ=sparse(r,c,UAngel(r)-UAngel(c)-Angle,Busnum,Busnum);
|
||||||
|
%% 赋初值,可以加快求解速度。
|
||||||
|
assign(Volt(:),1);
|
||||||
|
assign(UAngel(:),0);
|
||||||
|
assign(PD(:),PD0(:));
|
||||||
|
assign(QD(:),QD0(:));
|
||||||
|
% assign(PG(:),PG0(:));
|
||||||
|
% assign(QG(:),QG0(:));
|
||||||
%% YALMIP²¿·Ö
|
%% YALMIP²¿·Ö
|
||||||
Constraints = [AngleIJ==sparse(r,c,UAngel(r)-UAngel(c)-Angle,Busnum,Busnum), ...
|
dP=PG0-PD-diag(Volt)*Y.*cos( sparse(r,c,UAngel(r)-UAngel(c)-Angle,Busnum,Busnum) )*Volt;
|
||||||
PG-PD-diag(Volt)*Y.*cos(AngleIJ)*Volt==0, ...
|
dQ=QG0-QD-diag(Volt)*Y.*sin( sparse(r,c,UAngel(r)-UAngel(c)-Angle,Busnum,Busnum) )*Volt;
|
||||||
QG-QD-diag(Volt)*Y.*sin(AngleIJ)*Volt==0, ...
|
Loadi=PD0~=0 | QD0~=0;
|
||||||
|
Constraints = [%AngleIJ-sparse(r,c,UAngel(r)-UAngel(c)-Angle,Busnum,Busnum)==0, ...
|
||||||
|
dP==0, ...
|
||||||
|
dQ==0, ...
|
||||||
PG(setxor(1:Busnum, PGi) )==0, ...
|
PG(setxor(1:Busnum, PGi) )==0, ...
|
||||||
QG(setxor(1:Busnum, PVi) )==0, ...
|
QG(setxor(1:Busnum, PVi) )==0, ...
|
||||||
PD(PD0==0)==0, ...
|
PD(PD0==0)==0, ...
|
||||||
QD(QD0==0)==0, ...
|
QD(QD0==0)==0, ...
|
||||||
0.9*ones(Busnum,1)<=Volt<=1.1*ones(Busnum,1), ...
|
0.9*ones(Busnum,1)<=Volt<=1.1*ones(Busnum,1), ...
|
||||||
%GenL(:,2)<=PG(PGi)<=GenU(:,2), ...
|
-60*ones(Busnum,1)<=PG<=60*ones(Busnum,1), ...
|
||||||
%PVQL<=QG(PVi)<=PVQU, ...
|
-60*ones(Busnum,1)<=QG<=60*ones(Busnum,1)
|
||||||
-20*ones(Busnum,1)<=PG<=20*ones(Busnum,1), ...
|
|
||||||
-20*ones(Busnum,1)<=QG<=20*ones(Busnum,1)
|
|
||||||
];
|
];
|
||||||
options = sdpsettings('verbose',2,'showprogress',1,'debug',0,'solver','ipopt');
|
options = sdpsettings('verbose',2,'showprogress',1,'debug',0,'solver','ipopt','usex0','1');
|
||||||
sol = solvesdp(Constraints,Objective,options);
|
sol = solvesdp(Constraints,Objective,options);
|
||||||
if sol.problem == 0
|
if sol.problem == 0
|
||||||
fprintf('Volt\n');
|
fprintf('Volt\n');
|
||||||
dvolt=double(Volt);
|
dvolt=double(Volt)
|
||||||
fprintf('VoltAngle\n');
|
fprintf('VoltAngle\n');
|
||||||
dVangle=double(UAngel);
|
dVangle=double(UAngel)
|
||||||
fprintf('ojb\n');
|
fprintf('ojb\n');
|
||||||
optimalObj=double(Objective)
|
optimalObj=double(Objective)
|
||||||
double(PG)-PG0
|
double(PG)-PG0;
|
||||||
sol
|
sol
|
||||||
else
|
else
|
||||||
display('Hmm, something went wrong!');
|
display('Hmm, something went wrong!');
|
||||||
sol.info
|
sol.info
|
||||||
|
sol.solveroutput
|
||||||
yalmiperror(sol.problem)
|
yalmiperror(sol.problem)
|
||||||
end
|
end
|
||||||
|
|
||||||
toc
|
toc
|
||||||
|
|
|
||||||
42
SolveIt.asv
42
SolveIt.asv
|
|
@ -1,42 +0,0 @@
|
||||||
function XX=SolveIt(deltF,deltG,Init_L,Init_Z,Init_U,Init_W,deltdeltF,ddh,ddg,deltH,Init_Y,Ly,Lz,ContrlCount,Lw,Lul,Luu,RestraintCount,Lx,Balance,PVi,PGi,Busnum,Loadi)
|
|
||||||
LxComa=FormLxComa(deltF,deltG,deltH,Init_L,Luu,Lul,Init_Z,Init_Y,Lz,Init_U,Init_W,Lw,Lx);
|
|
||||||
H=-deltdeltF+ddh;%+ddg*(Init_Z'+Init_W');
|
|
||||||
t1=diag(Init_L.\Init_Z-Init_U.\Init_W);
|
|
||||||
t2=-deltG*( t1 )*deltG';
|
|
||||||
aa=[
|
|
||||||
(H+t2),deltH;
|
|
||||||
deltH',zeros(size(Init_Y,2));
|
|
||||||
];
|
|
||||||
yy=[LxComa;-Ly];
|
|
||||||
%% ƽºâ½Úµãµçѹ²»±ä
|
|
||||||
t=size(PVi,1)+size(PGi,1)+size(Loadi,1);
|
|
||||||
aa(t+Balance,:)=0;
|
|
||||||
aa(:,t+Balance)=0;
|
|
||||||
%aa(t+Balance,t+Balance)=1;
|
|
||||||
aa=aa+sparse(t+Balance,t+Balance,ones(length(Balance)),ContrlCount+2*Busnum,ContrlCount+2*Busnum);
|
|
||||||
deltG(t+Balance,:)=0;
|
|
||||||
%%
|
|
||||||
t=size(PVi,1)+size(PGi,1)+size(Loadi,1)+Busnum*1;
|
|
||||||
aa(t+Balance,:)=0;
|
|
||||||
aa(:,t+Balance)=0;
|
|
||||||
%aa(t+Balance,t+Balance)=1;
|
|
||||||
aa=aa+sparse(t+Balance,t+Balance,ones(length(Balance)),ContrlCount,ContrlCount);
|
|
||||||
deltG(t+Balance,:)=0;
|
|
||||||
%%
|
|
||||||
dxdy=aa\yy;
|
|
||||||
dX=dxdy(1:ContrlCount);
|
|
||||||
dY=dxdy(ContrlCount+1:ContrlCount+2*Busnum);
|
|
||||||
dL=Lz+deltG'*dX;
|
|
||||||
dU=-Lw-deltG'*dX;
|
|
||||||
dZ=-diag(Init_L)\Lul-diag(Init_L)\diag(Init_Z)*dL;
|
|
||||||
dW=-diag(Init_U)\Luu-diag(Init_U)\diag(Init_W)*dU;
|
|
||||||
XX=[
|
|
||||||
dX;
|
|
||||||
dY;
|
|
||||||
dZ;
|
|
||||||
dW;
|
|
||||||
dL;
|
|
||||||
dU;
|
|
||||||
|
|
||||||
];
|
|
||||||
end
|
|
||||||
51
SolveIt.m
51
SolveIt.m
|
|
@ -1,51 +0,0 @@
|
||||||
function XX=SolveIt(deltF,deltG,Init_L,Init_Z,Init_U,Init_W,deltdeltF,ddh,ddg,deltH,Init_Y,Ly,Lz,ContrlCount,Lw,Lul,Luu,RestraintCount,Lx,Balance,PVi,PGi,Busnum,Loadi)
|
|
||||||
LxComa=FormLxComa(deltF,deltG,deltH,Init_L,Luu,Lul,Init_Z,Init_Y,Lz,Init_U,Init_W,Lw,Lx);
|
|
||||||
H=-deltdeltF+ddh+ddg;
|
|
||||||
%t1=diag(Init_L.\Init_Z-Init_U.\Init_W);
|
|
||||||
t1=diag(Init_Z./Init_L-Init_W./Init_U);
|
|
||||||
t2=-deltG*( t1 )*deltG';
|
|
||||||
aa=[
|
|
||||||
(H+t2),deltH;
|
|
||||||
deltH',zeros(size(Init_Y,2));
|
|
||||||
];
|
|
||||||
yy=[LxComa;-Ly];
|
|
||||||
|
|
||||||
%% ƽºâ½Úµãµçѹ²»±ä
|
|
||||||
t=size(PVi,1)+size(PGi,1)+size(Loadi,1)*2;
|
|
||||||
yy(t+Balance)=0;
|
|
||||||
aa(t+Balance,:)=0;
|
|
||||||
aa(:,t+Balance)=0;
|
|
||||||
%aa(t+Balance,t+Balance)=1;
|
|
||||||
aa=aa+sparse(t+Balance,t+Balance,ones(length(Balance),1),ContrlCount+2*Busnum,ContrlCount+2*Busnum);
|
|
||||||
deltG(t+Balance,:)=0;
|
|
||||||
%%
|
|
||||||
t=size(PVi,1)+size(PGi,1)+size(Loadi,1)*2+Busnum*1;
|
|
||||||
yy(t+Balance)=0;
|
|
||||||
aa(t+Balance,:)=0;
|
|
||||||
aa(:,t+Balance)=0;
|
|
||||||
%aa(t+Balance,t+Balance)=1;
|
|
||||||
aa=aa+sparse(t+Balance,t+Balance,ones(length(Balance),1),ContrlCount+2*Busnum,ContrlCount+2*Busnum);
|
|
||||||
deltG(t+Balance,:)=0;
|
|
||||||
%%
|
|
||||||
|
|
||||||
dxdy=aa\yy;
|
|
||||||
%% KLU
|
|
||||||
%spy(aa)
|
|
||||||
%dxdy = klu(aa,'\',full(yy));
|
|
||||||
%%
|
|
||||||
dX=dxdy(1:ContrlCount);
|
|
||||||
dY=dxdy(ContrlCount+1:ContrlCount+2*Busnum);
|
|
||||||
dL=Lz+deltG'*dX;
|
|
||||||
dU=-Lw-deltG'*dX;
|
|
||||||
dZ=-diag(Init_L)\Lul-diag(Init_L)\diag(Init_Z)*dL;
|
|
||||||
dW=-diag(Init_U)\Luu-diag(Init_U)\diag(Init_W)*dU;
|
|
||||||
XX=[
|
|
||||||
dX;
|
|
||||||
dY;
|
|
||||||
dZ;
|
|
||||||
dW;
|
|
||||||
dL;
|
|
||||||
dU;
|
|
||||||
|
|
||||||
];
|
|
||||||
end
|
|
||||||
|
|
@ -1,31 +0,0 @@
|
||||||
function [G,B,GB,Y,r,c,Angle] = admmatrix(Busnum,Linei,Linej,Liner,Linex,Lineb,Transfori...
|
|
||||||
,Transforj,Transforr,Transforx,Transfork0,Branchi,Branchb)
|
|
||||||
%**************************************************************************
|
|
||||||
% 程序功能 : 子函数——形成节点导纳矩阵Y
|
|
||||||
% 编 者:
|
|
||||||
% 编制时间:2010.12
|
|
||||||
%**************************************************************************
|
|
||||||
%% 支路导纳计算
|
|
||||||
G = -sparse(Linei,Linej,Liner./(Liner.^2+Linex.^2),Busnum,Busnum) - sparse(Linej,Linei,Liner./(Liner.^2+Linex.^2),Busnum,Busnum);
|
|
||||||
G = G - sparse(1:Busnum,1:Busnum,sum(G,2)'); % 计算各线路支路电导
|
|
||||||
B = sparse(Linei,Linej,Linex./(Liner.^2+Linex.^2),Busnum,Busnum)+sparse(Linej,Linei,Linex./(Liner.^2+Linex.^2),Busnum,Busnum);
|
|
||||||
B = B - sparse(1:Busnum,1:Busnum,sum(B,2)')+sparse(Linei,Linei,Lineb,Busnum,Busnum)+sparse(Linej,Linej,Lineb,Busnum,Busnum);
|
|
||||||
%% 变压器支路计算
|
|
||||||
if Transfori>0
|
|
||||||
mr = Transforr./(Transforr.^2+Transforx.^2); % 计算变压器支路电导
|
|
||||||
mx = -Transforx./(Transforr.^2+Transforx.^2); % 计算变压器支路电纳
|
|
||||||
G = G-sparse(Transfori,Transforj,mr./Transfork0,Busnum,Busnum)-sparse(Transforj,Transfori,mr./Transfork0,Busnum,Busnum)...
|
|
||||||
+sparse(Transfori,Transfori,mr./Transfork0./Transfork0,Busnum,Busnum)+sparse(Transforj,Transforj,mr,Busnum,Busnum);
|
|
||||||
B = B-sparse(Transfori,Transforj,mx./Transfork0,Busnum,Busnum)-sparse(Transforj,Transfori,mx./Transfork0,Busnum,Busnum)...
|
|
||||||
+sparse(Transfori,Transfori,mx./Transfork0./Transfork0,Busnum,Busnum)+sparse(Transforj,Transforj,mx,Busnum,Busnum);
|
|
||||||
end
|
|
||||||
%% 接地支路计算
|
|
||||||
if Branchi>0 % 判断有无接地支路
|
|
||||||
B = B+sparse(Branchi,Branchi,Branchb,Busnum,Busnum);
|
|
||||||
end
|
|
||||||
%% 化作极坐标形式
|
|
||||||
GB = G+B.*1i; %将电导,电纳合并,写成复数形式
|
|
||||||
Y = abs(GB); %求节点导纳幅值
|
|
||||||
[r,c] = find(Y);
|
|
||||||
Angle = angle(GB(GB~=0)); %求节点导纳角度
|
|
||||||
%Angle=angle(GB);
|
|
||||||
176
func_ddf.m
176
func_ddf.m
|
|
@ -1,176 +0,0 @@
|
||||||
function ddf=func_ddh(AngleIJMat,GB,Volt,Init_Y,Busnum)
|
|
||||||
%% deltaPi/deltaThytai_deltaThytaj 非对角元素
|
|
||||||
t1=-Volt'*Volt;
|
|
||||||
t2=real(GB).*cos(AngleIJMat)+imag(GB).*sin(AngleIJMat);
|
|
||||||
dPdTidTj=t1.*t2; %%(保留了对角元素的)
|
|
||||||
dPidTidTj_2=dPdTidTj-diag(diag(dPdTidTj));%去掉了对角元素的
|
|
||||||
t3=repmat(Init_Y',1,size(Init_Y,2));%乘y的系数
|
|
||||||
dPidTidTj_2=dPidTidTj_2.*t3(1:2:size(Init_Y,2),1:2:size(Init_Y,2));
|
|
||||||
t2=real(GB).*cos(AngleIJMat)-imag(GB).*sin(AngleIJMat);
|
|
||||||
t3=t1.*t2;
|
|
||||||
dPjdTidTj=t3-diag(diag(t3));%去掉了对角元素的
|
|
||||||
t3=repmat(Init_Y,size(Init_Y,2),1);
|
|
||||||
dPjdTidTj=dPjdTidTj.*t3(1:2:size(Init_Y,2),1:2:size(Init_Y,2));%乘y的系数
|
|
||||||
dPdTidTj=dPidTidTj_2+dPjdTidTj;%最终非对角元素 @@
|
|
||||||
%% deltaP/deltaThyta_deltaThyta 对角元素
|
|
||||||
t1=sum(-dPidTidTj_2,2);
|
|
||||||
t2=diag(t1'.*Init_Y(1:2:size(Init_Y,2)));%乘y的系数
|
|
||||||
t3=sum(-dPidTidTj_2,1);
|
|
||||||
t4=diag(t3.*Init_Y(1:2:size(Init_Y,2)));%乘y的系数
|
|
||||||
dPdTidTi=t2+t4;%%最终对角元素 @@
|
|
||||||
%% deltaP/deltaThytai_dVi 对角元素
|
|
||||||
t1=ones(Busnum,1)*Volt.*(real(GB).*sin(AngleIJMat)-imag(GB).*cos(AngleIJMat));
|
|
||||||
t2=diag(t1);
|
|
||||||
t3=t1-diag(t2);%去掉了对角元素的
|
|
||||||
t4=sum(t3,2);
|
|
||||||
t4=t4'.*Init_Y(1:2:size(Init_Y,2));%乘y的系数
|
|
||||||
dPidTidVi=diag(t4);
|
|
||||||
dPidTjdVj=-t3;
|
|
||||||
t6=sum(dPidTjdVj,1);%乘y的系数
|
|
||||||
t6=t6.*Init_Y(1:2:size(Init_Y,2));%乘y的系数
|
|
||||||
dPdTidVi=dPidTidVi+diag(t6);%%最终对角元素 @@
|
|
||||||
%% deltaP/deltaThytai_dVj 非对角元素
|
|
||||||
t1=ones(Busnum,1)*Volt.*(real(GB).*sin(AngleIJMat)-imag(GB).*cos(AngleIJMat));
|
|
||||||
dPidTidVj=t1-diag(diag(t1));%%去掉对角元素的
|
|
||||||
t2=repmat(Init_Y',1,size(Init_Y,2));
|
|
||||||
dPidTidVj=dPidTidVj.*t2(1:2:size(Init_Y,2),1:2:size(Init_Y,2));%乘y的系数
|
|
||||||
t1=ones(Busnum,1)*Volt.*real(GB).*(sin(AngleIJMat)+imag(GB).*cos(AngleIJMat));
|
|
||||||
dPjdTidVj=t1-diag(diag(t1));%去掉对角元素的
|
|
||||||
t2=t2';
|
|
||||||
dPjdTidVj=dPjdTidVj.*t2(1:2:size(Init_Y,2),1:2:size(Init_Y,2));%乘y的系数
|
|
||||||
dPdTidVj=dPidTidVj+dPjdTidVj;%最终非对角元素 @@
|
|
||||||
%% deltaP/dVi_deltaThytaj 非对角元素
|
|
||||||
t1=-ones(Busnum,1)*Volt.*real(GB).*(sin(AngleIJMat)-imag(GB).*cos(AngleIJMat));
|
|
||||||
dPidVidTj=t1-diag(diag(t1));
|
|
||||||
t2=repmat(Init_Y',1,size(Init_Y,2));
|
|
||||||
dPidVidTj=dPidVidTj.*t2(1:2:size(Init_Y,2),1:2:size(Init_Y,2));%乘y的系数
|
|
||||||
t1=-ones(Busnum,1)*Volt.*(real(GB).*sin(AngleIJMat)+imag(GB).*cos(AngleIJMat));
|
|
||||||
dPjdVidTj=t1-diag(diag(t1));
|
|
||||||
t2=t2';
|
|
||||||
dPjdVidTj.*t2(1:2:size(Init_Y,2),1:2:size(Init_Y,2));%乘y的系数
|
|
||||||
dPdVidTj=dPidVidTj+dPjdVidTj;%最终非对角元素
|
|
||||||
%% deltaPi/dVi_deltaThyta 对角元素
|
|
||||||
dPdVidTi=dPdTidVi;%最终对角元素 @@
|
|
||||||
%% deltaP/dVi_dVj 非对角元素
|
|
||||||
t1=-(real(GB).*(cos(AngleIJMat)+imag(GB).*sin(AngleIJMat)));
|
|
||||||
dPidVidVj=t1-diag(diag(t1));%去掉对角元素的
|
|
||||||
t2=repmat(Init_Y',1,size(Init_Y,2));
|
|
||||||
dPidVidVj=dPidVidVj.*t2(1:2:size(Init_Y,2),1:2:size(Init_Y,2));%乘y的系数
|
|
||||||
t1=-(real(GB).*(cos(AngleIJMat)-imag(GB).*sin(AngleIJMat)));
|
|
||||||
dPjdVidVj=t1-diag(diag(t1));
|
|
||||||
t2=t2';
|
|
||||||
dPjdVidVj=dPjdVidVj.*t2(1:2:size(Init_Y,2),1:2:size(Init_Y,2));%乘y的系数
|
|
||||||
dPdVidVj=dPidVidVj+dPjdVidVj;%最终非对角元素 @@
|
|
||||||
%% deltaP/dVi_dVi 对角元素
|
|
||||||
t1=diag(real(GB));
|
|
||||||
t2=t1'.*Init_Y(1:2:size(Init_Y,2));%乘y的系数
|
|
||||||
dPidVidVi=-2*diag(t2);
|
|
||||||
dPidVjdVj=0;
|
|
||||||
dPdVidVi=dPidVidVi+dPidVjdVj;%最终对角元素 @@
|
|
||||||
%% 生成APi
|
|
||||||
APi=zeros(2*Busnum,2*Busnum);
|
|
||||||
APi(1:2:2*Busnum,1:2:2*Busnum)=dPdTidTj;%%非对角 TT
|
|
||||||
APi(1:2:2*Busnum,2:2:2*Busnum)=dPdTidVj;%%非对角 TV
|
|
||||||
APi(2:2:2*Busnum,1:2:2*Busnum)=dPdVidTj;%%非对角 VT
|
|
||||||
APi(2:2:2*Busnum,2:2:2*Busnum)=dPdVidVj;%%非对角 VV
|
|
||||||
APi(1:2:2*Busnum,1:2:2*Busnum)=dPdTidTi;%%对角
|
|
||||||
APi(1:2:2*Busnum,2:2:2*Busnum)=dPdTidVi;%%对角
|
|
||||||
APi(2:2:2*Busnum,1:2:2*Busnum)=dPdVidTi;%%对角
|
|
||||||
APi(2:2:2*Busnum,2:2:2*Busnum)=dPdVidVi;%%对角
|
|
||||||
%% deltaQ/deltaThyta_deltaThyta 非对角元素
|
|
||||||
t1=-Volt'*Volt;
|
|
||||||
t2=real(GB).*(sin(AngleIJMat)-imag(GB).*cos(AngleIJMat));
|
|
||||||
dQidTidTj=t1.*t2;%不去掉对角元素了,反正最后是要修正的
|
|
||||||
t3=repmat(Init_Y',1,size(Init_Y,2));
|
|
||||||
dQidTidTj=dQidTidTj.*t3(2:2:size(Init_Y,2),2:2:size(Init_Y,2));%乘y的系数
|
|
||||||
t2=real(GB).*(sin(AngleIJMat)+imag(GB).*cos(AngleIJMat));
|
|
||||||
dQjdTidTj=t1.*t2;
|
|
||||||
t3=t3';
|
|
||||||
dQjdTidTj=dQjdTidTj.*t3(2:2:size(Init_Y,2),2:2:size(Init_Y,2));%乘y的系数
|
|
||||||
dQdTidTj=dQidTidTj+dQjdTidTj;%最终非对角元素
|
|
||||||
%% deltaQ/deltaThyta_deltaThyta 对角元素
|
|
||||||
t1=dQidTidTj-diag(diag(dQidTidTj));
|
|
||||||
t2=sum(t1,2);
|
|
||||||
t3=t2'.*Init_Y(2:2:size(Init_Y,2));%乘y的系数
|
|
||||||
dQidTidTi=diag(t3);
|
|
||||||
t1=-Volt'*Volt;
|
|
||||||
t2=real(GB).*(sin(AngleIJMat)-imag(GB).*cos(AngleIJMat));
|
|
||||||
t3=t1.*t2;
|
|
||||||
t4=t3-diag(diag(t3));
|
|
||||||
t5=sum(t4,1);
|
|
||||||
t6=t5.*Init_Y(2:2:size(Init_Y,2));%乘y的系数
|
|
||||||
dQjdTidTi=diag(t6);
|
|
||||||
dQdTidTi=dQjdTidTi+dQidTidTi;%最终对角元素
|
|
||||||
%% deltaQ/deltaThyta_deltaV 非对角元素
|
|
||||||
t1=-Volt;
|
|
||||||
t2=real(GB).*(cos(AngleIJMat)+imag(GB).*sin(AngleIJMat));
|
|
||||||
t3=ones(Busnum,1)*t1.*t2;
|
|
||||||
t4=repmat(Init_Y',1,size(Init_Y,2));
|
|
||||||
dQidTidVj=t3.*t4(2:2:size(Init_Y,2),2:2:size(Init_Y,2));%乘y的系数
|
|
||||||
t2=real(GB).*(cos(AngleIJMat)-imag(GB).*sin(AngleIJMat));
|
|
||||||
t3=Volt'*ones(1,Busnum).*t2;
|
|
||||||
t4=t4';
|
|
||||||
dQjdTidVj=t2.*t3.*t4(2:2:size(Init_Y,2),2:2:size(Init_Y,2));%乘y的系数
|
|
||||||
dQdTidVj=dQidTidVj+dQjdTidVj;%最终非对角元素
|
|
||||||
%% deltaQ/deltaThyta_deltaV 对角元素
|
|
||||||
t1=sum(dQidTidVj,2)-diag(dQidTidVj);
|
|
||||||
t2=t1'.*Init_Y(2:2:size(Init_Y,2));%乘y的系数
|
|
||||||
dQidTidVi=diag(t2);
|
|
||||||
t1=-Volt'*ones(1,Busnum).*(real(GB).*cos(AngleIJMat)+imag(GB).*sin(AngleIJMat));
|
|
||||||
t2=t1-diag(diag(t1));
|
|
||||||
t3=sum(t2,1);
|
|
||||||
t4=t3.*Init_Y(2:2:size(Init_Y,2));%乘y的系数
|
|
||||||
dQjdTidVi=diag(t4);
|
|
||||||
dQdTidVi=dQidTidVi+dQjdTidVi;
|
|
||||||
%% deltaQ/deltaV_deltaV 非对角元素
|
|
||||||
t1=real(GB).*(sin(AngleIJMat)-imag(GB).*cos(AngleIJMat));
|
|
||||||
t2=repmat(Init_Y',1,size(Init_Y,2));
|
|
||||||
t3=t1.*t2(2:2:size(Init_Y,2),2:2:size(Init_Y,2));%乘y的系数
|
|
||||||
dQidVidVj=t3;
|
|
||||||
t1=-real(GB).*sin(AngleIJMat)-imag(GB).*cos(AngleIJMat);
|
|
||||||
t2=t2';
|
|
||||||
t3=t1.*t2(2:2:size(Init_Y,2),2:2:size(Init_Y,2));%乘y的系数
|
|
||||||
dQjdVidVj=t3;
|
|
||||||
dQdVidVj=dQidVidVj+dQjdVidVj;
|
|
||||||
%% deltaQ/deltaV_deltaV 对角元素
|
|
||||||
t1=2*(real(GB).*(sin(AngleIJMat)-imag(GB).*cos(AngleIJMat)));
|
|
||||||
t2=diag(t1);
|
|
||||||
t3=t2'.*Init_Y(2:2:size(Init_Y,2));%乘y的系数
|
|
||||||
dQidVidVi=diag(t3);
|
|
||||||
t1=-real(GB).*sin(AngleIJMat)-imag(GB).*cos(AngleIJMat);
|
|
||||||
t2=diag(t1);
|
|
||||||
t3=t1-diag(t2);
|
|
||||||
t4=sum(t3,1);
|
|
||||||
t5=t4.*Init_Y(2:2:size(Init_Y,2));%乘y的系数
|
|
||||||
dQjdVidVi=diag(t5);
|
|
||||||
dQdVidVi=dQidVidVi+dQjdVidVi;
|
|
||||||
%% deltaQ/deltaV_deltaThyta 非对角元素
|
|
||||||
t1=real(GB).*sin(AngleIJMat)+imag(GB).*cos(AngleIJMat);
|
|
||||||
t2=-ones(Busnum,1)*Volt.*(t1);
|
|
||||||
t3=repmat(Init_Y',1,size(Init_Y,2));
|
|
||||||
t4=t2'.*t3(2:2:size(Init_Y,2),2:2:size(Init_Y,2));%乘y的系数
|
|
||||||
dQidVidTj=t4;
|
|
||||||
t1=real(GB).*sin(AngleIJMat)-imag(GB).*cos(AngleIJMat);
|
|
||||||
t2=-ones(Busnum,1)*Volt.*(t1);
|
|
||||||
t3=t3';
|
|
||||||
t4=t2.*t3(2:2:size(Init_Y,2),2:2:size(Init_Y,2));%乘y的系数
|
|
||||||
dQjdVidTj=t4;
|
|
||||||
dQdVidTj=dQidVidTj+dQjdVidTj;
|
|
||||||
%% deltaQ/deltaV_deltaThyta 对角元素
|
|
||||||
dQdVidTi=dQdTidVi;
|
|
||||||
%% 生成AQi
|
|
||||||
AQi=zeros(2*Busnum,2*Busnum);
|
|
||||||
AQi(1:2:2*Busnum,1:2:2*Busnum)=dQdTidTj;%%非对角 TT
|
|
||||||
AQi(1:2:2*Busnum,2:2:2*Busnum)=dQdTidVj;%%非对角 TV
|
|
||||||
AQi(2:2:2*Busnum,1:2:2*Busnum)=dQdVidTj;%%非对角 VT
|
|
||||||
AQi(2:2:2*Busnum,2:2:2*Busnum)=dQdVidVj;%%非对角 VV
|
|
||||||
AQi(1:2:2*Busnum,1:2:2*Busnum)=dQdTidTi;%%对角
|
|
||||||
AQi(1:2:2*Busnum,2:2:2*Busnum)=dQdTidVi;%%对角
|
|
||||||
AQi(2:2:2*Busnum,1:2:2*Busnum)=dQdVidTi;%%对角
|
|
||||||
AQi(2:2:2*Busnum,2:2:2*Busnum)=dQdVidVi;%%对角
|
|
||||||
%% 生成ddf
|
|
||||||
t=[zeros(2*14-2*5,2*14);
|
|
||||||
zeros(2*5,2*14-2*5),AQi+APi;
|
|
||||||
];
|
|
||||||
ddf=t;
|
|
||||||
end
|
|
||||||
14
func_ddg.asv
14
func_ddg.asv
|
|
@ -1,14 +0,0 @@
|
||||||
function ddg=func_ddg(PGi,PVi,Busnum,ContrlCount,RestraintCount,Loadi,PD,QD)
|
|
||||||
% ¼ÆËãϵÊý
|
|
||||||
InitZW=Init_Z+InitW;
|
|
||||||
loc=InitZW(length(PGi)+length(PVi)+length(Loadi)*2:RestraintCount);
|
|
||||||
PD_=PD(Loadi);
|
|
||||||
QD_=QD(Loadi)
|
|
||||||
t1=sparse(length(PGi)+length(PVi),ContrlCount);
|
|
||||||
t21=sparse(length(PD_),length(PGi)+length(PVi));
|
|
||||||
t22=sparse(-2./PD_.^2-2*(QD_.^2-PD_.^2)./( PD_.^2+QD_.^2).^2)*diag(loc);%dF/dPdP
|
|
||||||
t23=sparse(4*PD_.*PQ_./(PD_.^2+PQ_.^2).^2); %dF/dPdQ
|
|
||||||
t24=sparse(length(PD_),2*Busnum);
|
|
||||||
t31=t21;
|
|
||||||
t32=sparse(-2.)
|
|
||||||
end
|
|
||||||
23
func_ddg.m
23
func_ddg.m
|
|
@ -1,23 +0,0 @@
|
||||||
function ddg=func_ddg(PGi,PVi,Busnum,ContrlCount,RestraintCount,Loadi,PD,QD,Init_Z,Init_W)
|
|
||||||
% ¼ÆËãϵÊý
|
|
||||||
InitZW=Init_Z+Init_W;
|
|
||||||
loc=InitZW(length(PGi)+length(PVi)+length(Loadi)*1+Busnum+1:RestraintCount);
|
|
||||||
PD_=PD(Loadi);
|
|
||||||
QD_=QD(Loadi);
|
|
||||||
t1=sparse(length(PGi)+length(PVi),ContrlCount);
|
|
||||||
t21=sparse(length(PD_),length(PGi)+length(PVi));
|
|
||||||
t22=diag(sparse(-2./PD_.^2-2*(QD_.^2-PD_.^2)./( PD_.^2+QD_.^2).^2))*diag(loc);%dF/dPdP
|
|
||||||
t23=diag(sparse(4*PD_.*QD_./(PD_.^2+QD_.^2).^2))*diag(loc); %dF/dPdQ
|
|
||||||
t24=sparse(length(PD_),2*Busnum);
|
|
||||||
t31=t21;
|
|
||||||
t32=t23;
|
|
||||||
t33=diag(sparse(-2*(PD_.^2-QD_.^2)./(PD_.^2+QD_.^2).^2 ))*diag(loc);
|
|
||||||
t34=t24;
|
|
||||||
t4=sparse(Busnum*2,ContrlCount);
|
|
||||||
ddg=[
|
|
||||||
t1;
|
|
||||||
t21,t22,t23,t24;
|
|
||||||
t31,t32,t33,t34;
|
|
||||||
t4;
|
|
||||||
];
|
|
||||||
end
|
|
||||||
58
func_ddh.m
58
func_ddh.m
|
|
@ -1,58 +0,0 @@
|
||||||
function ddh=func_ddh(Volt,Init_Y,Busnum,PVi,PGi,Y,UAngel,r,c,Angle,Loadi,ContrlCount)
|
|
||||||
%决定用循环重写
|
|
||||||
%ContrlCount=size(PVi,1)+size(PGi,1)+size(Loadi,1)+Busnum*2;
|
|
||||||
|
|
||||||
%AngleIJ=AngleIJMat-angle(GB);
|
|
||||||
mat_AngleIJ=sparse(r,c,UAngel(r)-UAngel(c)-Angle',Busnum,Busnum);
|
|
||||||
mat_INV_AngleIJ=mat_AngleIJ';
|
|
||||||
yP=Init_Y(1:size(Init_Y,2)/2);%暂时改这里 20111227
|
|
||||||
yQ=Init_Y(size(Init_Y,2)/2+1:size(Init_Y,2));%暂时改这里 20111227
|
|
||||||
t1=-sparse(1:Busnum,1:Busnum,Y.*cos(mat_INV_AngleIJ)*sparse(1:Busnum,1:Busnum,Volt,Busnum,Busnum)*yP',Busnum,Busnum);
|
|
||||||
t2=sparse(1:Busnum,1:Busnum,sparse(1:Busnum,1:Busnum,Volt,Busnum,Busnum)*yP',Busnum,Busnum)*Y.*cos(mat_AngleIJ);
|
|
||||||
t3=(t1+t2)*sparse(1:Busnum,1:Busnum,Volt,Busnum,Busnum);
|
|
||||||
t4=-(sparse(1:Busnum,1:Busnum,Y.*cos(mat_AngleIJ)*Volt',Busnum,Busnum) -sparse(1:Busnum,1:Busnum,Volt,Busnum,Busnum)*Y.*cos(mat_INV_AngleIJ) )*sparse(1:Busnum,1:Busnum,sparse(1:Busnum,1:Busnum,Volt,Busnum,Busnum)*yP',Busnum,Busnum);
|
|
||||||
ddPdTdT=t3+t4;%ok1
|
|
||||||
t1=(-sparse(1:Busnum,1:Busnum,Y.*sin(mat_AngleIJ)*Volt',Busnum,Busnum)+sparse(1:Busnum,1:Busnum,Volt,Busnum,Busnum)*Y.*sin(mat_INV_AngleIJ) )*sparse(1:Busnum,1:Busnum,yP,Busnum,Busnum);
|
|
||||||
t2= -sparse(1:Busnum,1:Busnum, sparse(1:Busnum,1:Busnum,Volt,Busnum,Busnum)*yP' ,Busnum,Busnum)*Y.*sin(mat_AngleIJ)+sparse(1:Busnum,1:Busnum,Y.*sin(mat_INV_AngleIJ)*sparse(1:Busnum,1:Busnum,Volt,Busnum,Busnum)*yP',Busnum,Busnum);
|
|
||||||
ddPdVdT=t1+t2;%ok1
|
|
||||||
t1=diag( Y.*sin(mat_INV_AngleIJ)*diag(Volt)*yP');
|
|
||||||
t2=diag(yP)*Y.*sin(mat_AngleIJ)*diag(Volt);
|
|
||||||
t3=-diag(yP)*diag(Y.*sin(mat_AngleIJ)*Volt');
|
|
||||||
t4=-Y.*sin(mat_INV_AngleIJ)*diag( diag(Volt)*yP' );
|
|
||||||
ddPdTdV=t1+t2+t3+t4;%存疑与我的不一样
|
|
||||||
t1=Y.*cos(mat_INV_AngleIJ)*diag(yP);
|
|
||||||
t2=diag(yP)*Y.*cos(mat_AngleIJ);
|
|
||||||
ddPdVdV=t1+t2;
|
|
||||||
t1=-diag(Y.*sin(mat_AngleIJ)*Volt');
|
|
||||||
t2=diag(Volt)*Y.*sin(mat_INV_AngleIJ);
|
|
||||||
t3=(t1+t2)*diag( diag(Volt)*yQ' );
|
|
||||||
t4=-diag( diag(Volt)*yQ' )*Y.*sin(mat_AngleIJ);
|
|
||||||
|
|
||||||
t5=diag(Y.*sin(mat_INV_AngleIJ)*diag(Volt)*yQ');
|
|
||||||
t6=-(t4+t5)*diag(Volt);
|
|
||||||
ddQdTdT=t3+t6;%ok1
|
|
||||||
t1=(diag(Y.*cos(mat_AngleIJ)*Volt')-diag(Volt)*Y.*cos(mat_INV_AngleIJ) )*diag(yQ);
|
|
||||||
t2=+diag( diag(Volt)*yQ' )*Y.*cos(mat_AngleIJ)-diag(Y.*cos(mat_INV_AngleIJ)*diag(Volt)*yQ');
|
|
||||||
ddQdVdT=t1+t2;
|
|
||||||
t1=Y.*cos(mat_INV_AngleIJ)*diag(diag(Volt)*yQ');
|
|
||||||
t2=diag(yQ)*diag(Y.*cos(mat_AngleIJ)*Volt');
|
|
||||||
t3=-diag(Y.*cos(mat_INV_AngleIJ)*diag(Volt)*yQ');
|
|
||||||
t4=-diag(yQ)*Y.*cos(mat_AngleIJ)*diag(Volt);
|
|
||||||
ddQdTdV=t1+t2+t3+t4;
|
|
||||||
t1=Y.*sin(mat_INV_AngleIJ)*diag(yQ);
|
|
||||||
t2=diag(yQ)*Y.*sin(mat_AngleIJ);
|
|
||||||
ddQdVdV=t1+t2;
|
|
||||||
% t=[ddPdVdV+ddQdVdV,ddPdTdV+ddQdTdV;
|
|
||||||
% ddPdVdT+ddQdVdT,ddPdTdT+ddQdTdT;
|
|
||||||
% ];
|
|
||||||
t=[ddPdVdV+ddQdVdV,ddPdVdT+ddQdVdT;
|
|
||||||
ddPdTdV+ddQdTdV,ddPdTdT+ddQdTdT;
|
|
||||||
];
|
|
||||||
sizePGi=size(PGi,1);
|
|
||||||
sizePVi=size(PVi,1);
|
|
||||||
sizeLoadi=size(Loadi,1)*2;
|
|
||||||
ddh=[
|
|
||||||
sparse(sizePGi+sizePVi+sizeLoadi,ContrlCount);
|
|
||||||
sparse(2*Busnum,sizePVi+sizePGi+sizeLoadi),-t;
|
|
||||||
];
|
|
||||||
end
|
|
||||||
|
|
@ -1,19 +0,0 @@
|
||||||
function deltF=func_deltF(PG,PVi,GenC,ContrlCount,PGi,wG,wD,PG0,PD0,PD,Busnum)
|
|
||||||
%t1=PG(setdiff(PVi,Balance));
|
|
||||||
% t2=Volt'*Volt;
|
|
||||||
% t3=real(GB).*cos(AngleIJMat)+imag(GB).*sin(AngleIJMat);
|
|
||||||
% t4=t2.*t3;
|
|
||||||
% t5=sum(t4,2);
|
|
||||||
% PBal=t5(Balance);
|
|
||||||
% PPG=([PQ(1),PBal])';%暂时用土办法处理一下
|
|
||||||
%%
|
|
||||||
t1=2
|
|
||||||
t2=2*wD.*(PD-PD0);
|
|
||||||
deltF=[
|
|
||||||
zeros(size(PGi));
|
|
||||||
zeros(size(PVi));
|
|
||||||
t2;
|
|
||||||
zeros(2*Busnum,1);
|
|
||||||
];
|
|
||||||
|
|
||||||
end
|
|
||||||
15
func_deltF.m
15
func_deltF.m
|
|
@ -1,15 +0,0 @@
|
||||||
function deltF=func_deltF(PG,QG,PVi,PGi,wPG,wQG,wPD,wQD,PG0,QG0,PD0,PD,QD0,QD,Busnum,Loadi)
|
|
||||||
|
|
||||||
t1=2*wPG.*(PG(PGi) - PG0(PGi) );
|
|
||||||
t2=2*wQG.*(QG(PVi) - QG0(PVi) );
|
|
||||||
t3=2*wPD.*(PD(Loadi)-PD0(Loadi));
|
|
||||||
t4=2*wQD.*(QD(Loadi)-QD0(Loadi));
|
|
||||||
deltF=[
|
|
||||||
sparse(t1);
|
|
||||||
sparse(t2);
|
|
||||||
sparse(t3);
|
|
||||||
sparse(length(Loadi),1);
|
|
||||||
sparse(2*Busnum,1);
|
|
||||||
];
|
|
||||||
|
|
||||||
end
|
|
||||||
|
|
@ -1,52 +0,0 @@
|
||||||
function deltG=func_deltG(Busnum,PVi,PGi,Loadi,PD,QD)
|
|
||||||
%%
|
|
||||||
PD_=PD(Loadi);
|
|
||||||
QD_=QD(Loadi);
|
|
||||||
sizePGi=size(PGi,1);
|
|
||||||
sizePVi=size(PVi,1);
|
|
||||||
sizeLoadi=size(Loadi,1);
|
|
||||||
%%
|
|
||||||
dg1_dPg=sparse(1:sizePGi,1:sizePGi,ones(sizePGi,1),sizePGi,sizePGi);
|
|
||||||
dg2_dPg=sparse(sizePGi,sizePVi);
|
|
||||||
dg3_dPg=sparse(sizePGi,sizeLoadi);
|
|
||||||
dg4_dPg=sparse(sizePGi,length(Loadi));
|
|
||||||
dg5_dPg=sparse(sizePGi,Busnum);
|
|
||||||
dg6_dPg=sparse(sizePGi,length(Loadi));
|
|
||||||
%%
|
|
||||||
dg1_dQr=sparse(sizePVi,sizePGi);
|
|
||||||
dg2_dQr=sparse(1:sizePVi,1:sizePVi,ones(sizePVi,1),sizePVi,sizePVi);
|
|
||||||
dg3_dQr=sparse(sizePVi,sizeLoadi);
|
|
||||||
dg4_dQr=sparse(sizePVi,length(Loadi));
|
|
||||||
dg5_dQr=sparse(sizePVi,Busnum);
|
|
||||||
dg6_dQr=sparse(sizePVi,length(Loadi));
|
|
||||||
%%
|
|
||||||
dg1_dPD=sparse(size(Loadi,1),size(PGi,1));
|
|
||||||
dg2_dPD=sparse(size(Loadi,1),size(PVi,1));
|
|
||||||
dg3_dPD=sparse(1:size(Loadi,1),1:size(Loadi,1),ones(size(Loadi,1),1),size(Loadi,1),size(Loadi,1));
|
|
||||||
dg4_dPD=sparse(size(Loadi,1),length(Loadi));
|
|
||||||
dg5_dPD=sparse(size(Loadi,1),Busnum);
|
|
||||||
dg6_dPD=diag(2./PD_-2*PD_./(PD_.^2+QD_.^2));
|
|
||||||
%%
|
|
||||||
dg1_dQD=sparse(size(Loadi,1),size(PGi,1));
|
|
||||||
dg2_dQD=sparse(size(Loadi,1),size(PVi,1));
|
|
||||||
dg3_dQD=sparse(length(Loadi),length(Loadi));
|
|
||||||
dg4_dQD=sparse(1:size(Loadi,1),1:size(Loadi,1),ones(size(Loadi,1),1),size(Loadi,1),size(Loadi,1));
|
|
||||||
dg5_dQD=sparse(size(Loadi,1),Busnum);
|
|
||||||
dg6_dQD=diag(-2*QD_./(PD_.^2+QD_.^2));
|
|
||||||
%%
|
|
||||||
dg1_dx=sparse(2*Busnum,sizePGi);
|
|
||||||
dg2_dx=sparse(2*Busnum,sizePVi);
|
|
||||||
dg3_dx=sparse(2*Busnum,sizeLoadi);
|
|
||||||
dg4_dx=sparse(2*Busnum,length(Loadi));
|
|
||||||
dg5_dx=[sparse(1:Busnum,1:Busnum,ones(Busnum,1),Busnum,Busnum);
|
|
||||||
sparse(Busnum,Busnum);
|
|
||||||
];
|
|
||||||
dg6_dx=sparse(2*Busnum,Busnum);
|
|
||||||
%%
|
|
||||||
deltG=[dg1_dPg,dg2_dPg,dg3_dPg,dg4_dPg,dg5_dPg,dg6_dPg;
|
|
||||||
dg1_dQr,dg2_dQr,dg3_dQr,dg4_dQr,dg5_dQr,dg6_dQr;
|
|
||||||
dg1_dPD,dg2_dPD,dg3_dPD,dg4_dPD,dg5_dPD,dg6_dPD;
|
|
||||||
dg1_dQD,dg2_dQD,dg3_dQD,dg4_dQD,dg5_dQD,dg6_dQD;
|
|
||||||
dg1_dx,dg2_dx,dg3_dx,dg4_dx,dg5_dx,dg6_dx;
|
|
||||||
];
|
|
||||||
end
|
|
||||||
58
func_deltG.m
58
func_deltG.m
|
|
@ -1,58 +0,0 @@
|
||||||
function deltG=func_deltG(Busnum,PVi,PGi,Loadi,PD,QD)
|
|
||||||
%%
|
|
||||||
PD_=PD(Loadi);
|
|
||||||
QD_=QD(Loadi);
|
|
||||||
sizePGi=size(PGi,1);
|
|
||||||
sizePVi=size(PVi,1);
|
|
||||||
sizeLoadi=size(Loadi,1);
|
|
||||||
%%
|
|
||||||
dg1_dPg=sparse(1:sizePGi,1:sizePGi,ones(sizePGi,1),sizePGi,sizePGi);
|
|
||||||
dg2_dPg=sparse(sizePGi,sizePVi);
|
|
||||||
dg3_dPg=sparse(sizePGi,sizeLoadi);
|
|
||||||
%dg4_dPg=sparse(sizePGi,length(Loadi));
|
|
||||||
dg5_dPg=sparse(sizePGi,Busnum);
|
|
||||||
dg6_dPg=sparse(sizePGi,length(Loadi));
|
|
||||||
%%
|
|
||||||
dg1_dQr=sparse(sizePVi,sizePGi);
|
|
||||||
dg2_dQr=sparse(1:sizePVi,1:sizePVi,ones(sizePVi,1),sizePVi,sizePVi);
|
|
||||||
dg3_dQr=sparse(sizePVi,sizeLoadi);
|
|
||||||
%dg4_dQr=sparse(sizePVi,length(Loadi));
|
|
||||||
dg5_dQr=sparse(sizePVi,Busnum);
|
|
||||||
dg6_dQr=sparse(sizePVi,length(Loadi));
|
|
||||||
%%
|
|
||||||
dg1_dPD=sparse(size(Loadi,1),size(PGi,1));
|
|
||||||
dg2_dPD=sparse(size(Loadi,1),size(PVi,1));
|
|
||||||
dg3_dPD=sparse(1:size(Loadi,1),1:size(Loadi,1),ones(size(Loadi,1),1),size(Loadi,1),size(Loadi,1));
|
|
||||||
%dg4_dPD=sparse(size(Loadi,1),length(Loadi));
|
|
||||||
dg5_dPD=sparse(size(Loadi,1),Busnum);
|
|
||||||
dg6_dPD=diag(2./PD_-2*PD_./(PD_.^2+QD_.^2));
|
|
||||||
%%
|
|
||||||
dg1_dQD=sparse(size(Loadi,1),size(PGi,1));
|
|
||||||
dg2_dQD=sparse(size(Loadi,1),size(PVi,1));
|
|
||||||
dg3_dQD=sparse(length(Loadi),length(Loadi));
|
|
||||||
%dg4_dQD=sparse(1:size(Loadi,1),1:size(Loadi,1),ones(size(Loadi,1),1),size(Loadi,1),size(Loadi,1));
|
|
||||||
dg5_dQD=sparse(size(Loadi,1),Busnum);
|
|
||||||
dg6_dQD=diag(-2*QD_./(PD_.^2+QD_.^2));
|
|
||||||
%%
|
|
||||||
dg1_dx=sparse(2*Busnum,sizePGi);
|
|
||||||
dg2_dx=sparse(2*Busnum,sizePVi);
|
|
||||||
dg3_dx=sparse(2*Busnum,sizeLoadi);
|
|
||||||
%dg4_dx=sparse(2*Busnum,length(Loadi));
|
|
||||||
dg5_dx=[sparse(1:Busnum,1:Busnum,ones(Busnum,1),Busnum,Busnum);
|
|
||||||
sparse(Busnum,Busnum);
|
|
||||||
];
|
|
||||||
dg6_dx=sparse(2*Busnum,length(Loadi));
|
|
||||||
%%
|
|
||||||
% deltG=[dg1_dPg,dg2_dPg,dg3_dPg,dg4_dPg,dg5_dPg,dg6_dPg;
|
|
||||||
% dg1_dQr,dg2_dQr,dg3_dQr,dg4_dQr,dg5_dQr,dg6_dQr;
|
|
||||||
% dg1_dPD,dg2_dPD,dg3_dPD,dg4_dPD,dg5_dPD,dg6_dPD;
|
|
||||||
% dg1_dQD,dg2_dQD,dg3_dQD,dg4_dQD,dg5_dQD,dg6_dQD;
|
|
||||||
% dg1_dx,dg2_dx,dg3_dx,dg4_dx,dg5_dx,dg6_dx;
|
|
||||||
% ];
|
|
||||||
deltG=[dg1_dPg,dg2_dPg,dg3_dPg,dg5_dPg,dg6_dPg;
|
|
||||||
dg1_dQr,dg2_dQr,dg3_dQr,dg5_dQr,dg6_dQr;
|
|
||||||
dg1_dPD,dg2_dPD,dg3_dPD,dg5_dPD,dg6_dPD;
|
|
||||||
dg1_dQD,dg2_dQD,dg3_dQD,dg5_dQD,dg6_dQD;
|
|
||||||
dg1_dx,dg2_dx,dg3_dx,dg5_dx,dg6_dx;
|
|
||||||
];
|
|
||||||
end
|
|
||||||
|
|
@ -1,8 +0,0 @@
|
||||||
function deltH=func_deltH(Busnum,Volt,PVi,Y,PGi,UAngel,r,c,Angle,Loadi)
|
|
||||||
|
|
||||||
dH_dPg=sparse(1:size(PGi,1),PGi,ones(size(PGi,1),1),size(PGi,1),2*Busnum);
|
|
||||||
dH_dQr=sparse(1:size(PVi,1),PVi+Busnum,ones(size(PVi,1),1),size(PVi,1),2*Busnum);
|
|
||||||
dH_dPD=[sparse(1:size(Loadi,1),Loadi,-ones(size(Loadi,1),1),size(Loadi,1),Busnum) sparse(size(Loadi,1),Busnum)];
|
|
||||||
dH_dx = jacobian_M3(Busnum,Volt,Y,Angle,UAngel,r,c); %形成雅克比矩阵
|
|
||||||
deltH=[dH_dPg;dH_dQr;dH_dPD;dH_dx'];
|
|
||||||
end
|
|
||||||
10
func_deltH.m
10
func_deltH.m
|
|
@ -1,10 +0,0 @@
|
||||||
function deltH=func_deltH(Busnum,Volt,PVi,Y,PGi,UAngel,r,c,Angle,Loadi)
|
|
||||||
|
|
||||||
dH_dPg=sparse(1:size(PGi,1),PGi,ones(size(PGi,1),1),size(PGi,1),2*Busnum);
|
|
||||||
dH_dQr=sparse(1:size(PVi,1),PVi+Busnum,ones(size(PVi,1),1),size(PVi,1),2*Busnum);
|
|
||||||
dH_dPD=[sparse(1:size(Loadi,1),Loadi,-ones(size(Loadi,1),1),size(Loadi,1),Busnum) sparse(size(Loadi,1),Busnum)];
|
|
||||||
dH_dQD=[sparse(size(Loadi,1),Busnum) sparse(1:size(Loadi,1),Loadi,-ones(size(Loadi,1),1),size(Loadi,1),Busnum)];
|
|
||||||
dH_dx = jacobian_M(Busnum,Volt,Y,Angle,UAngel,r,c); %ÐγÉÑſ˱ȾØÕó
|
|
||||||
deltH=[dH_dPg;dH_dQr;dH_dPD;dH_dQD;dH_dx'];
|
|
||||||
|
|
||||||
end
|
|
||||||
|
|
@ -1,11 +0,0 @@
|
||||||
function deltdeltF=func_deltdeltF(Busnum,PVi,PGi,wG,wD,ContrlCount)
|
|
||||||
|
|
||||||
ContrlCount=size(PVi,1)+size(PGi,1)+size(wD,1)+Busnum*2; %P,Q,Volt theta这些控制变量数
|
|
||||||
C=[wG' zeros(size(PVi))' wD'];
|
|
||||||
sizeC=size(C,2);
|
|
||||||
diagC=sparse(1:sizeC,1:sizeC,C,sizeC,sizeC);
|
|
||||||
deltdeltF=[
|
|
||||||
diagC*2,sparse(sizeC,ContrlCount-sizeC);
|
|
||||||
sparse(ContrlCount-sizeC,ContrlCount);
|
|
||||||
];
|
|
||||||
end
|
|
||||||
|
|
@ -1,11 +0,0 @@
|
||||||
function deltdeltF=func_deltdeltF(PVi,wPG,wQG,wPD,wQD,ContrlCount)
|
|
||||||
|
|
||||||
%ContrlCount=size(PVi,1)+size(PGi,1)+size(wD,1)+Busnum*2; %P,Q,Volt theta这些控制变量数
|
|
||||||
C=[wPG' wQG' wPD' wQD'];
|
|
||||||
sizeC=size(C,2);
|
|
||||||
diagC=sparse(1:sizeC,1:sizeC,C,sizeC,sizeC);
|
|
||||||
deltdeltF=[
|
|
||||||
diagC*2,sparse(sizeC,ContrlCount-sizeC);
|
|
||||||
sparse(ContrlCount-sizeC,ContrlCount);
|
|
||||||
];
|
|
||||||
end
|
|
||||||
20
jacobian_M.m
20
jacobian_M.m
|
|
@ -1,20 +0,0 @@
|
||||||
function [Jacob]=jacobian_M3(Busnum,Volt,Y,Angle,UAngel,r,c)
|
|
||||||
%**************************************************************************
|
|
||||||
% 程序功能 : 子函数——形成雅可比矩阵Jacobian
|
|
||||||
% 编 者:
|
|
||||||
% 编制时间:2010.12
|
|
||||||
%**************************************************************************
|
|
||||||
%%参照图书馆6楼的书编写
|
|
||||||
%% 分别求雅克比矩阵的子阵H,L,N,J及有功无功分量P,Q
|
|
||||||
AngleIJ=UAngel(r)-UAngel(c)-Angle';
|
|
||||||
mat_AngleIJ=sparse(r,c,AngleIJ,Busnum,Busnum);
|
|
||||||
mat_IvAngleIJ=mat_AngleIJ';
|
|
||||||
H=sparse(1:Busnum,1:Busnum,Volt,Busnum,Busnum)*Y.*sin(mat_IvAngleIJ)*sparse(1:Busnum,1:Busnum,Volt,Busnum,Busnum)-sparse(1:Busnum,1:Busnum,Y.*sin(mat_AngleIJ)*Volt',Busnum,Busnum)*sparse(1:Busnum,1:Busnum,Volt,Busnum,Busnum);
|
|
||||||
N=-sparse(1:Busnum,1:Busnum,Volt,Busnum,Busnum)*Y.*cos(mat_IvAngleIJ)*sparse(1:Busnum,1:Busnum,Volt,Busnum,Busnum)+sparse(1:Busnum,1:Busnum,Y.*cos(mat_AngleIJ)*Volt',Busnum,Busnum)*sparse(1:Busnum,1:Busnum,Volt,Busnum,Busnum);
|
|
||||||
J=sparse(1:Busnum,1:Busnum,Y.*cos(mat_AngleIJ)*Volt',Busnum,Busnum)+Y.*cos(mat_IvAngleIJ)*sparse(1:Busnum,1:Busnum,Volt,Busnum,Busnum);
|
|
||||||
L=sparse(1:Busnum,1:Busnum,Y.*sin(mat_AngleIJ)*Volt',Busnum,Busnum)+Y.*sin(mat_IvAngleIJ)*sparse(1:Busnum,1:Busnum,Volt,Busnum,Busnum);
|
|
||||||
t1=[J,L;
|
|
||||||
H,N;
|
|
||||||
]';
|
|
||||||
Jacob=-t1;
|
|
||||||
end
|
|
||||||
|
|
@ -1,18 +0,0 @@
|
||||||
function [new_G,new_B,GB,Y,r,c,Angle] = modifyadmmatrix(ii,jj,G,GG,B,BB)
|
|
||||||
%**************************************************************************
|
|
||||||
% 程序功能 : 子函数——形成节点导纳矩阵Y
|
|
||||||
% 编 者:
|
|
||||||
% 编制时间:2010.12
|
|
||||||
%**************************************************************************
|
|
||||||
%% 支路导纳计算
|
|
||||||
new_G=G;
|
|
||||||
new_G(ii,jj)=new_G(ii,jj)-G(ii,jj);
|
|
||||||
new_G(jj,ii)=new_G(jj,ii)-G(jj,ii);
|
|
||||||
new_G(ii,ii)=new_G(ii,ii)+G(ii,jj);
|
|
||||||
new_G(jj,jj)=new_G(jj,jj)+G(ii,jj);
|
|
||||||
|
|
||||||
%% 化作极坐标形式
|
|
||||||
GB = G+B.*1i; %将电导,电纳合并,写成复数形式
|
|
||||||
Y = abs(GB); %求节点导纳幅值
|
|
||||||
[r,c] = find(Y);
|
|
||||||
Angle = angle(GB(GB~=0)); %求节点导纳角度
|
|
||||||
|
|
@ -1,23 +0,0 @@
|
||||||
function [new_G,new_B,GB,Y,r,c,Angle] = modifyadmmatrix(ii,jj,G,B)
|
|
||||||
%**************************************************************************
|
|
||||||
% 程序功能 : 子函数——形成节点导纳矩阵Y
|
|
||||||
% 编 者:
|
|
||||||
% 编制时间:2010.12
|
|
||||||
%**************************************************************************
|
|
||||||
%% 支路导纳计算
|
|
||||||
new_G=G;
|
|
||||||
new_G(ii,jj)=new_G(ii,jj)-G(ii,jj);
|
|
||||||
new_G(jj,ii)=new_G(jj,ii)-G(jj,ii);
|
|
||||||
new_G(ii,ii)=new_G(ii,ii)+G(ii,jj);
|
|
||||||
new_G(jj,jj)=new_G(jj,jj)+G(ii,jj);
|
|
||||||
new_B=B;
|
|
||||||
new_B(ii,jj)=new_B(ii,jj)-B(ii,jj);
|
|
||||||
new_B(jj,ii)=new_B(jj,ii)-B(jj,ii);
|
|
||||||
new_B(ii,ii)=new_B(ii,ii)+B(ii,jj);
|
|
||||||
new_B(jj,jj)=new_B(jj,jj)+B(ii,jj);
|
|
||||||
|
|
||||||
%% 化作极坐标形式
|
|
||||||
GB = new_G+new_B.*1i; %将电导,电纳合并,写成复数形式
|
|
||||||
Y = abs(GB); %求节点导纳幅值
|
|
||||||
[r,c] = find(Y);
|
|
||||||
Angle = angle(GB(GB~=0)); %求节点导纳角度
|
|
||||||
34
pf.asv
34
pf.asv
|
|
@ -1,34 +0,0 @@
|
||||||
function [kmax,Precision,Uangle,U,Busnum,PVi,PVu,Balance,Y,Angle,P0,Q0,r,c,GB,Linei,Linej,Transfori,Transforj,GenU,GenL,GenC,PG,QG,PD,QD,CenterA,PGi,PVQU,PVQL]=pf(FileName)
|
|
||||||
%**************************************************************************
|
|
||||||
% 程序名称:电力系统潮流计算程序
|
|
||||||
% 程序算法:极坐标下的牛顿-拉夫逊法
|
|
||||||
% 程序功能:主函数
|
|
||||||
% 程序编者:
|
|
||||||
% 编制时间:2010.12
|
|
||||||
%**************************************************************************
|
|
||||||
clc;
|
|
||||||
tic;
|
|
||||||
%% 读取数据文件
|
|
||||||
[Busnum,Balance,PQstandard,Precision,Linei,Linej,Liner,Linex,Lineb,kmax,Transfori ,...
|
|
||||||
Transforj,Transforr,Transforx,Transfork0,Branchi,Branchb,Pointpoweri,PG,QG,PD,QD,PVi,PVu,GenU,GenL,GenC,CenterA,PGi,PVQU,PVQL]= openfile(FileName);
|
|
||||||
%% 形成节点导纳矩阵
|
|
||||||
[GB,Y,r,c,Angle] = admmatrix(Busnum,Linei,Linej,Liner,Linex,Lineb,Transfori,Transforj,Transforr,...
|
|
||||||
Transforx,Transfork0,Branchi,Branchb);
|
|
||||||
[P0,Q0,U,Uangle] = Initial(PG,PD,PQstandard,Pointpoweri,QG,QD,Busnum); %求功率不平衡量
|
|
||||||
disp('迭代次数i 最大不平衡量');
|
|
||||||
%% 循环体计算
|
|
||||||
for i = 0:kmax
|
|
||||||
[Jacob,PQ,U,Uangle] = jacobian(Busnum,Balance,PVi,PVu,U,Uangle,Y,Angle,P0,Q0,r,c); %形成雅克比矩阵
|
|
||||||
% disp('第一次雅克比');
|
|
||||||
m = max(abs(PQ));
|
|
||||||
m=full(m);
|
|
||||||
fprintf(' %u %.8f \n',i,m);
|
|
||||||
if m > Precision %判断不平衡量是否满足精度要求
|
|
||||||
[Uangle,U] = solvefun(Busnum,Jacob,PQ,Uangle,U); %求解修正方程,更新电压变量
|
|
||||||
else
|
|
||||||
disp(['收敛,迭代次数为',num2str(i),'次']);
|
|
||||||
break %若满足精度要求,则计算收敛
|
|
||||||
end
|
|
||||||
end
|
|
||||||
toc;
|
|
||||||
end
|
|
||||||
97
test.m
97
test.m
|
|
@ -1,97 +0,0 @@
|
||||||
%% Test
|
|
||||||
clear
|
|
||||||
arraymaxPDError=zeros(54,1);
|
|
||||||
for HH=0:53
|
|
||||||
|
|
||||||
arraymaxPDError(HH+1)=-100;
|
|
||||||
for CC=1:250
|
|
||||||
tic
|
|
||||||
[kmax,Precision,UAngel,Volt,Busnum,PVi,PVu,Balance,Y,Angle,P0,Q0,r,c,GB,Linei,Linej,Transfori,Transforj,GenU,GenL,GenC,PG,QG,PD,QD,CenterA,PGi,PVQU,PVQL]= ...
|
|
||||||
pf('D:\Project\最小化潮流\最小潮流算例\仙海919.txt');
|
|
||||||
%pf('c:/file31.txt');
|
|
||||||
%pf('ieee10471PG.dat');
|
|
||||||
|
|
||||||
%PVi电压节点序号
|
|
||||||
%PVu电压节点电压标幺值
|
|
||||||
Volt;
|
|
||||||
UAngel*180/3.1415926;
|
|
||||||
%Precision=Precision/1000000;
|
|
||||||
%Precision=Precision*1000;
|
|
||||||
%% 通过潮流计算PG
|
|
||||||
AngleIJ=sparse(r,c,UAngel(r)-UAngel(c)-Angle',Busnum,Busnum);
|
|
||||||
PGBal=diag(Volt)*Y.*cos(AngleIJ)*Volt';
|
|
||||||
|
|
||||||
%% 初值
|
|
||||||
PG0=PG;
|
|
||||||
PD0=PD;
|
|
||||||
%%
|
|
||||||
PG0(Balance)=PGBal(Balance);
|
|
||||||
%%
|
|
||||||
[Volt,UAngel,Init_Z,Init_W,Init_L,Init_U,Init_Y,PG,QG,RestraintCount,wG,wD,PD,PD0,randPDind,Loadi]=OPF_Init(Busnum,Balance,PG,QG,Volt,GenU,GenL,PVi,PGi,PVQU,PVQL,PD0,QD,HH);
|
|
||||||
Gap=(Init_L*Init_Z'-Init_U*Init_W');
|
|
||||||
KK=0;
|
|
||||||
plotGap=zeros(1,50);
|
|
||||||
ContrlCount=size(PVi,1)+size(PGi,1)+size(Loadi,1)+Busnum*2;
|
|
||||||
kmax=600;
|
|
||||||
%% 20120523 临时
|
|
||||||
QD_NON_ZERO=QD(PD==0 & QD~=0);
|
|
||||||
QD_NON_ZERO_IND=find(PD==0 & QD~=0);
|
|
||||||
%%
|
|
||||||
while(abs(Gap)>Precision)
|
|
||||||
if KK>kmax
|
|
||||||
break;
|
|
||||||
end
|
|
||||||
plotGap(KK+1)=Gap;
|
|
||||||
Init_u=Gap/2/RestraintCount*CenterA;
|
|
||||||
%AngleIJMat=0;
|
|
||||||
%% 开始计算OPF
|
|
||||||
%% 形成等式约束的雅克比
|
|
||||||
deltH=func_deltH(Busnum,Volt,PVi,Y,PGi,UAngel,r,c,Angle,Loadi);
|
|
||||||
%% 形成不等式约束的雅克比
|
|
||||||
deltG=func_deltG(Busnum,PVi,PGi,Loadi);
|
|
||||||
%%
|
|
||||||
L_1Z=diag(Init_Z./Init_L);
|
|
||||||
U_1W=diag(Init_W./Init_U);
|
|
||||||
%% 形成海森阵
|
|
||||||
deltdeltF=func_deltdeltF(PVi,wG,wD,ContrlCount);
|
|
||||||
%% 形成ddHy
|
|
||||||
ddh=func_ddh(Volt,Init_Y,Busnum,PVi,PGi,Y,UAngel,r,c,Angle,Loadi,ContrlCount);
|
|
||||||
%% 开始构建ddg
|
|
||||||
ddg=func_ddg(PGi,PVi,Busnum,RestraintCount,Loadi);
|
|
||||||
%% 开始构建deltF
|
|
||||||
deltF=func_deltF(PG,PVi,PGi,wG,wD,PG0,PD0,PD,Busnum,Loadi);
|
|
||||||
|
|
||||||
%% 形成方程矩阵
|
|
||||||
Luu=Init_U'.*Init_W'+Init_u*ones(RestraintCount,1);
|
|
||||||
Lul=Init_L'.*Init_Z'-Init_u*ones(RestraintCount,1);
|
|
||||||
Mat_G=FormG(Volt,PVi,PGi,PG,QG,PD,Loadi);
|
|
||||||
Mat_H=FormH(Busnum,Volt,PG,PD,QG,QD,Y,UAngel,r,c,Angle,QD_NON_ZERO,QD_NON_ZERO_IND);
|
|
||||||
Ly=Mat_H;
|
|
||||||
Lz=FormLz(Mat_G,Init_L,GenL,Busnum,PVQL,PD0,Loadi,KK);
|
|
||||||
Lw=FormLw(Mat_G,Init_U,GenU,Busnum,PVQU,PD0,Loadi,KK);
|
|
||||||
Lx=FormLx(deltF,deltH,Init_Y,deltG,Init_Z,Init_W);
|
|
||||||
%%YY=FormYY(Lul,Lz,Ly,Luu,Lw,Lx);
|
|
||||||
%% 开始解方程
|
|
||||||
XX=SolveIt(deltF,deltG,Init_L,Init_Z,Init_U,Init_W,deltdeltF,ddh,ddg,deltH,Init_Y,Ly,Lz,ContrlCount,Lw,Lul,Luu,RestraintCount,Lx,Balance,PVi,PGi,Busnum,Loadi);
|
|
||||||
%%取各分量
|
|
||||||
[deltZ,deltL,deltW,deltU,deltX,deltY]=AssignXX(XX,ContrlCount,RestraintCount,Busnum);
|
|
||||||
[Init_Z,Init_L,Init_W,Init_U,Init_Y,PG,QG,Volt,UAngel,PD]=Modification(Init_Z,Init_L,Init_W,Init_U,Init_Y,deltZ,deltL,deltW,deltU,deltX,deltY,PG,QG,Volt,UAngel,PVi,ContrlCount,Balance,Busnum,PGi,PD,Loadi);
|
|
||||||
Gap=(Init_L*Init_Z'-Init_U*Init_W');
|
|
||||||
KK=KK+1;
|
|
||||||
end
|
|
||||||
fprintf('迭代次数%d\n',KK);
|
|
||||||
ObjectiveFun(PG,PG0,PGi,PD,PD0,wG,wD,Loadi);
|
|
||||||
%DrawGap(plotGap);
|
|
||||||
%%
|
|
||||||
%Volt=full(Volt');
|
|
||||||
%PD=full(PD);
|
|
||||||
%% 统计PD误差
|
|
||||||
maxPDError=max(abs( (PD(Loadi)-PD0(Loadi))./PD0(Loadi) ));
|
|
||||||
if maxPDError>arraymaxPDError(HH+1)
|
|
||||||
arraymaxPDError(HH+1)=maxPDError;
|
|
||||||
end
|
|
||||||
toc
|
|
||||||
end
|
|
||||||
end
|
|
||||||
|
|
||||||
|
|
||||||
Loading…
Reference in New Issue