From d39075c0aa548d8260304e619828885f1919c1e3 Mon Sep 17 00:00:00 2001 From: facat Date: Sat, 17 Aug 2013 17:43:42 +0800 Subject: [PATCH] =?UTF-8?q?1.=E6=9B=B4=E6=96=B0=E5=8F=98=E5=8E=8B=E5=99=A8?= =?UTF-8?q?=E5=8A=9F=E7=8E=87=E7=9A=84=E5=87=BD=E6=95=B0=202.=E5=88=A0?= =?UTF-8?q?=E6=8E=89=E6=97=A0=E7=94=A8=E4=BB=A3=E7=A0=81=EF=BC=8C=E6=9B=B4?= =?UTF-8?q?=E6=98=93=E9=98=85=E8=AF=BB?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Signed-off-by: facat --- TransPower.m | 6 ++++-- TransReactivePower.m | 6 ++++-- run.m | 46 +++++++++----------------------------------- 3 files changed, 17 insertions(+), 41 deletions(-) diff --git a/TransPower.m b/TransPower.m index 58a5793..13f365b 100644 --- a/TransPower.m +++ b/TransPower.m @@ -1,11 +1,13 @@ -function [ output_args ] = TransPower( newwordParameter,cmpY,Volt,VAngel ) +function [ output_args ] = TransPower( newwordParameter,Volt,VAngel ) %TRANSPOWER Summary of this function goes here % Detailed explanation goes here transI=newwordParameter.trans.transI; transJ=newwordParameter.trans.transJ; +transR=newwordParameter.trans.transR; +transX=newwordParameter.trans.transX; cmpSEV=Volt.*exp(1j*VAngel); %复数电压 % cmpY=sparse(r,c,Y.*(1j*YAngel),length(Volt),length(Volt)); -cmpSETransI=(cmpSEV(transI)-cmpSEV(transJ)).*cmpY(sub2ind(size(cmpY),transI,transJ)); +cmpSETransI=(cmpSEV(transI)-cmpSEV(transJ)).*1./(transR+1j*transX); output_args=real(cmpSEV(transI).*conj(cmpSETransI)); end diff --git a/TransReactivePower.m b/TransReactivePower.m index 706d53b..b8a5ec1 100644 --- a/TransReactivePower.m +++ b/TransReactivePower.m @@ -1,11 +1,13 @@ -function [ output_args ] = TransReactivePower( newwordParameter,cmpY,Volt,VAngel ) +function [ output_args ] = TransReactivePower( newwordParameter,Volt,VAngel ) %TRANSREAVTIVEPOWER Summary of this function goes here % Detailed explanation goes here transI=newwordParameter.trans.transI; transJ=newwordParameter.trans.transJ; +transR=newwordParameter.trans.transR; +transX=newwordParameter.trans.transX; cmpSEV=Volt.*exp(1j*VAngel); %复数电压 % cmpY=sparse(r,c,Y.*(1j*YAngel),length(Volt),length(Volt)); -cmpSETransI=(cmpSEV(transI)-cmpSEV(transJ)).*cmpY(sub2ind(size(cmpY),transI,transJ)); +cmpSETransI=(cmpSEV(transI)-cmpSEV(transJ)).*1./(transR+1j*transX); output_args=imag(cmpSEV(transI).*conj(cmpSETransI)); end diff --git a/run.m b/run.m index 6d7155b..2040b8f 100644 --- a/run.m +++ b/run.m @@ -33,25 +33,14 @@ lineX=newwordParameter.line.lineX; lineB2=newwordParameter.line.lineB2; lineG=real(1./(lineR+1j*lineX)); lineB=imag(1./(lineR+1j*lineX)); -% 处理线路电阻或电抗为0的情况，即消除NaN -% zerosLRXInd=union(find(abs(lineR)<1e-5),find(abs(lineX)<1e-5)); cmpBranchI=BranchI( cmpV,lineI,lineJ,lineR,lineX );%复数支路电流 rBranchI=abs(cmpBranchI);% 支路电流幅值 mBranchI=rBranchI.*(normrnd(0,sigma,length(rBranchI),1)+1);%支路电流量测量 %% 支路功率 rBranchP= BranchP( cmpV,cmpBranchI,lineI,lineB2 ); -% rV=abs(cmpV); -% rA=angle(cmpV); -% a=rV(lineI).^2.*(lineG)-rV(lineI).*rV(lineJ).*(lineG.*cos(rA(lineI)-rA(lineJ) )+lineB.*sin(rA(lineI)-rA(lineJ) ) ); -% b=-rV(lineI).^2.*(lineB+lineB2)-rV(lineI).*rV(lineJ).*(lineG.*sin(rA(lineI)-rA(lineJ) )-lineB.*cos(rA(lineI)-rA(lineJ) ) ); -% rBranchP=rBranchP(abs(rBranchP)>1e-5); mBranchP=rBranchP.*(normrnd(0,sigma,length(rBranchP),1)+1);%支路功率量测量 rBranchQ=BranchQ( cmpV,cmpBranchI,lineI,lineB2 ); -% rBranchQ=rBranchQ(abs(rBranchQ)>1e-5); mBranchQ=rBranchQ.*(normrnd(0,sigma,length(rBranchQ),1)+1);%支路功率量测量 -% 处理线路电阻或电抗为0的情况，即消除NaN -% mBranchP=mBranchP(setdiff( 1:length(mBranchP),zerosLRXInd)); -% mBranchQ=mBranchQ(setdiff( 1:length(mBranchQ),zerosLRXInd)); %% 注入功率 rPD=PD; PDi=find(PD~=0); @@ -76,12 +65,8 @@ transR=newwordParameter.trans.transR; transX=newwordParameter.trans.transX; transG=real(1./(transR+1j*transX)); transB=imag(1./(transR+1j*transX)); -% 处理线路电阻或电抗为0的情况，即消除NaN -% zerosTRXInd=union(find(abs(transR)<1e-5),find(abs(transX)<1e-5)); -rTransP=TransPower( newwordParameter,cmpY,rVolt,rVAngel ); -rTransQ=TransReactivePower( newwordParameter,cmpY,rVolt,rVAngel ); -% rTransP=rTransP( setdiff( 1:length(rTransP),zerosTRXInd) ); -% rTransQ=rTransQ( setdiff( 1:length(rTransQ),zerosTRXInd) ); +rTransP=TransPower( newwordParameter,rVolt,rVAngel ); +rTransQ=TransReactivePower( newwordParameter,rVolt,rVAngel ); mTransP=rTransP.*(normrnd(0,sigma,length(rTransP),1)+1); mTransQ=rTransQ.*(normrnd(0,sigma,length(rTransQ),1)+1); %% 发电机注入功率 @@ -210,33 +195,20 @@ while maxD>1e-3; transG.*cos(SEVAngel(transI)-SEVAngel(transJ)) +transB.*sin(SEVAngel(transI)-SEVAngel(transJ))... ) ... ,length(transI),length(mVolt));%变压器 - % 处理线路电阻或电抗为0的情况，即消除NaN -% dLPij_dVi=dLPij_dVi(setdiff( 1:size(dLPij_dVi,1),zerosLRXInd),:); -% dLPij_dVj=dLPij_dVj(setdiff( 1:size(dLPij_dVj,1),zerosLRXInd),:); -% dLPij_dThetai=dLPij_dThetai(setdiff( 1:size(dLPij_dThetai,1),zerosLRXInd),:); -% dLPij_dThetaj=dLPij_dThetaj(setdiff( 1:size(dLPij_dThetaj,1),zerosLRXInd),:); -% dLQij_dVi=dLQij_dVi(setdiff( 1:size(dLQij_dVi,1),zerosLRXInd),:); -% dLQij_dVj=dLQij_dVj(setdiff( 1:size(dLQij_dVj,1),zerosLRXInd),:); -% dLQij_dThetai=dLQij_dThetai(setdiff( 1:size(dLQij_dThetai,1),zerosLRXInd),:); -% dLQij_dThetaj=dLQij_dThetaj(setdiff( 1:size(dLQij_dThetaj,1),zerosLRXInd),:); - % 对量测值做同样处理 - %% 进入迭代 % H=[dV_dV ,dV_dTyta; % dLPij_dVi+dLPij_dVj, dLPij_dThetai+dLPij_dThetaj; - % dLQij_dVi+dLQij_dVj, dLQij_dThetai+dLQij_dThetaj];%jacobi - + % dLQij_dVi+dLQij_dVj, dLQij_dThetai+dLQij_dThetaj];%jacobi H=[dV_dV,dV_dTyta; dLPij_dVi+dLPij_dVj,dLPij_dThetai+dLPij_dThetaj ; - dLQij_dVi+dLQij_dVj,dLQij_dThetai+dLQij_dThetaj ];%jacobi + dLQij_dVi+dLQij_dVj,dLQij_dThetai+dLQij_dThetaj ; + dTPij_dVi+dTPij_dVj,dTPij_dThetai+dTPij_dThetaj; + dTQij_dVi+dTQij_dVj,dTQij_dThetai+dTQij_dThetaj];%jacobi SEBranchI=BranchI( SEVolt.*exp(1j*SEVAngel),lineI,lineJ,lineR,lineX );%复数支路电流 SEBranchP=BranchP( SEVolt.*exp(1j*SEVAngel),SEBranchI,lineI,lineB2 ); SEBranchQ=BranchQ( SEVolt.*exp(1j*SEVAngel),SEBranchI,lineI,lineB2 ); -% SEBranchP=SEBranchP(setdiff( 1:length(SEBranchP),zerosLRXInd)); -% SEBranchQ=SEBranchQ(setdiff( 1:length(SEBranchQ),zerosLRXInd)); - % h=[SEVolt;SEBranchP;SEBranchQ]; - % z=[mVolt;mBranchP;mBranchQ]; + SETransP=TransPower( newwordParameter,cmpY,rVolt,rVAngel ); h=[SEVolt;SEBranchP;SEBranchQ]; z=[mVolt;mBranchP;mBranchP]; W=sparse(1:length(h),1:length(h),1/sigma.^2,length(h),length(h)); @@ -250,9 +222,9 @@ while maxD>1e-3; dX=G\-g; maxD=max(abs(dX)) % 更新变量 - % SEVolt=SEVolt+dX(1:length(mVolt)); + SEVolt=SEVolt+dX(1:length(mVolt)); Iteration=Iteration+1; - SEVAngel=SEVAngel+dX(1:end); + SEVAngel=SEVAngel+dX(length(mVolt)+1:end); end %% 输出结果 fprintf('迭代%d次\n',Iteration);