From 0bc157a3ba38ae738c50b3a6429bba492e91116d Mon Sep 17 00:00:00 2001
From: unknown <dmy@dmy-PC.(none)>
Date: Fri, 16 Nov 2012 22:10:18 +0800
Subject: [PATCH] =?UTF-8?q?=E8=BF=99=E6=98=AF=E4=B8=80=E4=B8=AA=E4=B8=8D?=
 =?UTF-8?q?=E6=8A=8A=E5=8A=9F=E7=8E=87=E5=9B=A0=E6=95=B0=E4=BD=9C=E4=B8=BA?=
 =?UTF-8?q?=E4=B8=8D=E7=AD=89=E5=BC=8F=E7=BA=A6=E6=9D=9F=EF=BC=8C=E8=BF=98?=
 =?UTF-8?q?=E6=98=AF=E5=9C=A8QD=E7=9A=84=E4=B8=8D=E7=AD=89=E5=BC=8F?=
 =?UTF-8?q?=E7=BA=A6=E6=9D=9F=E4=B8=AD=E6=9C=9F=E4=BD=9C=E7=94=A8=E7=9A=84?=
 =?UTF-8?q?=E7=89=88=E6=9C=AC=E3=80=82?=
MIME-Version: 1.0
Content-Type: text/plain; charset=UTF-8
Content-Transfer-Encoding: 8bit

Signed-off-by: unknown <dmy@dmy-PC.(none)>
---
 FormLw.m         |  8 +++++---
 FormLz.m         |  8 +++++---
 OPF.m            | 21 ++++++++++++---------
 OPF_Init.m       |  5 +++--
 ObjectiveFun.m   | 11 +++++++----
 func_deltF.m     |  4 ++--
 func_deltdeltF.m |  4 ++--
 pf.m             |  6 +++++-
 8 files changed, 41 insertions(+), 26 deletions(-)

diff --git a/FormLw.m b/FormLw.m
index 86e4298..d2778a2 100644
--- a/FormLw.m
+++ b/FormLw.m
@@ -1,8 +1,10 @@
 function Lw=FormLw(Mat_G,Init_U,GenU,Busnum,PVQU,PD0,QD0,Loadi,KK)
 KK=999;
-PU=GenU(:,2);%发电机有功上界
-QU=PVQU(:,1);%发电机无功上界
-VoltU=(1.2)*ones(1,Busnum);
+%PU=GenU(:,2);%发电机有功上界
+PU=5*ones(length(GenU(:,2)),1);
+%QU=PVQU(:,1);%发电机无功上界
+QU=5*ones(length(PVQU(:,1)),1);
+VoltU=(1.1)*ones(1,Busnum);
 %VoltU=10*ones(1,Busnum);
 PDU=PD0(Loadi);
 PDU(PDU>0)=1.200*PDU(PDU>0);
diff --git a/FormLz.m b/FormLz.m
index cda108f..f76380f 100644
--- a/FormLz.m
+++ b/FormLz.m
@@ -1,8 +1,10 @@
 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);
+%PL=GenL(:,2);%发电机有功下界
+PL=-5*ones(length(GenL(:,2)),1);
+%QL=PVQL(:,1);%发电机无功下界
+QL=-5*ones(length(PVQL(:,1)),1);
+VoltL=(0.9)*ones(1,Busnum);
 %VoltL=-10*ones(1,Busnum);
 PDL=PD0(Loadi);
 PDL(PDL>0)=0.800*PDL(PDL>0);
diff --git a/OPF.m b/OPF.m
index 7234876..6ba6041 100644
--- a/OPF.m
+++ b/OPF.m
@@ -2,7 +2,7 @@ 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('ieee10471.dat');
+pf('ieee118.dat');
 %pf('D:\Project\青秀降损项目\最小化潮流\最小潮流算例\原始\津头站津视922（3-1）_0.5_120%.txt'); 
 %pf('D:\Project\最小化潮流\最小潮流算例\仙海919.txt'); 
 %pf('c:/file31.txt'); 
@@ -14,8 +14,8 @@ 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';
+PGBal=PD+diag(Volt)*Y.*cos(AngleIJ)*Volt';
+QGBal=QD+diag(Volt)*Y.*sin(AngleIJ)*Volt';
 %% 初值-即测量值
 PG0=PG;
 QG0=QG;
@@ -25,9 +25,12 @@ PDReal=PD;%
 %PD0(12)=PD0(12)+0.001;
 %%
 PG0(Balance)=PGBal(Balance);
+PG(Balance)=PGBal(Balance);
 QG0(Balance)=QGBal(Balance);
+QG0(PVi)=QGBal(PVi);
+QG(PVi)=QGBal(PVi);
 %%
-[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);
+[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);
@@ -53,13 +56,13 @@ while(abs(Gap)>Precision)
     L_1Z=diag(Init_Z./Init_L);
     U_1W=diag(Init_W./Init_U);
     %% 形成海森阵
-    deltdeltF=func_deltdeltF(PVi,wPG,wQG,wD,ContrlCount);
+    deltdeltF=func_deltdeltF(PVi,wPG,wQG,wPD,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,PD,QD);
     %%  开始构建deltF
-    deltF=func_deltF(PG,QG,PVi,PGi,wPG,wQG,wD,PG0,QG0,PD0,PD,Busnum,Loadi);
+    deltF=func_deltF(PG,QG,PVi,PGi,wPG,wQG,wPD,PG0,QG0,PD0,PD,Busnum,Loadi);
     
     %% 形成方程矩阵
     Luu=Init_U'.*Init_W'+Init_u*ones(RestraintCount,1);
@@ -81,7 +84,7 @@ while(abs(Gap)>Precision)
     KK=KK+1;
 end
 fprintf('迭代次数%d\n',KK);
-%fprintf('目标值%f\n',full(ObjectiveFun(PG,PG0,PGi,PD,PD0,wG,wD,Loadi)));
+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');
@@ -89,13 +92,13 @@ DrawGap(plotGap);
 %% 统计PD误差
 % absPDLoad=abs(  (PD(Loadi)-PD0(Loadi))./PD0(Loadi)  );
 absPDLoad=abs(  (PD(Loadi)-PDReal(Loadi))./PDReal(Loadi)  );
-maxPDError=max(absPDLoad);
+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);
+%Lineloss(Linei,Linej,Liner,Linex,Lineb,Transfori,Transforj,Transforr,Transforx,Transfork0,Volt,UAngel);
 toc  
 
diff --git a/OPF_Init.m b/OPF_Init.m
index cdacdf9..25d3ccb 100644
--- a/OPF_Init.m
+++ b/OPF_Init.m
@@ -1,4 +1,4 @@
-function [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)
+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]=OPF_Init(Busnum,Balance,PG,QG,Volt,GenU,GenL,PVi,PGi,PVQU,PVQL,PD0,QD0,QD,PD)
 Loadi=find(QD~=0 | PD~=0);
 %Loadi=[1:Busnum]';
 RestraintCount=size(PVi,1)+size(PGi,1)+size(Loadi,1)*2+Busnum*1; %约束条件数,放开所有QD
@@ -22,7 +22,8 @@ wQG=1*ones(size(PVi,1),1);
 %randInt=randperm(size(Loadi,1));
 %randPDind=randInt(1:10);
 randPDind=0;
-wD=1*ones(size(Loadi,1),1);
+wPD=1*ones(size(Loadi,1),1);
+wQD=1*ones(size(Loadi,1),1);
 %wD(randPDind)=0;%一些负荷不约束
 %wD(7)=0;
 % wD(11)=0;
diff --git a/ObjectiveFun.m b/ObjectiveFun.m
index 967e682..590d571 100644
--- a/ObjectiveFun.m
+++ b/ObjectiveFun.m
@@ -1,6 +1,9 @@
-function [out_arg]=ObjectiveFun(PG,PG0,PGi,PD,PD0,wG,wD,Loadi)
+function [out_arg]=ObjectiveFun(PG,PG0,PGi,QG,QG0,PVi,PD,PD0,QD,QD0,wPG,wQG,wPD,wQD,Loadi)
 t1=PG(PGi)-PG0(PGi);
-t2=wG.*t1.*t1;
-t3=wD.*((PD(Loadi)-PD0(Loadi)).*(PD(Loadi)-PD0(Loadi)));
-out_arg= sum(t2)+sum(t3);
+t2=wPG.*t1.*t1;
+t1=QG(PVi)-QG0(PVi);
+t3=wQG.*t1.*t1;
+t4=wPD.*((PD(Loadi)-PD0(Loadi)).*(PD(Loadi)-PD0(Loadi)));
+t5=wQD.*((QD(Loadi)-QD0(Loadi)).*(QD(Loadi)-QD0(Loadi)));
+out_arg= sum(t2)+sum(t3)+sum(t4)+sum(t5);
 end
\ No newline at end of file
diff --git a/func_deltF.m b/func_deltF.m
index e579389..5a22b53 100644
--- a/func_deltF.m
+++ b/func_deltF.m
@@ -1,8 +1,8 @@
-function deltF=func_deltF(PG,QG,PVi,PGi,wPG,wQG,wD,PG0,QG0,PD0,PD,Busnum,Loadi)
+function deltF=func_deltF(PG,QG,PVi,PGi,wPG,wQG,wPD,PG0,QG0,PD0,PD,Busnum,Loadi)
 
 t1=2*wPG.*(PG(PGi) - PG0(PGi) );
 t2=2*wQG.*(QG(PVi) - QG0(PVi) );
-t3=2*wD.*(PD(Loadi)-PD0(Loadi));
+t3=2*wPD.*(PD(Loadi)-PD0(Loadi));
 deltF=[
     sparse(t1);
     sparse(t2);
diff --git a/func_deltdeltF.m b/func_deltdeltF.m
index 2f5f97c..d08f22b 100644
--- a/func_deltdeltF.m
+++ b/func_deltdeltF.m
@@ -1,7 +1,7 @@
-function deltdeltF=func_deltdeltF(PVi,wPG,wQG,wD,ContrlCount)
+function deltdeltF=func_deltdeltF(PVi,wPG,wQG,wPD,ContrlCount)
 
 %ContrlCount=size(PVi,1)+size(PGi,1)+size(wD,1)+Busnum*2; %P,Q,Volt theta这些控制变量数
-C=[wPG' wQG' wD'];
+C=[wPG' wQG' wPD'];
 sizeC=size(C,2);
 diagC=sparse(1:sizeC,1:sizeC,C,sizeC,sizeC);
 deltdeltF=[
diff --git a/pf.m b/pf.m
index b585f71..c3df35a 100644
--- a/pf.m
+++ b/pf.m
@@ -31,4 +31,8 @@ for i = 0:kmax
     end
 end
 toc;
-end
\ No newline at end of file
+PG=PG/PQstandard;
+QG=QG/PQstandard;
+PD=PD/PQstandard;
+QD=QD/PQstandard;
+end