From 1723ffec75c86f6498fbd1c03de87d962e383203 Mon Sep 17 00:00:00 2001 From: "facat@lab.com" Date: Mon, 3 Nov 2014 20:15:57 +0800 Subject: [PATCH] =?UTF-8?q?1.=E6=8A=8A=E6=95=B0=E6=8D=AE=E6=96=87=E4=BB=B6?= =?UTF-8?q?=E7=9A=84=E7=94=B5=E5=AE=B9=E5=88=A0=E6=8E=89=E4=BA=86=E3=80=82?= =?UTF-8?q?=202.=E8=B4=9F=E8=8D=B7=E4=BC=B0=E8=AE=A1=E5=8F=98=E9=87=8F?= =?UTF-8?q?=E7=9A=84=E7=BB=B4=E6=95=B0=E4=B8=8E=E8=8A=82=E7=82=B9=E6=95=B0?= =?UTF-8?q?=E4=B8=8D=E4=B8=80=E8=87=B4=EF=BC=8C=E5=8F=AA=E7=AD=89=E4=BA=8E?= =?UTF-8?q?=E8=B4=9F=E8=8D=B7=E6=95=B0=E3=80=82=203.=E5=8A=A0=E4=BA=86?= =?UTF-8?q?=E6=96=B9=E5=B7=AE=E3=80=82?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Signed-off-by: facat@lab.com --- FormH.m | 6 +++--- Modification.m | 4 ++-- OPF.m | 51 +++++++++++++++++++++++++++++++++++----------- feeder13/data1.txt | 2 +- func_deltdeltF.m | 2 +- 5 files changed, 46 insertions(+), 19 deletions(-) diff --git a/FormH.m b/FormH.m index 1d558d4..b0a55a8 100644 --- a/FormH.m +++ b/FormH.m @@ -1,4 +1,4 @@ -function Mat_H=FormH(Busnum,Volt,PG,PD,QG,QD,Y,UAngel,r,c,Angle) +function Mat_H=FormH(Busnum,Volt,PG,PD,QG,QD,Y,UAngel,r,c,Angle,Loadi) %% %QDcos=textread('D:\Project\最小化潮流\最小潮流算例\仙海919PDQDglys.txt'); %QD(QD~=0)=PD(QD~=0)./tan(QDcos); @@ -6,8 +6,8 @@ function Mat_H=FormH(Busnum,Volt,PG,PD,QG,QD,Y,UAngel,r,c,Angle) %% %PD(Loadi)=QD(Loadi)./tan(acos(0.98)); AngleIJ=sparse(r,c,UAngel(r)-UAngel(c)-Angle,Busnum*3,Busnum*3); -dP=PG-PD-diag(Volt)*Y.*cos(AngleIJ)*Volt; -dQ=QG-QD-diag(Volt)*Y.*sin(AngleIJ)*Volt; +dP=PG-sparse(Loadi,1,PD,Busnum*3,1)-diag(Volt)*Y.*cos(AngleIJ)*Volt; +dQ=QG-sparse(Loadi,1,QD,Busnum*3,1)-diag(Volt)*Y.*sin(AngleIJ)*Volt; Mat_H=[dP;dQ;]; diff --git a/Modification.m b/Modification.m index 646b1cc..8eabcf8 100644 --- a/Modification.m +++ b/Modification.m @@ -14,8 +14,8 @@ Init_Y=Init_Y+AlphaD*deltY'; %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(1:size(Loadi,1)*2); -PD(Loadi)=PD(Loadi)+AlphaP*t(1:length(Loadi)); -QD(Loadi)=QD(Loadi)+AlphaP*t(length(Loadi)+1:length(Loadi)*2); +PD=PD+AlphaP*t(1:length(Loadi)); +QD=QD+AlphaP*t(length(Loadi)+1:length(Loadi)*2); t=deltX(size(Loadi,1)*2+1:ContrlCount)'; t(Busnum*3+(Balance-1)*3+1)=0; t(Busnum*3+(Balance-1)*3+2)=0; diff --git a/OPF.m b/OPF.m index c63fc32..00524cc 100644 --- a/OPF.m +++ b/OPF.m @@ -123,7 +123,6 @@ ub=checkSSatisfied(Balance,phaseABCY,VoltpABC, ... fprintf('最大不平衡量为%f\n\n',full(max(abs(ub)))) %% 潮流计算end - busNum=length(nodeNum); Busnum=busNum; PQi=setxor(nodeNum,Balance); @@ -139,14 +138,17 @@ PQi3P=zeros(length(PQi)*3,1); PQi3P(1:3:end)=(PQi-1)*3+1; PQi3P(2:3:end)=(PQi-1)*3+2; PQi3P(3:3:end)=(PQi-1)*3+3; -PD3P=sparse(busNum*3,1); -QD3P=sparse(busNum*3,1); +Loadi=PQi3P; +PD3P=sparse(Busnum*3,1); +QD3P=sparse(Busnum*3,1); PD3P(1:3:end)=phaseASpotLoadP; PD3P(2:3:end)=phaseBSpotLoadP; PD3P(3:3:end)=phaseCSpotLoadP; QD3P(1:3:end)=phaseASpotLoadQ; QD3P(2:3:end)=phaseBSpotLoadQ; QD3P(3:3:end)=phaseCSpotLoadQ; +PD3P=PD3P(Loadi); +QD3P=QD3P(Loadi); QGi3P=zeros(length(QGi)*3,1); QGi3P(1:3:end)=(QGi-1)*3+1; QGi3P(2:3:end)=(QGi-1)*3+2; @@ -168,13 +170,38 @@ Balance3P(3:3:end)=(Balance-1)*3+3; Vp3a((Balance-1)*3+1)=0; Vp3a((Balance-1)*3+2)=-120/180*pi; Vp3a((Balance-1)*3+3)=+120/180*pi; -Loadi=PQi3P; -wVolt=ones(busNum*3,1); -wPD=ones(length(Loadi),1); -wQD=ones(length(Loadi),1); -mVolt=Vp3m; -PD0=PD3P; -QD0=QD3P; +%准备量测量和方差 +%真实值 +%三相电压幅值 +rVoltABCV=zeros(busNum*3,1); +rVoltABCV(1:3:end)=abs(VoltpABC(1,:)); +rVoltABCV(2:3:end)=abs(VoltpABC(2,:)); +rVoltABCV(3:3:end)=abs(VoltpABC(3,:)); +%三相电压相角 +rVoltABCA=zeros(busNum*3,1); +rVoltABCA(1:3:end)=angle(VoltpABC(1,:)); +rVoltABCA(2:3:end)=angle(VoltpABC(2,:)); +rVoltABCA(3:3:end)=angle(VoltpABC(3,:)); +%三相负荷 +rPD3P=zeros(busNum*3,1); +rPD3P(1:3:end)=phaseASpotLoadP; +rPD3P(2:3:end)=phaseBSpotLoadP; +rPD3P(3:3:end)=phaseCSpotLoadP; +rQD3P=zeros(busNum*3,1); +rQD3P(1:3:end)=phaseASpotLoadQ; +rQD3P(2:3:end)=phaseBSpotLoadQ; +rQD3P(3:3:end)=phaseCSpotLoadQ; +rPD3P=rPD3P(Loadi); +rQD3P=rQD3P(Loadi); +%量测量 +sigma=0.03; +mVoltABCV=rVoltABCV.*(1+normrnd(0,sigma,length(rVoltABCV),1)); +mPD3P=rPD3P.*(1+normrnd(0,sigma,length(rPD3P),1)); +mQD3P=rPD3P.*(1+normrnd(0,sigma,length(rQD3P),1)); +%量测方差 +wVolt=1./(abs(mVoltABCV*sigma).^2); +wPD=1./(abs(mPD3P*sigma).^2); +wQD=1./(abs(mQD3P*sigma).^2); %% RestraintCount=size(Loadi,1)*2; %约束条件数,放开所有QD Init_Z=sparse(ones(1,RestraintCount)); @@ -215,13 +242,13 @@ while(abs(Gap)>Precision) % ddg=func_ddg(PGi,PVi,Busnum,RestraintCount,Loadi,PD,QD); ddg=0; %% 开始构建deltF - deltF=func_deltF(wVolt,wPD,wQD,PD0,PD3P,QD3P,QD0,Volt,mVolt,Busnum,Loadi); + deltF=func_deltF(wVolt,wPD,wQD,mPD3P,PD3P,QD3P,mQD3P,Volt,mVoltABCV,Busnum,Loadi); % deltF=0; %% 形成方程矩阵 Luu=Init_U'.*Init_W'+Init_u*ones(RestraintCount,1); Lul=Init_L'.*Init_Z'-Init_u*ones(RestraintCount,1); Mat_G=FormG(PD3P,QD3P,Loadi); - Mat_H=FormH(Busnum,Volt,PG3P,PD3P,QG3P,QD3P,Y,UAngel,r,c,Yangle); + Mat_H=FormH(Busnum,Volt,PG3P,PD3P,QG3P,QD3P,Y,UAngel,r,c,Yangle,Loadi); Ly=Mat_H; Lz=FormLz(Mat_G,Init_L,Loadi); Lw=FormLw(Mat_G,Init_U,Loadi); diff --git a/feeder13/data1.txt b/feeder13/data1.txt index a4d2f06..a873457 100644 --- a/feeder13/data1.txt +++ b/feeder13/data1.txt @@ -4,7 +4,7 @@ 602 632 631 1800 602 631 630 1700 0 -632 30 0 0 +632 0 0 0 0 630 80 70 60 50 40 30 631 40 10 11 50 30 30 diff --git a/func_deltdeltF.m b/func_deltdeltF.m index 77f6b9c..bae5e34 100644 --- a/func_deltdeltF.m +++ b/func_deltdeltF.m @@ -3,7 +3,7 @@ function deltdeltF=func_deltdeltF(wVolt,wPD,wQD,ContrlCount) %ContrlCount=size(PVi,1)+size(PGi,1)+size(wD,1)+Busnum*2; %P,Q,Volt theta这些控制变量数 C=[wPD' wQD']; CVolt=wVolt; -sizeC=size(C,2); +sizeC=length(C); diagC=sparse(1:sizeC,1:sizeC,C,sizeC,sizeC); deltdeltF=[ diagC*2,sparse(sizeC,ContrlCount-sizeC);