facat
/
distributionnetwork-power2current
1
0
Fork 0
Code Issues Pull Requests Releases Wiki Activity

还在找错误，没调好. 

Signed-off-by: dugg@lab-desk <dugg@lab-desk>
This commit is contained in:
dugg@lab-desk 2014-10-18 22:58:30 +08:00
parent e861147774
commit 81e04f0748
5 changed files with 35 additions and 19 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

12
IPMLoop.m
View File

@ -1,4 +1,4 @@
function [ V1r,V1i,I1r,I1i ] = IPMLoop(measurement,BalI1r,BalI1i,busNum,Loadi,fsY1,Balance )
function [ V1r,V1i,I1r,I1i ] = IPMLoop(measurement,BalI1r,BalI1i,busNum,Loadi,fsY1,Balance,Vref )
%把每个序的循环写在这个函数中。其实也就是内点法循环。
V1r=1*ones(busNum,1);
V1i=0*ones(busNum,1);
@ -29,13 +29,13 @@ while Gap>1e-5 && KK<20
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);
% 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();
ddg=func_ddg(busNum,Loadi,Init_Z,Init_W);
deltF=func_deltF(measurement,busNum,fsY1,Loadi,V1r,V1i,I1r,I1i);
% deltF=0;
Luu=Init_U.*Init_W+Init_u*ones(RestraintCount,1);
@ -46,13 +46,13 @@ while Gap>1e-5 && KK<20
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);
% 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);
[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,Vref);
Gap=(Init_L'*Init_Z-Init_U'*Init_W);
end

4
Modification.m
View File

@ -1,5 +1,5 @@
function [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)
V1r,V1i,I1r,I1i,ContrlCount,Balance,busNum,Loadi,Vref)
AlphaP=FormAlphaP(Init_L,deltL,Init_U,deltU);
AlphaD=FormAlphaD(Init_Z,deltZ,Init_W,deltW);
Init_Z=Init_Z+AlphaD*deltZ;
@ -9,7 +9,7 @@ Init_U=Init_U+AlphaP*deltU;
Init_Y=Init_Y+AlphaD*deltY;
t=deltX(1:busNum*2);
V1r=V1r+AlphaP*t(1:busNum);
V1r(Balance)=1;
V1r(Balance)=Vref;
V1i=V1i+AlphaP*t(busNum+1:end);
V1i(Balance)=0;
t=deltX(busNum*2+1:ContrlCount);

2
SolveIt.m
View File

@ -1,6 +1,6 @@
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,Lx,Balance,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;
H=-deltdeltF+ddh+ddg;
t1=diag(Init_Z./Init_L-Init_W./Init_U);
t2=-deltG*( t1 )*deltG';
aa=[

13
func_ddg.m
View File

@ -1,4 +1,11 @@
function ddg=func_ddg()
ddg=0;
function ddg=func_ddg(busNum,Loadi,Init_Z,Init_W)
% ddg=0;
Init_ZW=Init_Z+Init_W;
t=[2*diag( ones(length(Loadi),1).*Init_ZW),zeros(length(Loadi));
zeros(length(Loadi)),2*diag( ones(length(Loadi),1).*Init_ZW);
];
ddg=[
zeros(busNum*2,busNum*2+length(Loadi)*2);
zeros(length(Loadi)*2,busNum*2), t;
];
end

23
run.m
View File

@ -200,6 +200,7 @@ fprintf('
% VoltpC=Vp3m(3:3:end).*exp(1j*Vp3a(3:3:end));
%% 用牛顿法求解end
%% 开始进入状态估计
clear PD QD PG QG;
%准备量测量
iterPhaseASpotLoadP=phaseASpotLoadP;
iterPhaseBSpotLoadP=phaseBSpotLoadP;
@ -220,21 +221,30 @@ If2=conj(f012(3,:)');
%试着算一下正序电流
% fsY11*V1;
%形成负荷序电流的测量值
mIf0=If0;
mIf0=-If0;
mIf1=-If1;
% mIf1(3)=-mIf1(2);
mIf2=If2;
mIf2=-If2;
% fsY11(:,Balance)=0;
% fsY11(Balance,:)=0;
% fsY11=fsY11+sparse(Balance,Balance,ones(length(Balance),1),busNum,busNum);
% mIf1(3)=1;
%% 先算正序的
%平衡节点电流
BalI1r=real(-sum(mIf1));
BalI1i=imag(-sum(mIf1));
inv(fsY11)*(mIf1);
% inv(fsY11)*(mIf1);
measurement=-mIf1(Loadi);
clear PD QD PG QG;
[ V1r,V1i,I1r,I1i ]=IPMLoop(measurement,BalI1r,BalI1i,busNum,Loadi,fsY1,Balance,1 );
f=sum(([real(measurement);imag(measurement)]-[-I1r;-I1i]).^2);
%% 算负序的
BalI2r=real(-sum(mIf2));
BalI2i=imag(-sum(mIf2));
measurement=-mIf2(Loadi);
[ V2r,V2i,I2r,I2i ]=IPMLoop(measurement,BalI2r,BalI2i,busNum,Loadi,fsY2,Balance,1 );
f=sum([real(measurement);imag(measurement)]-[-I2r;-I2i]);
%状态量
% SEVoltpA=sparse(ones(busNum,1));
% SEVoltpB=sparse(ones(busNum,1)).*exp(1j*-120/180*pi);
@ -250,6 +260,5 @@ clear PD QD PG QG;
% SEVaf1=sparse(zeros(busNum,1));
% SEPD=sparse(zeros(busNum,1));
% SEQD=sparse(zeros(busNum,1));
[ V1r,V1i,I1r,I1i ]=IPMLoop(measurement,BalI1r,BalI1i,busNum,Loadi,fsY1,Balance );
%检查目标函数
f=sum([real(measurement);imag(measurement)]-[-I1r;-I1i]);