From 5e07471ce34825ad289edc4f9850dcadfeea3a6e Mon Sep 17 00:00:00 2001
From: facat <facat@facat.cn>
Date: Thu, 15 Aug 2013 21:51:31 +0800
Subject: [PATCH] =?UTF-8?q?=E6=8A=8A=E7=BA=BF=E8=B7=AF=E7=94=B5=E6=B5=81?=
 =?UTF-8?q?=EF=BC=8C=E5=8A=9F=E7=8E=87=E7=9A=84=E8=AE=A1=E7=AE=97=E5=86=99?=
 =?UTF-8?q?=E6=88=90=E5=87=BD=E6=95=B0=E5=BD=A2=E5=BC=8F?=
MIME-Version: 1.0
Content-Type: text/plain; charset=UTF-8
Content-Transfer-Encoding: 8bit

Signed-off-by: facat <facat@facat.cn>
---
 BranchI.m |  5 +++++
 BranchP.m |  4 ++++
 BranchQ.m |  4 ++++
 run.m     | 13 ++++++++-----
 4 files changed, 21 insertions(+), 5 deletions(-)
 create mode 100644 BranchI.m
 create mode 100644 BranchP.m
 create mode 100644 BranchQ.m

diff --git a/BranchI.m b/BranchI.m
new file mode 100644
index 0000000..8e53f36
--- /dev/null
+++ b/BranchI.m
@@ -0,0 +1,5 @@
+function [ output_args ] = BranchI( V,lineI,lineJ,lineR,lineX )
+output_args=(V(lineI)-V(lineJ))./(lineR+1j*lineX);
+
+end
+
diff --git a/BranchP.m b/BranchP.m
new file mode 100644
index 0000000..5028e94
--- /dev/null
+++ b/BranchP.m
@@ -0,0 +1,4 @@
+function [ output_args ] = BranchP( V,I,lineI,lineB2 )
+output_args=real((V(lineI)).*conj(I))+real(V(lineI) .*conj(1j*lineB2.*V(lineI) ) );
+end
+
diff --git a/BranchQ.m b/BranchQ.m
new file mode 100644
index 0000000..5aafde9
--- /dev/null
+++ b/BranchQ.m
@@ -0,0 +1,4 @@
+function [ output_args ] = BranchQ( V,I,lineI,lineB2 )
+output_args=imag((V(lineI)).*conj(I))+imag(V(lineI) .*conj(1j*lineB2.*V(lineI) ) );
+end
+
diff --git a/run.m b/run.m
index 75ada0c..578125f 100644
--- a/run.m
+++ b/run.m
@@ -33,18 +33,18 @@ lineX=newwordParameter.line.lineX;
 lineB2=newwordParameter.line.lineB2;
 lineG=real(1./(lineR+1j*lineX));
 lineB=imag(1./(lineR+1j*lineX));
-cmpBranchI=(cmpV(lineI)-cmpV(lineJ))./(lineR+1j*lineX);%复数支路电流
+cmpBranchI=BranchI( cmpV,lineI,lineJ,lineR,lineX );%复数支路电流
 rBranchI=abs(cmpBranchI);% 支路电流幅值
 mBranchI=rBranchI.*(normrnd(0,sigma,length(rBranchI),1)+1);%支路电流量测量
 %% 支路功率
-rBranchP=real((cmpV(lineI)).*conj(cmpBranchI))+real(cmpV(lineI) .*conj(1j*lineB2.*cmpV(lineI) ) );
+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=imag((cmpV(lineI)).*conj(cmpBranchI))+imag(cmpV(lineI) .*conj(1j*lineB2.*cmpV(lineI) ) );
+rBranchQ=BranchQ( cmpV,cmpBranchI,lineI,lineB2 );
 % rBranchQ=rBranchQ(abs(rBranchQ)>1e-5);
 mBranchQ=rBranchQ.*(normrnd(0,sigma,length(rBranchQ),1)+1);%支路功率量测量
 %% 注入功率
@@ -159,13 +159,16 @@ dLQij_dVi=dLQij_dVi(setdiff( 1:size(dLQij_dVi,1),zerosRXInd),:);
 dLQij_dVj=dLQij_dVj(setdiff( 1:size(dLQij_dVj,1),zerosRXInd),:);
 dLQij_dThetai=dLQij_dThetai(setdiff( 1:size(dLQij_dThetai,1),zerosRXInd),:);
 dLQij_dThetaj=dLQij_dThetaj(setdiff( 1:size(dLQij_dThetaj,1),zerosRXInd),:);
+% 对量测值做同样处理
+mBranchP=mBranchP(setdiff( 1:length(mBranchP),zerosRXInd));
+mBranchQ=mBranchQ(setdiff( 1:length(mBranchQ),zerosRXInd));
 %% 进入迭代
 H=[dV_dV ,dV_dTyta;
   dLPij_dVi+dLPij_dVj, dLPij_dThetai+dLPij_dThetaj;
   dLQij_dVi+dLQij_dVj, dLQij_dThetai+dLQij_dThetaj];%jacobi
 % H=[dV_dV];
-h=SEVolt;
-z=mVolt;
+h=[SEVolt;mBranchP;mBranchQ];
+z=[mVolt;mBranchP;mBranchQ];
 W=sparse(1:length(mVolt),1:length(mVolt),1/sigma.^2,length(mVolt),length(mVolt));
 G=H'*W*H;
 g=-H'*W*(z-h);