将电压加入目标函数中。

OPF.m

@@ -42,7 +42,7 @@ for II=1:1
     RealQD=QD0;
     %%
     %PGi=zeros(1,1);
-    [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,RealPD,RealQD,QD,PD);
+    [Volt,UAngel,Init_Z,Init_W,Init_L,Init_U,Init_Y,PG,QG,RestraintCount,wPG,wQG,wPD,wQD,wVolt,PD,PD0,QD,randPDind,Loadi]=OPF_Init(Busnum,Balance,PG,QG,Volt,GenU,GenL,PVi,PGi,PVQU,PVQL,RealPD,RealQD,QD,PD);
     Gap=(Init_L*Init_Z'-Init_U*Init_W');
     KK=0;
     plotGap=zeros(1,60);
@@ -80,13 +80,13 @@ for II=1:1
         L_1Z=diag(Init_Z./Init_L);
         U_1W=diag(Init_W./Init_U);
         %% 形成海森阵
-        deltdeltF=func_deltdeltF(PVi,wPG,wQG,wPD,wQD,ContrlCount);
+        deltdeltF=func_deltdeltF(PVi,wPG,wQG,wVolt,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,RestraintCount,Loadi,PD,QD);
         %%  开始构建deltF
-        deltF=func_deltF(PG,QG,PVi,PGi,wPG,wQG,wPD,wQD,PG0,QG0,PD0,PD,QD,QD0,Busnum,Loadi);
+        deltF=func_deltF(Volt,PG,QG,PVi,PGi,wPG,wQG,wPD,wQD,wVolt,PG0,QG0,PD0,PD,QD,QD0,Volt0,Busnum,Loadi);
         %% 形成方程矩阵
         Luu=Init_U'.*Init_W'+Init_u*ones(RestraintCount,1);
         Lul=Init_L'.*Init_Z'-Init_u*ones(RestraintCount,1);

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,wPD,wQD,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,wVolt,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(Loadi,1)*2+Busnum*1; %约束条件数,放开所有QD
@@ -20,10 +20,10 @@ wQG=0;
 %randInt=randperm(size(Loadi,1));
 %randPDind=randInt(1:10);
 randPDind=0;
-wPD=1/.05^2*zeros(size(Loadi,1),1);
-wPD(1:2:end)=0;
-wQD=1/.05^2*zeros(size(Loadi,1),1);
-wQD(1:2:end)=0;
+wPD=1/.05^2*ones(size(Loadi,1),1);
+% wPD(1:2:end)=0;
+wQD=1/.05^2*ones(size(Loadi,1),1);
+% wQD(1:2:end)=0;
 %wD(randPDind)=0;%一些负荷不约束
 %wD(7)=0;
 % wD(11)=0;
@@ -31,4 +31,5 @@ PD=0.5*PD0;
 %powerFacter=0.98;
 %QD=PD*sqrt((1-powerFacter^2)/powerFacter^2);
 QD=.5*QD0;
+wVolt=1/.05^2*ones(size(Busnum,1),1);
 end

func_deltF.m

@@ -1,9 +1,11 @@
-function deltF=func_deltF(PG,QG,PVi,PGi,wPG,wQG,wPD,wQD,PG0,QG0,PD0,PD,QD,QD0,Busnum,Loadi)
+function deltF=func_deltF(Volt,PG,QG,PVi,PGi,wPG,wQG,wPD,wQD,wVolt,PG0,QG0,PD0,PD,QD,QD0,Volt0,Busnum,Loadi)
 t3=2*wPD.*(PD(Loadi)-PD0(Loadi));
 t4=2*wQD.*(QD(Loadi)-QD0(Loadi));
+t5=2*wVolt.*(Volt-Volt0);
 deltF=[sparse(t3);
     sparse(t4);
-    sparse(2*Busnum,1);
+    sparse(t5');
+    sparse(1*Busnum,1);
 ];
 
 end

func_deltdeltF.m

@@ -1,7 +1,7 @@
-function deltdeltF=func_deltdeltF(PVi,wPG,wQG,wPD,wQD,ContrlCount)
+function deltdeltF=func_deltdeltF(PVi,wPG,wQG,wVolt,wPD,wQD,ContrlCount)
 
 %ContrlCount=size(PVi,1)+size(PGi,1)+size(wD,1)+Busnum*2; %P,Q,Volt theta这些控制变量数
-C=[wPD' wQD'];
+C=[wPD' wQD' wVolt];
 sizeC=size(C,2);
 diagC=sparse(1:sizeC,1:sizeC,C,sizeC,sizeC);
 deltdeltF=[