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

之前搞错了，其实可以用量测量做负荷曲线。之前是把三序顺序搞错了。 

Signed-off-by: dmy@lab <dmy@lab.lab>
This commit is contained in:
dmy@lab 2015-04-02 17:15:39 +08:00
parent b6a1f8058c
commit fa03721b18
4 changed files with 29 additions and 18 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

6
FormLwstate1.m
View File

@ -1,9 +1,9 @@
function Lw=FormLwstate1(Loadi,Mat_G,Init_U,I1measurement,dI_F,flag,rIf1) function Lw=FormLwstate1(Loadi,Mat_G,Init_U,I1measurement,dI_F,flag)
Ir=real(rIf1(Loadi)); Ir=real(I1measurement);
pIr=find(Ir>0); pIr=find(Ir>0);
nIr=find(Ir<0); nIr=find(Ir<0);
Ii=imag(rIf1(Loadi)); Ii=imag(I1measurement);
pIi=find(Ii>0); pIi=find(Ii>0);
nIi=find(Ii<0); nIi=find(Ii<0);

6
FormLzstate1.m
View File

@ -1,8 +1,8 @@
function Lz=FormLzstate1(Loadi,Mat_G,Init_L,I1measurement,dI_F,flag,rIf1) function Lz=FormLzstate1(Loadi,Mat_G,Init_L,I1measurement,dI_F,flag)
Ir=real(rIf1(Loadi)); Ir=real(I1measurement);
pIr=find(Ir>0); pIr=find(Ir>0);
nIr=find(Ir<0); nIr=find(Ir<0);
Ii=imag(rIf1(Loadi)); Ii=imag(I1measurement);
pIi=find(Ii>0); pIi=find(Ii>0);
nIi=find(Ii<0); nIi=find(Ii<0);

8
IPMLoop.m
View File

@ -1,4 +1,4 @@
function [ V1r,V1i,I1r,I1i,isConverged ] = IPMLoop(V1measurement,wV1r,wV1i,I1measurement,wI1r,wI1i,BalI1r,BalI1i,busNum,Loadi,fsY1,Balance,Vref,dI_F,flag,rIf1 ) function [ V1r,V1i,I1r,I1i,isConverged ] = IPMLoop(V1measurement,wV1r,wV1i,I1measurement,wI1r,wI1i,BalI1r,BalI1i,busNum,Loadi,fsY1,Balance,Vref,dI_F,flag )
%把每个序的循环写在这个函数中。其实也就是内点法循环。 %把每个序的循环写在这个函数中。其实也就是内点法循环。
V1r=Vref*ones(busNum,1); V1r=Vref*ones(busNum,1);
V1i=0*ones(busNum,1); V1i=0*ones(busNum,1);
@ -17,7 +17,7 @@ Init_L=1*sparse(ones(RestraintCount,1));
Init_U=1*sparse(ones(RestraintCount,1)); Init_U=1*sparse(ones(RestraintCount,1));
Init_Y=sparse(2*busNum,1);%等式约束乘子 Init_Y=sparse(2*busNum,1);%等式约束乘子
Gap=(Init_L'*Init_Z-Init_U'*Init_W); Gap=(Init_L'*Init_Z-Init_U'*Init_W);
kmax=20; kmax=10;
while(abs(Gap)>0.00001) while(abs(Gap)>0.00001)
if KK>=kmax if KK>=kmax
break; break;
@ -49,8 +49,8 @@ while(abs(Gap)>0.00001)
Mat_G=FormGstate1(I1r,I1i); Mat_G=FormGstate1(I1r,I1i);
Mat_H=FormH(fsY1,Loadi,V1r,V1i,I1r,I1i,BalI1r,BalI1i,Balance); Mat_H=FormH(fsY1,Loadi,V1r,V1i,I1r,I1i,BalI1r,BalI1i,Balance);
Ly=Mat_H; Ly=Mat_H;
Lz=FormLzstate1(Loadi,Mat_G,Init_L,I1measurement,dI_F,flag,rIf1); Lz=FormLzstate1(Loadi,Mat_G,Init_L,I1measurement,dI_F,flag);
Lw=FormLwstate1(Loadi,Mat_G,Init_U,I1measurement,dI_F,flag,rIf1); Lw=FormLwstate1(Loadi,Mat_G,Init_U,I1measurement,dI_F,flag);
Lx=FormLx(deltF,deltH,Init_Y,deltG,Init_Z,Init_W); 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);
%% 开始解方程 %% 开始解方程

27
run.m
View File

@ -328,9 +328,9 @@ wV0i=wV0i(Loadi);
% dI_F=conj(dI_F)'; % dI_F=conj(dI_F)';
% dI_F=dI_F(Loadi,:); % dI_F=dI_F(Loadi,:);
% 换一个方式来计算 % 换一个方式来计算
deltmCurpA=conj((phaseASpotLoadP+1j*phaseASpotLoadQ)./rVoltpA)*0.2; deltmCurpA=conj((mphaseASpotLoadP+1j*mphaseASpotLoadQ)./mVoltpA)*0.2;
deltmCurpB=conj((phaseBSpotLoadP+1j*phaseBSpotLoadQ)./rVoltpB)*0.2; deltmCurpB=conj((mphaseBSpotLoadP+1j*mphaseBSpotLoadQ)./mVoltpB)*0.2;
deltmCurpC=conj((phaseCSpotLoadP+1j*phaseCSpotLoadQ)./rVoltpC)*0.2; deltmCurpC=conj((mphaseCSpotLoadP+1j*mphaseCSpotLoadQ)./mVoltpC)*0.2;
%转换为序电流 %转换为序电流
deltF012=Tp2f*conj([deltmCurpA';deltmCurpB';deltmCurpC']); deltF012=Tp2f*conj([deltmCurpA';deltmCurpB';deltmCurpC']);
dI_F=conj(deltF012)'; dI_F=conj(deltF012)';
@ -353,6 +353,17 @@ dI_F=dI_F(Loadi,:);
% dI_F=conj(deltF012)'; % dI_F=conj(deltF012)';
% dI_F=dI_F(Loadi,:); % dI_F=dI_F(Loadi,:);
%考虑负荷曲线不准确的情况
% curveIf1=rIf1.*(1+unifrnd(-0.2,0.2,length(rIf1),1));
% curveIf2=rIf2.*(1+unifrnd(-0.2,0.2,length(rIf2),1));
% curveIf0=rIf0.*(1+unifrnd(-0.2,0.2,length(rIf0),1));
%
% curveCurpA=conj((phaseASpotLoadP.*(1+unifrnd(-0.2,0.2,length(phaseASpotLoadP),1))+1j*phaseASpotLoadQ.*(1+unifrnd(-0.2,0.2,length(phaseASpotLoadQ),1)))./rVoltpA)*0.2;
% curveCurpB=conj((phaseBSpotLoadP.*(1+unifrnd(-0.2,0.2,length(phaseBSpotLoadP),1))+1j*phaseBSpotLoadQ.*(1+unifrnd(-0.2,0.2,length(phaseBSpotLoadQ),1)))./rVoltpB)*0.2;
% curveCurpC=conj((phaseCSpotLoadP.*(1+unifrnd(-0.2,0.2,length(phaseCSpotLoadP),1))+1j*phaseCSpotLoadQ.*(1+unifrnd(-0.2,0.2,length(phaseCSpotLoadQ),1)))./rVoltpC)*0.2;
% deltF012=Tp2f*conj([curveCurpA';curveCurpB';curveCurpC']);
tic tic
for II=1:3 for II=1:3
@ -360,19 +371,19 @@ for II=1:3
if II==1 if II==1
fprintf('正序\n'); fprintf('正序\n');
tic tic
[ V1r,V1i,I1r,I1i,isConverged1 ]=IPMLoop(V1measurement,wV1r,wV1i,I1measurement,wI1r,wI1i,BalI1r,BalI1i,busNum,Loadi,fsY11,Balance,1,dI_F,II,rIf1 ); [ V1r,V1i,I1r,I1i,isConverged1 ]=IPMLoop(V1measurement,wV1r,wV1i,I1measurement,wI1r,wI1i,BalI1r,BalI1i,busNum,Loadi,fsY11,Balance,1,dI_F,2 );
toc toc
end end
if II==2 if II==2
fprintf('负序\n'); fprintf('负序\n');
tic tic
[ V2r,V2i,I2r,I2i,isConverged2 ]=IPMLoop(V2measurement,wV2r,wV2i,I2measurement,wI2r,wI2i,BalI2r,BalI2i,busNum,Loadi,fsY22,Balance,0,dI_F,II,rIf2 ); [ V2r,V2i,I2r,I2i,isConverged2 ]=IPMLoop(V2measurement,wV2r,wV2i,I2measurement,wI2r,wI2i,BalI2r,BalI2i,busNum,Loadi,fsY22,Balance,0,dI_F,3 );
toc toc
end end
if II==3 if II==3
fprintf('零序\n'); fprintf('零序\n');
tic tic
[ V0r,V0i,I0r,I0i,isConverged0 ]=IPMLoop(V0measurement,wV0r,wV0i,I0measurement,wI0r,wI0i,BalI0r,BalI0i,busNum,Loadi,fsY00,Balance,0,dI_F,II,rIf0 ); [ V0r,V0i,I0r,I0i,isConverged0 ]=IPMLoop(V0measurement,wV0r,wV0i,I0measurement,wI0r,wI0i,BalI0r,BalI0i,busNum,Loadi,fsY00,Balance,0,dI_F,1 );
toc toc
end end
@ -393,7 +404,7 @@ SEVoltpABC=Tf2p*conj([ SEVf0'; SEVf1'; SEVf2']);
SEIf0=I0r+1j*I0i; SEIf0=I0r+1j*I0i;
SEIf1=I1r+1j*I1i; SEIf1=I1r+1j*I1i;
SEIf2=I2r+1j*I2i; SEIf2=I2r+1j*I2i;
SEIpABC=full(Tf2p*conj([SEIf0';SEIf1(1:3)';SEIf2'])); SEIpABC=full(Tf2p*conj([SEIf0';SEIf1';SEIf2']));
%看一下差多少 %看一下差多少
VError=(full(abs(VoltpABC))-abs(SEVoltpABC))./abs(VoltpABC)*100; VError=(full(abs(VoltpABC))-abs(SEVoltpABC))./abs(VoltpABC)*100;
VError=reshape(VError,size(VError,1)*size(VError,2),1); VError=reshape(VError,size(VError,1)*size(VError,2),1);
@ -404,7 +415,7 @@ rThreeLoad=[ phaseASpotLoadP'+1j*phaseASpotLoadQ';
phaseCSpotLoadP'+1j*phaseCSpotLoadQ'; phaseCSpotLoadP'+1j*phaseCSpotLoadQ';
]; ];
% rThreeLoad=rThreeLoad(:,setxor(1:size(SEVoltpABC,2),Balance)); % rThreeLoad=rThreeLoad(:,setxor(1:size(SEVoltpABC,2),Balance));
SEThreeLoad=SEVoltpABC(:,setxor(1:size(SEVoltpABC,2),Balance)).*conj(-SEIpABC); SEThreeLoad=SEVoltpABC(:,Loadi).*conj(-SEIpABC);
% phaseLoadPError=real(rThreeLoad-SEThreeLoad)./real(rThreeLoad)*100; % phaseLoadPError=real(rThreeLoad-SEThreeLoad)./real(rThreeLoad)*100;
% phaseLoadQError=imag(rThreeLoad-SEThreeLoad)./imag(rThreeLoad)*100; % phaseLoadQError=imag(rThreeLoad-SEThreeLoad)./imag(rThreeLoad)*100;
% phaseLoadPError=reshape(phaseLoadPError,size(phaseLoadPError,1)*size(phaseLoadPError,2),1); % phaseLoadPError=reshape(phaseLoadPError,size(phaseLoadPError,1)*size(phaseLoadPError,2),1);