From e8611477743efa83051f1a92d1a08689d1b1af05 Mon Sep 17 00:00:00 2001 From: "dugg@lab-desk" Date: Sat, 18 Oct 2014 21:24:17 +0800 Subject: [PATCH] =?UTF-8?q?=E6=8A=8A=E5=86=85=E7=82=B9=E6=B3=95=E5=BE=AA?= =?UTF-8?q?=E7=8E=AF=E5=8D=95=E7=8B=AC=E6=94=BE=E5=88=B0=E4=B8=80=E4=B8=AA?= =?UTF-8?q?=E5=87=BD=E6=95=B0=E6=96=87=E4=BB=B6=E4=B8=AD=E3=80=82?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Signed-off-by: dugg@lab-desk --- IPMLoop.m | 61 +++++++++++++++++++++++++++++++++++++++++++++++++++++++ run.m | 57 ++------------------------------------------------- 2 files changed, 63 insertions(+), 55 deletions(-) create mode 100644 IPMLoop.m diff --git a/IPMLoop.m b/IPMLoop.m new file mode 100644 index 0000000..c3654e5 --- /dev/null +++ b/IPMLoop.m @@ -0,0 +1,61 @@ +function [ V1r,V1i,I1r,I1i ] = IPMLoop(measurement,BalI1r,BalI1i,busNum,Loadi,fsY1,Balance ) +%把每个序的循环写在这个函数中。其实也就是内点法循环。 +V1r=1*ones(busNum,1); +V1i=0*ones(busNum,1); +I1r=0.01*ones(length(Loadi),1);%注入电流，相当于放电机电流，与负荷电流负荷相反 +I1i=0.01*ones(length(Loadi),1);%注入电流，相当于放电机电流，与负荷电流负荷相反 +KK=0; +plotGap=zeros(1,60); +%初始化 +%状态量为 SEPD SEQD SEVmf1 SEVaf1 +RestraintCount=length(Loadi)*1; +ContrlCount=busNum*2+length(Loadi)*2; +CenterA=0.1; +Init_Z=sparse(ones(RestraintCount,1)); +Init_W=sparse(-1*ones(RestraintCount,1)); +Init_L=1*sparse(ones(RestraintCount,1)); +Init_U=1*sparse(ones(RestraintCount,1)); +Init_Y=sparse(2*busNum,1);%等式约束乘子 +Gap=(Init_L'*Init_Z-Init_U'*Init_W); +% PG=sparse(busNum,1); +% PG(Balance)=0.1105; +% QG=sparse(busNum,1); +% QG(Balance)=0.0984; +% SEPD(2)=0.1105; +% SEQD(2)=0.0984; +while Gap>1e-5 && KK<20 + KK=KK+1; + Init_u=Gap/2/RestraintCount*CenterA; + deltH=func_deltH(busNum,fsY1,Loadi,Balance); + deltG=func_deltG(busNum,Loadi,I1r,I1i); +% deltG=0; + L_1Z=diag(Init_Z./Init_L); + U_1W=diag(Init_W./Init_U); + deltdeltF=func_deltdeltF(busNum,fsY1,Loadi); +% deltdeltF=0; +% ddh=func_ddh(SEVmf1,Init_Y,busNum,fsY1amp,SEVaf1,r,c,fsY1ang,Loadi,ContrlCount); + ddh=0; + ddg=func_ddg(); + deltF=func_deltF(measurement,busNum,fsY1,Loadi,V1r,V1i,I1r,I1i); +% 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(I1r,I1i); + Mat_H=FormH(fsY1,Loadi,V1r,V1i,I1r,I1i,BalI1r,BalI1i,Balance); + Ly=Mat_H; + Lz=FormLz(Loadi,Mat_G,Init_L); + Lw=FormLw(Loadi,Mat_G,Init_U); + Lx=FormLx(deltF,deltH,Init_Y,deltG,Init_Z,Init_W); + YY=FormYY(Lul,Lz,Ly,Luu,Lw,Lx); + %% 开始解方程 + plotGap(KK)=Gap; + fprintf('迭代次数 %d Gap %f\n',KK,plotGap(KK)); + XX=SolveIt(deltF,deltG,Init_L,Init_Z,Init_U,Init_W,deltdeltF,ddh,ddg,deltH,Init_Y,Ly,Lz,ContrlCount,Lw,Lul,Luu,Lx,Balance,busNum,Loadi); + [deltZ,deltL,deltW,deltU,deltX,deltY]=AssignXX(XX,ContrlCount,RestraintCount,busNum); + [Init_Z,Init_L,Init_W,Init_U,Init_Y,V1r,V1i,I1r,I1i]=Modification(Init_Z,Init_L,Init_W,Init_U,Init_Y,deltZ,deltL,deltW,deltU,deltX,deltY,V1r,V1i,I1r,I1i,ContrlCount,Balance,busNum,Loadi); + Gap=(Init_L'*Init_Z-Init_U'*Init_W); +end + + +end + diff --git a/run.m b/run.m index e7a657c..f85eb23 100644 --- a/run.m +++ b/run.m @@ -232,10 +232,7 @@ mIf2=If2; BalI1r=real(-sum(mIf1)); BalI1i=imag(-sum(mIf1)); inv(fsY11)*(mIf1); -V1r=1*ones(busNum,1); -V1i=0*ones(busNum,1); -I1r=0.01*ones(length(Loadi),1);%注入电流，相当于放电机电流，与负荷电流负荷相反 -I1i=0.01*ones(length(Loadi),1);%注入电流，相当于放电机电流，与负荷电流负荷相反 + measurement=-mIf1(Loadi); clear PD QD PG QG; %状态量 @@ -253,56 +250,6 @@ clear PD QD PG QG; % SEVaf1=sparse(zeros(busNum,1)); % SEPD=sparse(zeros(busNum,1)); % SEQD=sparse(zeros(busNum,1)); -KK=0; -plotGap=zeros(1,60); -%初始化 -%状态量为 SEPD SEQD SEVmf1 SEVaf1 -RestraintCount=length(Loadi)*1; -ContrlCount=busNum*2+length(Loadi)*2; -CenterA=0.1; -Init_Z=sparse(ones(RestraintCount,1)); -Init_W=sparse(-1*ones(RestraintCount,1)); -Init_L=1*sparse(ones(RestraintCount,1)); -Init_U=1*sparse(ones(RestraintCount,1)); -Init_Y=sparse(2*busNum,1);%等式约束乘子 -Gap=(Init_L'*Init_Z-Init_U'*Init_W); -% PG=sparse(busNum,1); -% PG(Balance)=0.1105; -% QG=sparse(busNum,1); -% QG(Balance)=0.0984; -% SEPD(2)=0.1105; -% SEQD(2)=0.0984; -while Gap>1e-5 && KK<20 - KK=KK+1; - Init_u=Gap/2/RestraintCount*CenterA; - deltH=func_deltH(busNum,fsY1,Loadi,Balance); - deltG=func_deltG(busNum,Loadi,I1r,I1i); -% deltG=0; - L_1Z=diag(Init_Z./Init_L); - U_1W=diag(Init_W./Init_U); - deltdeltF=func_deltdeltF(busNum,fsY1,Loadi); -% deltdeltF=0; -% ddh=func_ddh(SEVmf1,Init_Y,busNum,fsY1amp,SEVaf1,r,c,fsY1ang,Loadi,ContrlCount); - ddh=0; - ddg=func_ddg(); - deltF=func_deltF(measurement,busNum,fsY1,Loadi,V1r,V1i,I1r,I1i); -% 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(I1r,I1i); - Mat_H=FormH(fsY1,Loadi,V1r,V1i,I1r,I1i,BalI1r,BalI1i,Balance); - Ly=Mat_H; - Lz=FormLz(Loadi,Mat_G,Init_L); - Lw=FormLw(Loadi,Mat_G,Init_U); - Lx=FormLx(deltF,deltH,Init_Y,deltG,Init_Z,Init_W); - YY=FormYY(Lul,Lz,Ly,Luu,Lw,Lx); - %% 开始解方程 - plotGap(KK)=Gap; - fprintf('迭代次数 %d Gap %f\n',KK,plotGap(KK)); - XX=SolveIt(deltF,deltG,Init_L,Init_Z,Init_U,Init_W,deltdeltF,ddh,ddg,deltH,Init_Y,Ly,Lz,ContrlCount,Lw,Lul,Luu,Lx,Balance,busNum,Loadi); - [deltZ,deltL,deltW,deltU,deltX,deltY]=AssignXX(XX,ContrlCount,RestraintCount,busNum); - [Init_Z,Init_L,Init_W,Init_U,Init_Y,V1r,V1i,I1r,I1i]=Modification(Init_Z,Init_L,Init_W,Init_U,Init_Y,deltZ,deltL,deltW,deltU,deltX,deltY,V1r,V1i,I1r,I1i,ContrlCount,Balance,busNum,Loadi); - Gap=(Init_L'*Init_Z-Init_U'*Init_W); -end +[ V1r,V1i,I1r,I1i ]=IPMLoop(measurement,BalI1r,BalI1i,busNum,Loadi,fsY1,Balance ); %检查目标函数 f=sum([real(measurement);imag(measurement)]-[-I1r;-I1i]); \ No newline at end of file