diff --git a/OPF.m b/OPF.m index 7505c94..0a17982 100644 --- a/OPF.m +++ b/OPF.m @@ -5,10 +5,8 @@ arrayA=zeros(21,10); sumCaseA_SE=0; sumCaseB_SE=0; sumCaseC_SE=0; -badDataResult=zeros(10,33); -badDataLocation=zeros(34,10); -nodeMaxDVolt=zeros(33,1); -nodeMaxDVAngle=zeros(33,1); +VoltAAE=0; +VAngleAAE=0; for badDataNode=1:1 loopN=1; maxDVolt=0; @@ -24,9 +22,16 @@ for badDataNode=1:1 Loadi=find(PD~=0); PD0=sparse(Busnum,1); QD0=sparse(Busnum,1); - PD0(Loadi)=RealPD(Loadi).*(1+normrnd(0,sigma,length(Loadi),1)); - QD0(Loadi)=RealQD(Loadi).*(1+normrnd(0,sigma,length(Loadi),1)); - mVolt=rVolt.*(1+normrnd(0,sigma,length(rVolt),1))'; + %先试一下均匀分布 + PD0(Loadi)=RealPD(Loadi).*(1+unifrnd(-RealPD(Loadi)*sigma,RealPD(Loadi)*sigma,length(Loadi),1)); + QD0(Loadi)=RealQD(Loadi).*(1+unifrnd(-RealQD(Loadi)*sigma,RealQD(Loadi)*sigma,length(Loadi),1)); + mVolt=rVolt.*(1+unifrnd(-rVolt*sigma,rVolt*sigma,1,length(rVolt))); + + + +% PD0(Loadi)=RealPD(Loadi).*(1+normrnd(0,sigma,length(Loadi),1)); +% QD0(Loadi)=RealQD(Loadi).*(1+normrnd(0,sigma,length(Loadi),1)); +% mVolt=rVolt.*(1+normrnd(0,sigma,length(rVolt),1))'; %加载保存的变量 % PD0=load('PD0'); @@ -40,39 +45,21 @@ for badDataNode=1:1 %% 画Case A的图 % figure('Color',[1 1 1]); [Busnum,Loadi,Volt,PD,QD,rVolt,UAngel,RealPD,RealQD,rUAngel,Vbi,PDbi,QDbi,plotGapA,isConverge]=subOPF([],PD0,QD0,mVolt,sigma);%全部有 - % badDataResult(I,badDataNode)=sum(Vbi); - % badDataLocation(1:33,I)=Vbi; - % badDataLocation(34,I)=sum(abs((rVolt-Volt)./rVolt./length(rVolt)))+sum(abs( (UAngel(2:33)-rUAngel(2:33))./rUAngel(2:33)./length(rUAngel(2:33)))); -% break; if isConverge==0 continue; end - % 加最大偏差 - maxDVolt_=max(abs((rVolt-Volt))); - if maxDVolt_>maxDVolt - maxDVolt=maxDVolt_; - end - maxDVAngle_=max(abs((UAngel(2:33)-rUAngel(2:33)))); - if maxDVAngle_>maxDVAngle - maxDVAngle=maxDVAngle_; - end - nodeMaxDVolt_t=abs((rVolt-Volt))'; - nodeMaxDVAngle_t=abs((UAngel-rUAngel))'; - nodeMaxDVolt_t([1:18,23:33])=0; - nodeMaxDVAngle_t([1:18,23:33])=0; - nodeMaxDVolt(nodeMaxDVolt=1 -% nodeMaxDVolt(badDataNode)=maxDVolt; -% nodeMaxDVAngle(badDataNode)=maxDVAngle; + VoltAAE=VoltAAE+sum(abs((Volt-rVolt)./rVolt)); + VAngleAAE=VAngleAAE+sum(abs((UAngel(2:end)-rUAngel(2:end))./rUAngel(2:end))); + loopN=loopN+1; + if loopN>=500 break; end - - loopN=loopN+1; + end end - +VoltAAE=VoltAAE/(loopN*length(Volt))*100; +VAngleAAE=VAngleAAE/(loopN*length(UAngel(2:end)))*100; %% PLOTING % % end