diff --git a/FormLwstate1.m b/FormLwstate1.m index 2dbb78e..ae3459d 100644 --- a/FormLwstate1.m +++ b/FormLwstate1.m @@ -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); nIr=find(Ir<0); -Ii=imag(rIf1(Loadi)); +Ii=imag(I1measurement); pIi=find(Ii>0); nIi=find(Ii<0); diff --git a/FormLzstate1.m b/FormLzstate1.m index fac89fc..de28d2a 100644 --- a/FormLzstate1.m +++ b/FormLzstate1.m @@ -1,8 +1,8 @@ -function Lz=FormLzstate1(Loadi,Mat_G,Init_L,I1measurement,dI_F,flag,rIf1) -Ir=real(rIf1(Loadi)); +function Lz=FormLzstate1(Loadi,Mat_G,Init_L,I1measurement,dI_F,flag) +Ir=real(I1measurement); pIr=find(Ir>0); nIr=find(Ir<0); -Ii=imag(rIf1(Loadi)); +Ii=imag(I1measurement); pIi=find(Ii>0); nIi=find(Ii<0); diff --git a/IPMLoop.m b/IPMLoop.m index 5dd0e51..8138c96 100644 --- a/IPMLoop.m +++ b/IPMLoop.m @@ -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); V1i=0*ones(busNum,1); @@ -17,7 +17,7 @@ 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); -kmax=20; +kmax=10; while(abs(Gap)>0.00001) if KK>=kmax break; @@ -49,8 +49,8 @@ while(abs(Gap)>0.00001) Mat_G=FormGstate1(I1r,I1i); Mat_H=FormH(fsY1,Loadi,V1r,V1i,I1r,I1i,BalI1r,BalI1i,Balance); Ly=Mat_H; - Lz=FormLzstate1(Loadi,Mat_G,Init_L,I1measurement,dI_F,flag,rIf1); - Lw=FormLwstate1(Loadi,Mat_G,Init_U,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); Lx=FormLx(deltF,deltH,Init_Y,deltG,Init_Z,Init_W); YY=FormYY(Lul,Lz,Ly,Luu,Lw,Lx); %% 开始解方程 diff --git a/run.m b/run.m index 165c84c..da27642 100644 --- a/run.m +++ b/run.m @@ -328,9 +328,9 @@ wV0i=wV0i(Loadi); % dI_F=conj(dI_F)'; % dI_F=dI_F(Loadi,:); % 换一个方式来计算 -deltmCurpA=conj((phaseASpotLoadP+1j*phaseASpotLoadQ)./rVoltpA)*0.2; -deltmCurpB=conj((phaseBSpotLoadP+1j*phaseBSpotLoadQ)./rVoltpB)*0.2; -deltmCurpC=conj((phaseCSpotLoadP+1j*phaseCSpotLoadQ)./rVoltpC)*0.2; +deltmCurpA=conj((mphaseASpotLoadP+1j*mphaseASpotLoadQ)./mVoltpA)*0.2; +deltmCurpB=conj((mphaseBSpotLoadP+1j*mphaseBSpotLoadQ)./mVoltpB)*0.2; +deltmCurpC=conj((mphaseCSpotLoadP+1j*mphaseCSpotLoadQ)./mVoltpC)*0.2; %转换为序电流 deltF012=Tp2f*conj([deltmCurpA';deltmCurpB';deltmCurpC']); dI_F=conj(deltF012)'; @@ -353,6 +353,17 @@ dI_F=dI_F(Loadi,:); % dI_F=conj(deltF012)'; % 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 for II=1:3 @@ -360,19 +371,19 @@ for II=1:3 if II==1 fprintf('正序\n'); 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 end if II==2 fprintf('负序\n'); 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 end if II==3 fprintf('零序\n'); 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 end @@ -393,7 +404,7 @@ SEVoltpABC=Tf2p*conj([ SEVf0'; SEVf1'; SEVf2']); SEIf0=I0r+1j*I0i; SEIf1=I1r+1j*I1i; 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=reshape(VError,size(VError,1)*size(VError,2),1); @@ -404,7 +415,7 @@ rThreeLoad=[ phaseASpotLoadP'+1j*phaseASpotLoadQ'; phaseCSpotLoadP'+1j*phaseCSpotLoadQ'; ]; % 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; % phaseLoadQError=imag(rThreeLoad-SEThreeLoad)./imag(rThreeLoad)*100; % phaseLoadPError=reshape(phaseLoadPError,size(phaseLoadPError,1)*size(phaseLoadPError,2),1);