From 006adb9fc008b44e5ef3aacad18638f31597db3a Mon Sep 17 00:00:00 2001 From: facat Date: Wed, 14 Aug 2013 10:46:45 +0800 Subject: [PATCH] =?UTF-8?q?=E5=87=86=E5=A4=87=E8=87=AA=E5=B7=B1=E5=86=99?= =?UTF-8?q?=E4=BB=A3=E7=A0=81=E6=B1=82=E5=AF=BC?= MIME-Version: 1.0 Content-Type: text/plain; charset=UTF-8 Content-Transfer-Encoding: 8bit Signed-off-by: facat --- @SEOpti/SEOpti.m | 4 ++++ @SEOpti/equ.m | 4 +++- @SEOpti/fun.asv | 31 +++++++++++++++++++++++-------- @SEOpti/fun.m | 25 +++++++++++++++++++------ @SEOpti/init.m | 7 +++++-- mBranchI.mat | Bin 374 -> 374 bytes mBranchP.mat | Bin 357 -> 357 bytes mBranchQ.mat | Bin 373 -> 373 bytes mPG.mat | Bin 196 -> 196 bytes mQG.mat | Bin 218 -> 217 bytes mTransP.mat | Bin 205 -> 205 bytes mTransQ.mat | Bin 206 -> 206 bytes mVolt.mat | Bin 332 -> 334 bytes run.m | 25 ++++++++++++++----------- 14 files changed, 68 insertions(+), 28 deletions(-) diff --git a/@SEOpti/SEOpti.m b/@SEOpti/SEOpti.m index d28e6d5..1963266 100644 --- a/@SEOpti/SEOpti.m +++ b/@SEOpti/SEOpti.m @@ -25,6 +25,10 @@ classdef SEOpti transI=NaN; transJ=NaN; newwordParameter=NaN; + BalanceVolt=NaN; + %真实值 + rVolt=NaN; + rVAngel end methods diff --git a/@SEOpti/equ.m b/@SEOpti/equ.m index 3e286bb..3b5d007 100644 --- a/@SEOpti/equ.m +++ b/@SEOpti/equ.m @@ -11,7 +11,9 @@ PQ=diag(cmpSEV)*conj(this.cmpY*cmpSEV); out_arg=[real(PQ(this.zerosInjectionIndex));imag(PQ(this.zerosInjectionIndex));]; % out_arg=[PQ1(this.zerosInjectionIndex);PQ2(this.zerosInjectionIndex);]; % out_arg=[out_arg;x(this.Balance+this.Busnum)]; -out_arg=[out_arg;x(this.Balance+this.Busnum)]; +out_arg=[out_arg;x(this.Balance+this.Busnum)];%相角等于0 +BalanceVolt=this.BalanceVolt; +out_arg=[out_arg;SEVolt(this.Balance)-BalanceVolt];%平衡节点电压 out_arg=full(out_arg); this.mnle=zeros(length(out_arg),1); this.mnlrhs=this.mnle; diff --git a/@SEOpti/fun.asv b/@SEOpti/fun.asv index e3682ed..c41e226 100644 --- a/@SEOpti/fun.asv +++ b/@SEOpti/fun.asv @@ -4,28 +4,43 @@ function [ Objective ] = fun(this, x ) Objective=0; SEVolt=x(1:this.Busnum); SEVAngel=x(this.Busnum+1:2*this.Busnum); -Objective=sum(((SEVolt-this.mVolt)./this.sigma).^2);%电压 +rVolt=this.rVolt; +sigma=this.sigma; +VoltSigma=rVolt*sigma; +Objective=sum(((SEVolt-this.mVolt)./VoltSigma).^2);%电压 %% 支路电流 cmpSEV=SEVolt.*exp(1j*SEVAngel); %复数电压 cmpSEBranchI=(cmpSEV(this.lineI)-cmpSEV(this.lineJ))./(this.lineR+1j*this.lineX);%复数支路电流 SEBranchI=abs(cmpSEBranchI);% 支路电流幅值 -Objective=Objective+sum((SEBranchI-this.mBranchI.^2).^2./this.sigma^2);%%电流 +% (SEBranchI-this.mBranchI).^2./(this.sigma^2) +% Objective=Objective+sum((SEBranchI-this.mBranchI).^2./(this.sigma^2));%%电流 %% 线路功率 +rVAngel=this.rVAngel; +cmpV=rVolt.*exp(1j*rVAngel); +rBranchPQ=(cmpV(this.lineI)-cmpV(this.lineJ)).*conj( (cmpV(this.lineI)-cmpV(this.lineJ))./(this.lineR+1j*this.lineX) ); +SEBranchPSigma=real(rBranchPQ)*sigma; +indBranchPS=find(SEBranchPSigma>1e-5); +SEBranchPSigma=SEBranchPSigma(indBranchPS); +SEBranchQSigma=imag(rBranchPQ)*sigma; +indBranchPS=find(SEBranchPSigma>1e-5); +SEBranchQSigma=SEBranchQSigma(abs(SEBranchQSigma)>1e-5); SEBranchP=real((cmpSEV(this.lineI)-cmpSEV(this.lineJ)).*conj(cmpSEBranchI)); SEBranchQ=imag((cmpSEV(this.lineI)-cmpSEV(this.lineJ)).*conj(cmpSEBranchI)); -Objective=Objective+sum((SEBranchP-this.mBranchP).^2/this.sigma^2); -Objective=Objective+sum((SEBranchQ-this.mBranchQ).^2/this.sigma^2); +Objective=Objective+sum((SEBranchP-this.mBranchP).^2./(SEBranchPSigma.^2)); +% Objective=Objective+sum((SEBranchQ-this.mBranchQ).^2./(SEBranchQSigma.^2)); + + %% 变压器 newwordParameter=this.newwordParameter; cmpY=this.cmpY; transP=TransPower( newwordParameter,cmpY,SEVolt,SEVAngel ); -Objective=Objective+sum((transP-this.mTransP).^2/this.sigma^2); +Objective=Objective+sum((transP-this.mTransP).^2/(this.sigma^2)); transQ=TransReactivePower( newwordParameter,cmpY,SEVolt,SEVAngel ); -Objective=Objective+sum((transQ-this.mTransQ).^2/this.sigma^2); +Objective=Objective+sum((transQ-this.mTransQ).^2/(this.sigma^2)); %% 0注入节点 PQ=diag(cmpSEV)*conj(this.cmpY*cmpSEV); %% 发电机注入功率 -Objective=Objective+(this.mPG(this.onlyPG)-real(PQ(this.onlyPG)))/this.sigmas.onlyPG)-real(PQ(this.onlyPG))); -% Objective=Objective+(this.mQG(this.onlyQG)-imag(PQ(this.onlyQG)))'*(1./this.sigma^2*eye(length(this.mQG(this.onlyQG))))*(this.mQG(this.onlyQG)-imag(PQ(this.onlyQG))); +% Objective=Objective+sum(((this.mPG(this.onlyPG)-real(PQ(this.onlyPG)))/this.sigma).^2); +% Objective=Objective+sum(((this.mQG(this.onlyQG)-imag(PQ(this.onlyQG)))./this.sigma).^2); end diff --git a/@SEOpti/fun.m b/@SEOpti/fun.m index ebe5058..9422f88 100644 --- a/@SEOpti/fun.m +++ b/@SEOpti/fun.m @@ -4,26 +4,39 @@ function [ Objective ] = fun(this, x ) Objective=0; SEVolt=x(1:this.Busnum); SEVAngel=x(this.Busnum+1:2*this.Busnum); -Objective=sum(((SEVolt-this.mVolt)./this.sigma).^2);%电压 +rVolt=this.rVolt; +sigma=this.sigma; +VoltSigma=rVolt*sigma; +Objective=sum(((SEVolt-this.mVolt)./VoltSigma).^2);%电压 %% 支路电流 cmpSEV=SEVolt.*exp(1j*SEVAngel); %复数电压 cmpSEBranchI=(cmpSEV(this.lineI)-cmpSEV(this.lineJ))./(this.lineR+1j*this.lineX);%复数支路电流 SEBranchI=abs(cmpSEBranchI);% 支路电流幅值 % (SEBranchI-this.mBranchI).^2./(this.sigma^2) -Objective=Objective+sum((SEBranchI-this.mBranchI).^2./(this.sigma^2));%%电流 +% Objective=Objective+sum((SEBranchI-this.mBranchI).^2./(this.sigma^2));%%电流 %% 线路功率 +rVAngel=this.rVAngel; +cmpV=rVolt.*exp(1j*rVAngel); +rBranchPQ=(cmpV(this.lineI)-cmpV(this.lineJ)).*conj( (cmpV(this.lineI)-cmpV(this.lineJ))./(this.lineR+1j*this.lineX) ); +SEBranchPSigma=real(rBranchPQ)*sigma; +indBranchPS=find(real(rBranchPQ)>1e-5);%不考虑太小的数据,以免权重太大。例如除以0 +SEBranchPSigma=SEBranchPSigma(indBranchPS); +SEBranchQSigma=imag(rBranchPQ)*sigma; +indBranchQS=find(imag(rBranchPQ)>1e-5);%不考虑太小的数据,以免权重太大。例如除以0 +SEBranchQSigma=SEBranchQSigma(indBranchQS); SEBranchP=real((cmpSEV(this.lineI)-cmpSEV(this.lineJ)).*conj(cmpSEBranchI)); +SEBranchP=SEBranchP(indBranchPS); SEBranchQ=imag((cmpSEV(this.lineI)-cmpSEV(this.lineJ)).*conj(cmpSEBranchI)); -% (SEBranchP-this.mBranchP).^2/(this.sigma^2)./(this.mBranchP.^2); -Objective=Objective+sum((SEBranchP-this.mBranchP).^2/(this.sigma^2)); -Objective=Objective+sum((SEBranchQ-this.mBranchQ).^2/(this.sigma^2)); +SEBranchQ=SEBranchQ(indBranchQS); +Objective=Objective+sum((SEBranchP-this.mBranchP).^2./(SEBranchPSigma.^2)); +Objective=Objective+sum((SEBranchQ-this.mBranchQ).^2./(SEBranchQSigma.^2)); %% 变压器 newwordParameter=this.newwordParameter; cmpY=this.cmpY; transP=TransPower( newwordParameter,cmpY,SEVolt,SEVAngel ); Objective=Objective+sum((transP-this.mTransP).^2/(this.sigma^2)); transQ=TransReactivePower( newwordParameter,cmpY,SEVolt,SEVAngel ); -Objective=Objective+sum((transQ-this.mTransQ).^2/(this.sigma^2)); +% Objective=Objective+sum((transQ-this.mTransQ).^2/(this.sigma^2)); %% 0注入节点 PQ=diag(cmpSEV)*conj(this.cmpY*cmpSEV); %% 发电机注入功率 diff --git a/@SEOpti/init.m b/@SEOpti/init.m index 6b7f752..87a38d1 100644 --- a/@SEOpti/init.m +++ b/@SEOpti/init.m @@ -1,4 +1,4 @@ -function [ this ] = init(this,Busnum, mVolt,sigma,newwordParameter,zerosInjectionIndex,cmpY,onlyPG,onlyQG,mPG,mQG,Balance,mBranchI,mBranchP,mBranchQ,mTransP,mTransQ) +function [ this ] = init(this,Busnum, mVolt,sigma,newwordParameter,zerosInjectionIndex,cmpY,onlyPG,onlyQG,mPG,mQG,Balance,mBranchI,mBranchP,mBranchQ,mTransP,mTransQ,BalanceVolt,rVolt,rVAngel) %INIT Summary of this function goes here % Detailed explanation goes here this.mVolt=mVolt; @@ -22,7 +22,10 @@ function [ this ] = init(this,Busnum, mVolt,sigma,newwordParameter,zerosInjectio this.mTransP=mTransP; this.mTransQ=mTransQ; this.newwordParameter=newwordParameter; - + this.BalanceVolt=BalanceVolt; + %真实值 + this.rVolt=rVolt; + this.rVAngel=rVAngel; % this.Y=Y; % this.Yangle=Yangle; % this.r=r; diff --git a/mBranchI.mat b/mBranchI.mat index 502780e6b452ae7aaa26dfab474cd7d66db37d60..561bacf2486bf09e864b1d9aabc0f2ba61bb4ed2 100644 GIT binary patch delta 179 zcmeyy^o?nP9gmrnfuWV5iGq=Vq4C5(?THEO6KhiIkFQ;ms{Oup`}u%hn{tZ382`5I zx_(bKsU}bVvCaF6iVx1!hBy0p#rB!ME#)g`RgbsLKR>fdr+ofh;os}~_pM*M^A!6{^nYdisCyxBTZa^D~31kNi8lehC=UGJsuNf&iD ozqxOc`l2p!@3lIyeYvNYXRfqPxq^{_q4C5(?THEO6KhiI+q{qeyK(#eTlJDJ^I!h_bmrfy z#cJP{v7Ucp-J-j@)y`(acJu8qj|HnP#}#Z#{qbyN`dfL=x=VYDUu~>?5!_j8z1r;d z^Y<(EswHeXyYI`>ec$)A9)AYL-*|QD>7;6g{yZN;0PLGz&;S4c diff --git a/mBranchP.mat b/mBranchP.mat index 1ba197c6eed986cdc51d29c1975ca1a0e694e8a0..1cd943b2f983c48295ae1b60a180ee5e330e36b0 100644 GIT binary patch delta 170 zcmaFL^pt6W9gmrnfuWV5iGq=Vq4C5(?THEO6Ki7Y3okm|yzpYd(r=>D^&jgd&Pv{& z{kQS`zL}Tb-rTur*8V^96YN#jRp0Dl_Wia!BJFC(&*wM2wY>c~?!+d^D}H*MsFKr%`)eL_=+9v=2fX7*H delta 170 zcmaFL^pt6W9gl&Psil>Pxq^{_q4C5(?THEO6Ki7YH=h0;e(lbo({J9JpXDyfpOZX4 ze%Jqp2V*b4U99`9I4#R>%lZES+YYYZcyeFf*?-=GY?a??v`=mSl0N0u>tkyFAH=-R z$TzC|x2{0Z{@`u+sIR(-3SyP4IykEBO!wW09 zjQZ^ztJ9CJQTV-f;q_(v;^(E$U-tfkaY5!n+tYWBervk@`pxc`WG(Y)GjCs?`Lazl pCTQFL|LK?b*VP~VCbnI7yXcx6-)-VkV)stazj8v3K}X5m4*+jAT^;}c delta 179 zcmey$^p$CX9gl&Psil>Pxq^{_q4C5(?THEO6KhiHcXWTgl6(9o`-xo^UMnBCr7@n?T*eOA-nmsNaczq>cT zQP!J3r~Fg?gxZ1!Rp&FyKfUGdnRjn)Ty)BM%lrQtlH8jAHhsC1ADH*-TYUZg!t?hR nr^sLW))kl*I$^Jv&voO!wf&Rcii*CI+um)uiJc+hx|<>ZJiuIo diff --git a/mPG.mat b/mPG.mat index c08e25ae860c207f3373f6f7670a97af53100ee2..020a2afd656fa87b13c0497458bc5450eed54588 100644 GIT binary patch delta 55 zcmX@Yc!Y6+9gmrnfuWV5iGq=Vq4C5(?THEO6KnLOeoXrB_w{MG{7vqaOV6AZ+M9ab MA(??8U>)ZL02{p(f&c&j delta 55 zcmX@Yc!Y6+9gl&Psil>Pxq^{_q4C5(?THEO6KnLOcKp4%DKza`|J$jn*OzOG-r17A Mz>tBV;{;bM03xUrjsO4v diff --git a/mQG.mat b/mQG.mat index 6df8a65e80ae33c1e590890109d445da29deee9c..f1f11f6c6eb96e8b3e2f624646e4c45efeda5f9b 100644 GIT binary patch delta 79 zcmcb`c$0C09gmrnfuWV5iGq=Vq4C5(?THEO6Kgma11Gks>;C#{GyA=F-p%%n=Tf5V key_1*meQC|yG}^3P(V+w)GT@b;rYoDuQ(VP#1w@D0ACLq=Kufz delta 80 zcmcb~c#Cm@9gl&Psil>Pxq^{_q4C5(?THEO6KgmagC@4B>wRua>fT+Rck_L>{qlcz me$KIFmb!2=X!je}6LmXwe4JLhz3z(4;Q)6I28QP?BIg12aU(Pxq^{_q4C5(?THEO6Kjm+0)NXDRdrb3u<^dR{dUE=J-svj Q9(rHR!koh(Cc)bb0BOn=fdBvi diff --git a/mTransQ.mat b/mTransQ.mat index f7a3b7a40a434259d4b66d4edf066ecb9f8eeaae..3e328057a4973d0f78a4680214ff8c404f2d61a3 100644 GIT binary patch delta 60 zcmX@dc#d&`9gmrnfuWV5iGq=Vq4C5(?THEO6Kjm+uT}fk2}T+u@4qC}`26x~wKwr* R)wBAhyngkU!9qbm0RVQQ7f%2H delta 60 zcmX@dc#d&`9gl&Psil>Pxq^{_q4C5(?THEO6Kjm+HNS`1mCm$GuJ`kpyX)W5e}B`$ R@9WK}FTYyHU>_m45&(I37pMRL diff --git a/mVolt.mat b/mVolt.mat index 492b7862a05315f429d6608f975952be3719c19f..63417b73c405dc439a11d0c6a43a3a75c7ab2903 100644 GIT binary patch delta 158 zcmX@ZbdG6)9gmrnfuWV5iGq=Vq4C5(?THEO6Kgmak4tUp`x+-cGa_v*j* zR{x(jy?#r1`aks(?>=1I_hx1J#GQB8f7@?Pxq^{_q4C5(?THEO6Kgmak4$X!t-3KMfB%iZd;Q;nzpq=d z{q?Ve|Bvonsy1%2zv0bl5&KW-=fS$&R;!n}-edoLQ#V!kRjb`C^>aUdm^MFN{P+0w z)ALpHer4!?t6Q0xaGL!$>+RFwb}!jh*=I~D`}dvwH_!i5_k+*XY1e-5); mBranchP=rBranchP.*(normrnd(0,sigma,length(rBranchP),1)+1);%支路功率量测量 rBranchQ=imag((cmpV(lineI)-cmpV(lineJ)).*conj(cmpBranchI)); +rBranchQ=rBranchQ(abs(rBranchQ)>1e-5); mBranchQ=rBranchQ.*(normrnd(0,sigma,length(rBranchQ),1)+1);%支路功率量测量 %% 注入功率 rPD=PD; @@ -82,15 +85,15 @@ seOpti=SEOpti(); % save('mTransP','mTransP'); % save('mTransQ','mTransQ'); %% load -load('mVolt'); -load('mPG'); -load('mQG'); -load('mBranchI'); -load('mBranchP'); -load('mBranchQ'); -load('mTransP'); -load('mTransQ'); -seOpti=seOpti.init(Busnum, mVolt,sigma,newwordParameter,zerosInjectionIndex,cmpY,onlyPG,onlyQG,mPG,mQG,Balance,mBranchI,mBranchP,mBranchQ,mTransP,mTransQ); +% load('mVolt'); +% load('mPG'); +% load('mQG'); +% load('mBranchI'); +% load('mBranchP'); +% load('mBranchQ'); +% load('mTransP'); +% load('mTransQ'); +seOpti=seOpti.init(Busnum, mVolt,sigma,newwordParameter,zerosInjectionIndex,cmpY,onlyPG,onlyQG,mPG,mQG,Balance,mBranchI,mBranchP,mBranchQ,mTransP,mTransQ,BalanceVolt,rVolt,rVAngel); opts = optiset('solver','ipopt'); opts.maxiter=85500; opts.maxtime=3000; @@ -100,7 +103,7 @@ opts.tolrfun=1e-4; opts.tolafun=1e-4; opts.warnings='all'; opts.display='off'; -x0=[ones(Busnum,1);-0.5*ones(Busnum,1)]; +x0=[ones(Busnum,1);-0.2*ones(Busnum,1)]; % x0=[rVolt;rVAngel]; [~,seOpti]=seOpti.equ(x0); nlrhs=seOpti.nlrhs();