一、简介
参照何仰赞电力系统分析（下）的牛拉法程序流程框图，将程序分为节点导纳矩阵的形成、判敛、不平衡量计算、雅可比矩阵计算、节点电压修正、计算功率等几个部分。
二、源代码
clear
%先求节点导纳矩阵
N1=30; %节点数
L1=43; %支路数
%N1=5;L1=7;
B1=xlsread('XIANLU2.xlsx','B2:F44');
%B1=xlsread('XIANLU.xlsx','B2:F44');
G=zeros(N1);
B=zeros(N1);
%每行存储一条支路
%第一列存储支路的一个端点I,变压器所在端加个负号，励磁支路为0
%第二列存储另一个端点J，变压器所在端加个负号，励磁支路为0
%第三列存储支路的电阻R
%第四列存储支路电抗X
%第五列存储支路对地电纳或变压器变比k,k在首端为正
for m1=1:L1
I=B1(m1,1);J=B1(m1,2);R=B1(m1,3);X=B1(m1,4);k=B1(m1,5);
if I*J==0 %励磁支路？
if I==0
G(J,J)=G(J,J)+R;
B(J,J)=B(J,J)+X;
end
if J==0
G(I,I)=G(I,I)+R;
B(I,I)=B(I,I)+X;
end
else
if I*J>0 %线路支路？
B(I,I)= B(I,I)+ k;
B(J,J)= B(J,J)+ k;
k=1;
else
if I<0
t=I;
I=J;
J=t;
end
J=abs(J);
if k<0
k=-1/k;
end
end
G(I,J)=-(1/k)*R/(R^2+X^2);
G(J,I)=G(I,J);
B(I,J)=(1/k)*X/(R^2+X^2);
B(J,I)=B(I,J);
G(I,I)=G(I,I)+R/(R^2+X^2);
B(I,I)=B(I,I)-X/(R^2+X^2);
G(J,J)=G(J,J)+1/(k^2)*R/(R^2+X^2);
B(J,J)=B(J,J)-1/(k^2)*X/(R^2+X^2);
end
end
Y=G+B*1i;
%求差.
b=xlsread('JIEDIAN2.xlsx','C2:F31');
%b=xlsread('JIEDIAN.xlsx','C2:F31');
precision=1; %误差
t=0; %迭代次数
pq=24; %pq节点数
%pq=4;
U=b(pq+2:N1,1);
%开始牛拉法
while precision>0.00001
P1=zeros(pq,1);
Q=zeros(pq,1);
P2=zeros(N1-pq-1,1);
U2=zeros(N1-pq-1,1);
deltpqu=zeros(2*(N1-1),1);
deltu2=zeros(N1-pq-1,1);
for i=1:pq
for j=1:N1
P1(i,1)=P1(i,1)+b(i+1,1)*(G(i+1,j)*b(j,1)-B(i+1,j)*b(j,2))+b(i+1,2)*(G(i+1,j)*b(j,2)+B(i+1,j)*b(j,1));
Q(i,1) =Q(i,1) +b(i+1,2)*(G(i+1,j)*b(j,1)-B(i+1,j)*b(j,2))-b(i+1,1)*(G(i+1,j)*b(j,2)+B(i+1,j)*b(j,1));
end
end
for i=1:N1-pq-1
for j=1:N1
P2(i,1)=P2(i,1)+b(pq+i+1,1)*(G(pq+i+1,j)*b(j,1)-B(pq+i+1,j)*b(j,2))+b(pq+i+1,2)*(G(pq+i+1,j)*b(j,2)+B(pq+i+1,j)*b(j,1));
end
U2(i,1)=b(pq+i+1,1)^2+b(pq+i+1,2)^2;
end
deltp1=b(2:pq+1,3)-P1(1:pq,1);
deltq =b(2:pq+1,4)-Q(1:pq,1);
deltp2=b(pq+2:N1,3)-P2(1:N1-pq-1,1);
for i=1:N1-pq-1
deltu2(i,1)=U(i,1)^2-U2(i,1);
end
三、文件
四、备注
版本：2014a
