Matlab潮流计算程序改进.docx
《Matlab潮流计算程序改进.docx》由会员分享，可在线阅读，更多相关《Matlab潮流计算程序改进.docx（16页珍藏版）》请在冰点文库上搜索。
Matlab潮流计算程序改进
Matlab 潮流计算程序（改进）
MATlab 潮流计算
花了很多时间做的
感觉很有意义
从早上7点到现在调试的都没吃饭
饿了一天了
有点成就感；
但是ieee118点的我做的电压前边117个收敛的还可，但为什么最后一个相差很大，而且相角度也有很大的误差；
但是同样的程序对于5ieee,14ieee,30ieee都运行的不错。
这还得进一步研究一下。
得吃饭去了，饿的慌，眼睛都冒星星了.............
估计赵老师还是不会满意，所以得还得改改
%==========================================================================
%==========================================================================
%==========================================================================
%潮流计算MATLAB粗略程序 C.zhou 2009.3.27
%==========================================================================
%==========================================================================
%==========================================================================
%creatanew_data
t=0;
s=0;
r=0;
w=0;
number=input（'Howmanynodearethere='）;
%ConvertPqtoanewarray
forii=1:
number
ifdata（ii,4）==1
t=t+1;
forjj=1:
14
new_data1（t,jj）=data（ii,jj）;
end;
a（1,t）=ii;
s=s+1; %recordthenumberofthePQnode
end;
end;
%Convertpvtoanewarray
forii=1:
number
ifdata（ii,4）==2
t=t+1;
forjj=1:
14
new_data1（t,jj）=data（ii,jj）;
end;
a（1,t）=ii;
r=r+1; %recordthenumberofthePVnode
end;
end;
%Convertset_vtoanewarray
forii=1:
number
ifdata（ii,4）==3
t=t+1;
forjj=1:
14
new_data1（t,jj）=data（ii,jj）;
end;
a（1,t）=ii;
w=w+1;
end;
end;
%creatanew_data2
[x,y]=size（data2）
forii=1:
x
forjj=1:
2
formm=1:
number
ifdata2（ii,jj）==a（1,mm）
new_data2（ii,jj）=mm;
end;
end;
end;
end;
forii=1:
x
forjj=3:
14
new_data2（ii,jj）=data2（ii,jj）;
end;
end;
%creataY
Y=zeros（number,number）;
YY=zeros（number,number）;
yy=zeros（number,number）;
forii=1:
x
%forjj=1:
14
iii=new_data2（ii,1）;
jjj=new_data2（ii,2）;
ifnew_data2（ii,5）==2
sub=new_data2（ii,6）./（new_data2（ii,7）.*new_data2（ii,7）+new_data2（ii,6）.*new_data2（ii,6））-new_data2（ii,7）./（new_data2（ii,7）.*new_data2（ii,7）+new_data2（ii,6）.*new_data2（ii,6））*i;
Y（iii,jjj）=-sub./new_data2（ii,14）;
YY（iii,jjj）=sub./new_data2（ii,14）;
Y（jjj,iii）=-sub/new_data2（ii,14）;
YY（jjj,iii）=sub./new_data2（ii,14）;
yy（iii,jjj）=（1.-new_data2（ii,14））./（new_data2（ii,14）.*new_data2（ii,14））.*sub;
yy（jjj,iii）=（new_data2（ii,14）-1）./（new_data2（ii,14））.*sub;
else
Y（iii,jjj）=-new_data2（ii,6）./（new_data2（ii,7）.*new_data2（ii,7）+new_data2（ii,6）.*new_data2（ii,6））+new_data2（ii,7）./（new_data2（ii,7）.*new_data2（ii,7）+new_data2（ii,6）.*new_data2（ii,6））*i;
YY（iii,jjj）=new_data2（ii,6）./（new_data2（ii,7）.*new_data2（ii,7）+new_data2（ii,6）.*new_data2（ii,6））-new_data2（ii,7）./（new_data2（ii,7）.*new_data2（ii,7）+new_data2（ii,6）.*new_data2（ii,6））*i;
Y（jjj,iii）=-new_data2（ii,6）./（new_data2（ii,7）.*new_data2（ii,7）+new_data2（ii,6）.*new_data2（ii,6））+new_data2（ii,7）./（new_data2（ii,7）.*new_data2（ii,7）+new_data2（ii,6）.*new_data2（ii,6））*i;
YY（jjj,iii）=new_data2（ii,6）./（new_data2（ii,7）.*new_data2（ii,7）+new_data2（ii,6）.*new_data2（ii,6））-new_data2（ii,7）./（new_data2（ii,7）.*new_data2（ii,7）+new_data2（ii,6）.*new_data2（ii,6））*i;
yy（iii,jjj）=new_data2（ii,8）./2.*i;
yy（jjj,iii）=new_data2（ii,8）./2.*i;
end;
%end;
end;
foriii=1:
number
Y（iii,iii）=0;
end;
%forii=1:
x
% forjj=1:
14
for iii=1:
number
forjj=1:
number
%ifiii~=jj
Y（iii,iii）=Y（iii,iii）+YY（iii,jj）+yy（iii,jj）;
%end;
end;
end;
%creatB,G
forii=1:
number
forjj=1:
number
G（ii,jj）=real（Y（ii,jj））;
B（ii,jj）=imag（Y（ii,jj））;
end;
end;
%creatInitial_PInitial_QInitial_V
forii=1:
（s+r）
set_P（ii,1）=（new_data1（ii,9）-new_data1（ii,7））./100;
end;
forii=1:
s;
set_Q（ii,1）=（new_data1（ii,10）-new_data1（ii,8））./100;
end;
forii=1:
r
set_V（ii,1）=new_data1（ii+s,12）.*new_data1（ii+s,12）;%trytomodifyforsikeofcorrecting
end;
Initial_p_q_v=[set_P;set_Q;set_V];
disp（Initial_p_q_v）;
%creatInitial_e,Initial_f
forii=1:
number-1
e（ii,1）=1;
f（ii,1）=0.0;%changeftotestusedtobe1.0
end;
e（number,1）=new_data1（number,12）;
f（number,1）=0;
%e（64,1）=0.88;%test118ieee
%f（64,1）=0.39395826829394;
%f（14,1）=0;
%e（10,1）=1.045;
%e（11,1）=1.01;
%e（12,1）=1.07;
%e（13,1）=1.09;
%//////////////////////////////////////////////////////////////////////////
%/////////////////////////////////////////////////////////////////////////
%//////////////////////////////////////////////////////////////////////////
%//////////////////////////////////////////////////////////////////////////
%StartNEWTOWNCALULATION
fortry_time=1:
25
%CreateverynodeconsumePQandU
n=s;
m=r;
forii=1:
（n+m）
sum1=0;
forjj=1:
（n+m+1）
sum1=sum1+e（ii,1）.*（G（ii,jj）.*e（jj,1）-B（ii,jj）.*f（jj,1））+f（ii,1）.*（G（ii,jj）.*f（jj,1）+B（ii,jj）.*e（jj,1））;
end;
p（ii,1）=sum1;
end;
forii=1:
n
sum2=0;
forjj=1:
（n+m+1）
sum2=sum2+f（ii,1）.*（G（ii,jj）.*e（jj,1）-B（ii,jj）.*f（jj,1））-e（ii,1）.*（G（ii,jj）.*f（jj,1）+B（ii,jj）.*e（jj,1））;
end;
q（ii,1）=sum2;
end;
disp（'q='）;
disp（q）;
u=zeros（（n+m）,1）;
forii=（n+1）:
（n+m）
u（ii,1）=e（ii,1）.*e（ii,1）+f（ii,1）.*f（ii,1）;
end;
forii=n+1:
（n+m）
extra_u（（ii-n）,1）=u（ii,1）;
end;
disp（'extra_u='）;
disp（extra_u）;
sum=[p;q;extra_u];
disp（sum）
disp（s）;
disp（p）;
%creatJacobian
disp（n）;
disp（m）;
forii=1:
（n+m）
forjj=1:
（n+m）
if（ii~=jj）
PF（ii,jj）=B（ii,jj）.*e（ii,1）-G（ii,jj）.*f（ii,1）;
PE（ii,jj）=-G（ii,jj）.*e（ii,1）-B（ii,jj）.*f（ii,1）;
else
ss=0;
qq=0;
fornum=1:
（n+m+1）
ss=ss+G（ii,num）.*f（num,1）+B（ii,num）.*e（num,1）;
qq=qq+G（ii,num）.*e（num,1）-B（ii,num）.*f（num,1）;
end;
PF（ii,jj）=-ss+B（ii,jj）.*e（ii,1）-G（ii,jj）.*f（ii,1）;%TEST+1
PE（ii,jj）=-qq-G（ii,jj）.*e（ii,1）-B（ii,jj）.*f（ii,1）;%TEST+1
end;
end;
end;
copy=3.14159;
disp（'================copy================'）
forii=1:
n
forjj=1:
m+n
if（ii~=jj）
QE（ii,jj）=B（ii,jj）.*e（ii,1）-G（ii,jj）.*f（ii,1）;%TEST+1
QF（ii,jj）=G（ii,jj）.*e（ii,1）+B（ii,jj）.*f（ii,1）;%TEST+1
else
ss=0;
qq=0;
fornum=1:
（n+m+1）
ss=ss+G（ii,num）.*f（num,1）+B（ii,num）.*e（num,1）;
qq=qq+G（ii,num）.*e（num,1）-B（ii,num）.*f（num,1）;
end;
QF（ii,jj）=-qq+G（ii,jj）.*e（ii,1）+B（ii,jj）.*f（ii,1）;%TEST+1
QE（ii,jj）=ss+B（ii,jj）.*e（ii,1）-G（ii,jj）.*f（ii,1）;%TEST+1
end;
end;
end;
%disp（'QF'）;
%disp（QF）;
%disp（'QE'）;
%disp（QE）;
UE=zeros（（n+m）,（n+m））;
UF=zeros（（n+m）,（n+m））;
forii=n+1:
n+m
forjj=1:
（n+m）
if（ii~=jj）
UE（ii,jj）=0;
UF（ii,jj）=0;
else
ss=0;
qq=0;
fornum=1:
（n+m+1）
ss=ss+G（ii,num）.*f（num,1）+B（ii,num）.*e（num,1）;
qq=qq+G（ii,num）.*e（num,1）-B（ii,num）.*f（num,1）;
end;
UF（ii,jj）=-2.*f（ii,1）;
UE（ii,jj）=-2.*e（ii,1）;
end;
end;
end;
forii=（n+1）:
（n+m）
forjj=1:
（n+m）
extra_UE（（ii-n）,jj）=UE（ii,jj）;
extra_UF（（ii-n）,jj）=UF（ii,jj）;
end;
end;
%disp（'extra_UE'）;
%disp（extra_UE）;
%disp（'extra_Uf'）;
%disp（extra_UF）;
Jacobian=[PF,PE;QF,QE;extra_UF,extra_UE];
%disp（'Jacobian='）;
%disp（Jacobian）;
%creatsubstractresult
substract_result=Initial_p_q_v-sum;
%disp（'substract_result'）;
%disp（substract_result）;
%calculatedelta_f_e
delta_f_e=-inv（Jacobian）*substract_result;
%disp（delta_f_e）;
forii=1:
number-1;
f（ii,1）=f（ii,1）+delta_f_e（ii,1）;
e（ii,1）=e（ii,1）+delta_f_e（ii+number-1,1）;
end;
ifmax（substract_result）<1e-4
break;
end;
end;
%disp（'substract_result'）;
%disp（substract_result）;
%disp（'e='）;
%disp（e）;
%disp（'f='）;
%disp（f）;
forii=1:
number
uuu（ii,1）=e（ii,1）.*e（ii,1）+f（ii,1）.*f（ii,1）;
U_RESULT（ii,1）=sqrt（uuu（ii,1））;
end;
forii=1:
number
for jj=1:
number
ifii==a（1,jj）
Old_Uresult（ii,1）=U_RESULT（jj,1）
end;
end;
end;
forii=1:
number
Old_Uresult（ii,2）=ii;
end;
%disp（'U_result'）;
%disp（U_RESULT）;
disp（'====================================='）;
disp（'Thelastresultis:
'）
disp（'===========U===================BUS-NO.'）;
disp（'U='）
disp（Old_Uresult）;
%calculatetheangle
PI=3.141592
forii=1:
number
Angle（ii,1）=atan（f（ii,1）./e（ii,1））./PI*180;
end;
forii=1:
number
for jj=1:
number
ifii==a（1,jj）
Old_Angle（ii,1）=Angle（jj,1）;
Old_Angle（ii,2）=ii;
end;
end;
end;
disp（'=======Angle===================BUS-NO.'）;
disp（'Angle='）;
disp（Old_Angle）;
disp（'=====Try-times======================='）
disp（'Try-times='）
disp（try_time）;
%disp（'p===================='）;
%disp（p）;
%for jj=1:
number
%ifa（1,jj）==118
% man=jj
% end;
%end;
%disp（'man========='）;
%disp（man）
sum4=0;
forjj=1:
number
Y_conj（number,jj）=conj（Y（number,jj））;
sum4=sum4+Y_conj（number,jj）.*（e（jj,1）-f（jj,1）*i）;
end;
%sum4=sum4*e（number,1）;
disp（'===============BalancePQ=========BUS-NO'）;
%disp（sum4）;
Blance_Q（1,1）=imag（sum4）*100;
Blance_Q（1,2）=a（1,number）;
Blance_P（1,1）=real（sum4）*100;
Blan
升级会员 