matlab潮流计算程序287401

A=[ 1	2	0.01938	0.05917	0.0528	0	0	0	0	0	1	-360	360;
	1	5	0.05403	0.22304	0.0492	0	0	0	0	0	1	-360	360;
	2	3	0.04699	0.19797	0.0438	0	0	0	0	0	1	-360	360;
	2	4	0.05811	0.17632	0.034	0	0	0	0	0	1	-360	360;
	2	5	0.05695	0.17388	0.0346	0	0	0	0	0	1	-360	360;
	3	4	0.06701	0.17103	0.0128	0	0	0	0	0	1	-360	360;
	4	5	0.01335	0.04211	0	0	0	0	0	0	1	-360	360;
	4	7	0	0.20912	0	0	0	0	0.978	0	1	-360	360;
	4	9	0	0.55618	0	0	0	0	0.969	0	1	-360	360;
	5	6	0	0.25202	0	0	0	0	0.932	0	1	-360	360;
	6	11	0.09498	0.1989	0	0	0	0	0	0	1	-360	360;
	6	12	0.12291	0.25581	0	0	0	0	0	0	1	-360	360;
	6	13	0.06615	0.13027	0	0	0	0	0	0	1	-360	360;
	7	8	0	0.17615	0	0	0	0	0	0	1	-360	360;
	7	9	0	0.11001	0	0	0	0	0	0	1	-360	360;
	9	10	0.03181	0.0845	0	0	0	0	0	0	1	-360	360;
	9	14	0.12711	0.27038	0	0	0	0	0	0	1	-360	360;
	10	11	0.08205	0.19207	0	0	0	0	0	0	1	-360	360;
	12	13	0.22092	0.19988	0	0	0	0	0	0	1	-360	360;
	13	14	0.17093	0.34802	0	0	0	0	0	0	1	-360	360;
];
Bus=[
	1	3	0	0	0	0	1	1.06	0	0	1	1.06	0.94;
	2	2	21.7	12.7	0	0	1	1.045	-4.98	0	1	1.06	0.94;
	3	2	94.2	19	0	0	1	1.01	-12.72	0	1	1.06	0.94;
	4	1	47.8	-3.9	0	0	1	1.019	-10.33	0	1	1.06	0.94;
	5	1	7.6	1.6	0	0	1	1.02	-8.78	0	1	1.06	0.94;
	6	2	11.2	7.5	0	0	1	1.07	-14.22	0	1	1.06	0.94;
	7	1	0	0	0	0	1	1.062	-13.37	0	1	1.06	0.94;
	8	2	0	0	0	0	1	1.09	-13.36	0	1	1.06	0.94;
	9	1	29.5	16.6	0	19	1	1.056	-14.94	0	1	1.06	0.94;
	10	1	9	5.8	0	0	1	1.051	-15.1	0	1	1.06	0.94;
	11	1	3.5	1.8	0	0	1	1.057	-14.79	0	1	1.06	0.94;
	12	1	6.1	1.6	0	0	1	1.055	-15.07	0	1	1.06	0.94;
	13	1	13.5	5.8	0	0	1	1.05	-15.16	0	1	1.06	0.94;
	14	1	14.9	5	0	0	1	1.036	-16.04	0	1	1.06	0.94;
];
n=size（Bus1）;%求系统节点数
m=size（A1）;%确定支路数
%%
%做一些基础工作
%按PQ PV排列的节点信息矩阵
Bus2=sortrows（Bus2）;
list=zeros（n1）;
 for i=1:n;
     a=Bus2（i1）;
     list（i）=a;
 end
 PQn=length（find（Bus2（:2）==1））;%PQ节点的个数
 PVn=length（find（Bus2（:2）==2））;%PV节点的个数
 %求标幺值
 sb=100;
 for i=1:n
     Bus2（i3）=Bus2（i3）/sb;
     Bus2（i4）=Bus2（i4）/sb;
 end
%%
%节点导纳矩阵的求解
Y=zeros（n）;%形成一个n阶的零矩阵
M=zeros（n）;%形成一个n阶的零矩阵用于储存自导的临时结果
K=zeros（n）;%形成一个n阶的零矩阵用于储存自导的临时结果
for i=1:m;
    Y（A（i1）A（i2））=（A（i3）+j*A（i4））^（-1）;
    a=A（i1）;
    b=A（i2）;
    r=A（i3）;
    x=A（i4）;
    Y（ab）=-1/（r+j*x）;
    Y（ba）=Y（ab）;
end%求互导
for i=1:m;
  a=A（i1）;
  M（aa）=j*A（i5）/2+M（aa）;
end%扫描第一列求对地导纳
for i=1:m;
  b=A（i2）;
  K（bb）=A（b5）*j/2+K（bb）;
end%扫描第二列求对地导纳
for i=1:n;
Y（ii）=-sum（Y（i:））+M（ii）+K（ii）;
end%求每个节点的自导
G=real（Y）;%电导矩阵
B=imag（Y）;%电纳矩阵
G=G（list:）;%按照PQ、PV顺序的电导矩阵
B=B（list:）;%按照PQ、PV顺序的电纳矩阵
%%
%迭代开始
COUNT=0;
Tmax=100;%迭代次数和最大迭代次数。
V=10^（-6）;%设置精度。
%设电压初始值
e=ones（n1）;
f=zeros（n1）;
%求不平衡量dp、dq、dv
Pi=zeros（n1）;
Qi=zeros（n1）;
Vi=zeros（n1）;
for i=1:Tmax;
for i=1:n;
    Pi（i1）=Bus2（i3