From 276db0e00e30cffac38d8879d4e06a3cf702c904 Mon Sep 17 00:00:00 2001 From: facat Date: Wed, 23 Apr 2014 11:33:34 +0800 Subject: [PATCH] =?UTF-8?q?1.=E9=87=8D=E6=96=B0=E6=8E=92=E5=88=97=E4=BA=86?= =?UTF-8?q?=E4=BB=A3=E7=A0=81=202.=E8=80=83=E8=99=91=E4=BA=86=E6=95=B0?= =?UTF-8?q?=E6=8D=AE=E6=96=87=E4=BB=B6=E4=B8=AD=E6=B2=A1=E6=9C=89=E7=94=B5?= =?UTF-8?q?=E5=AE=B9=E5=99=A8=E9=A1=B9=E7=9A=84=E6=83=85=E5=86=B5?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Signed-off-by: facat --- dataRead.m | 6 ++++++ run.m | 44 ++++++++++++++++++++++---------------------- 2 files changed, 28 insertions(+), 22 deletions(-) diff --git a/dataRead.m b/dataRead.m index c370ee2..9769b1e 100644 --- a/dataRead.m +++ b/dataRead.m @@ -45,6 +45,7 @@ phaseBSpotLoadQ( ismember(setIJ,spotloads(:,1) ) )=spotloads(:,5); phaseCSpotLoadQ( ismember(setIJ,spotloads(:,1) ) )=spotloads(:,7); %% 补偿电容 cap=data(zeroEntries(2)+1:zeroEntries(3)-1,:); +if ~isempty(cap) capNode=nodeNum(cap(:,1)==setIJ); %数据文件中的补偿电容给的是补偿容量，要将其转换为导纳。 capB=cap(:,2:4)/1000; @@ -58,5 +59,10 @@ offSet=(capNode-1)*3+1; phaseABCY=phaseABCY+sparse(offSet,offSet,1j*capB(:,2),busNum*3,busNum*3); offSet=(capNode-1)*3+1; phaseABCY=phaseABCY+sparse(offSet,offSet,1j*capB(:,3),busNum*3,busNum*3); +else + cap=struct(); + cap.capNode=[]; + cap.capB=sparse(1,3); +end end diff --git a/run.m b/run.m index caf2d3f..40a0f0a 100644 --- a/run.m +++ b/run.m @@ -1,9 +1,9 @@ clc clear -lineZ=readLineZ('.\feeder13\lineParameter.txt'); +lineZ=readLineZ('D:\share\feeder123\lineParameter.txt'); [ fsY0, fsY1, fsY2,phaseASpotLoadP,phaseBSpotLoadP,phaseCSpotLoadP ... phaseASpotLoadQ,phaseBSpotLoadQ,phaseCSpotLoadQ,setIJ,nodeNum,Balance,phaseABCY ... - cap]=dataRead(lineZ,'.\feeder13\data1.txt'); + cap]=dataRead(lineZ,'D:\share\feeder123\data.txt'); a=exp(1j*2*pi/3); Tp2f=1/3*[1 1 1; 1 a a^2; @@ -55,18 +55,18 @@ VoltpB=sparse(ones(busNum,1)).*exp(1j*-120/180*pi); VoltpC=sparse(ones(busNum,1)).*exp(1j*+120/180*pi); while(k<=kmax && maxD> EPS) k=k+1; -%把补偿电容看作负荷 -SA=VoltpA.*conj(VoltpA.*sparse(cap.capNode,1,1j*cap.capB(:,1),busNum,1)); -SB=VoltpB.*conj(VoltpA.*sparse(cap.capNode,1,1j*cap.capB(:,2),busNum,1)); -SC=VoltpC.*conj(VoltpA.*sparse(cap.capNode,1,1j*cap.capB(:,3),busNum,1)); -iterPD=PD+real(SA+SB+SC)/3; -iterQD=QD+imag(SA+SB+SC)/3; -iterPhaseASpotLoadP=phaseASpotLoadP+real(SA); -iterPhaseBSpotLoadP=phaseBSpotLoadP+real(SB); -iterPhaseCSpotLoadP=phaseCSpotLoadP+real(SC); -iterPhaseASpotLoadQ=phaseASpotLoadQ+imag(SA); -iterPhaseBSpotLoadQ=phaseBSpotLoadQ+imag(SB); -iterPhaseCSpotLoadQ=phaseCSpotLoadQ+imag(SC); + %把补偿电容看作负荷 + SA=VoltpA.*conj(VoltpA.*sparse(cap.capNode,1,1j*cap.capB(:,1),busNum,1)); + SB=VoltpB.*conj(VoltpA.*sparse(cap.capNode,1,1j*cap.capB(:,2),busNum,1)); + SC=VoltpC.*conj(VoltpA.*sparse(cap.capNode,1,1j*cap.capB(:,3),busNum,1)); + iterPD=PD+real(SA+SB+SC)/3; + iterQD=QD+imag(SA+SB+SC)/3; + iterPhaseASpotLoadP=phaseASpotLoadP+real(SA); + iterPhaseBSpotLoadP=phaseBSpotLoadP+real(SB); + iterPhaseCSpotLoadP=phaseCSpotLoadP+real(SC); + iterPhaseASpotLoadQ=phaseASpotLoadQ+imag(SA); + iterPhaseBSpotLoadQ=phaseBSpotLoadQ+imag(SB); + iterPhaseCSpotLoadQ=phaseCSpotLoadQ+imag(SC); [dP, dQ, YdotSinVolt, YdotCosVolt, diag_Volt_YdotSin, diag_Volt_YdotCos]=Unbalance(Balance,busNum, ... PQi,PG,QG,QGi,iterPD,iterQD,Vmf1,Vaf1,fsY1amp,fsY1ang,r,c,Vf2,If2,Vf0,If0);%不平衡量 maxD=max(abs([dP;dQ;])); @@ -154,14 +154,14 @@ maxD=10000; tic while(k<=kmax && maxD> EPS) k=k+1; -[dP, dQ, YdotSinVolt, YdotCosVolt, diag_Volt_YdotSin, diag_Volt_YdotCos]=Unbalance(Balance3P,busNum*3, ... - PQi3P,PG,QG,QGi3P,PD3P,QD3P,Vp3m,Vp3a,Y,Yangle,r,c,0,0,0,0); -maxD=max(abs([dP;dQ;])); -jaco=Jacobi(Balance3P,busNum*3,QGi3P,Vp3m,YdotSinVolt,YdotCosVolt,diag_Volt_YdotSin,diag_Volt_YdotCos);%雅克比矩阵 -[dV, dVangle]=Solv(busNum*3,jaco,dP,dQ);%解出修正量 -[Vp3m, Vp3a]=Modify(Vp3m,Vp3a,dV,dVangle,1); -fprintf('第 %d 次迭代， 最大不平衡量为 %f\n',k,full(maxD)); -fprintf('迭代时间%f\n',toc); + [dP, dQ, YdotSinVolt, YdotCosVolt, diag_Volt_YdotSin, diag_Volt_YdotCos]=Unbalance(Balance3P,busNum*3, ... + PQi3P,PG,QG,QGi3P,PD3P,QD3P,Vp3m,Vp3a,Y,Yangle,r,c,0,0,0,0); + maxD=max(abs([dP;dQ;])); + jaco=Jacobi(Balance3P,busNum*3,QGi3P,Vp3m,YdotSinVolt,YdotCosVolt,diag_Volt_YdotSin,diag_Volt_YdotCos);%雅克比矩阵 + [dV, dVangle]=Solv(busNum*3,jaco,dP,dQ);%解出修正量 + [Vp3m, Vp3a]=Modify(Vp3m,Vp3a,dV,dVangle,1); + fprintf('第 %d 次迭代， 最大不平衡量为 %f\n',k,full(maxD)); + fprintf('迭代时间%f\n',toc); end NewtonToc=toc; fprintf('牛顿法计算时间 %f\n',NewtonToc);