1 主要内容

该部分资料是牛拉法和PQ分解法两种潮流计算方法的代码和对应的资料，程序针对5节点系统，也可以自行修改节点和线路参数改成其他节点系统，程序通用性较强，注释清晰，方便学习。

2 部分代码

n=5;                      %输入节点数
nl=5;                     %输入支路数
isb=1;                    %输入平衡母线节点号  
pr=0.00001;               %输入误差精度
B1=[1 2 0.03i 0 1.05 0;2 3 0.08+0.3i 0.5i 1 0;2 4 0.1+0.35i 0 1 0;3 4 0.04+0.25i 0.5i 1 0;3 5 0.015i 0 1.05 1];           %输入由支路参数形成的矩阵
B2=[0 0 1.05 1.05 0 1;0 3.7+1.3i 1.05 0 0 2;0 2+1i 1.05 0 0 2;0 1.6+0.8i 1.05 0 0 2;5 0 1.05 1.05 0 3];                      %输入由节点参数形成的矩阵
X=[1 0;2 0;3 0;4 0;5 0];     %输入由节点号及其对地阻抗形成的矩阵
Y=zeros(n);e=zeros(1,n);f=zeros(1,n);V=zeros(1,n);O=zeros(1,n);S1=zeros(nl);O1=zeros(1,n);
%形成导纳矩阵
for i=1:n
    if X(i,2)~=0                    %判断节点对地是否存在阻抗
        p=X(i,1);                   %如果该节点对地存在阻抗，则将该节点编号赋值给p
        Y(p,p)=Y(p,p)+1./X(i,2);      %计算对地阻抗对于节点自阻抗的影响
    end
end
for i=1:nl
    if B1(i,6)==0                   %判断该支路参数是否已折算到变压器高压侧
        p=B1(i,1);q=B1(i,2);         %如果该支路参数是否已折算到变压器高压侧，则将B1矩阵第一列元素赋值给p，将B2矩阵第二列元素赋值给q 
    else
        p=B1(i,2);q=B1(i,1);        %如果该支路参数是否已折算到变压器高压侧，则将B1矩阵第一列元素赋值给q，将B2矩阵第二列元素赋值给p
    end
    Y(p,q)=Y(p,q)-1./(B1(i,3)*B1(i,5));     %计算节点间互导纳
    Y(q,p)=Y(p,q);                          %由于节点导纳矩阵为对称矩阵，因此Y(q,p)=Y(p,q)
    Y(q,q)=Y(q,q)+1./(B1(i,3)*B1(i,5)^2)+B1(i,4)./2;   %计算q节点自导纳
    Y(p,p)=Y(p,p)+1./B1(i,3)+B1(i,4)./2;             %计算p节点自导纳
end

%形成新的导纳矩阵(将PQ,PV和平衡节点进行分类)
node_sign=zeros(1,n);                    %记录新的导纳矩阵节点编号
PQ_count=0;                           %初始化PQ节点计数器
PV_count=0;                           %初始化PV节点计数器
for i=1:n
    if B2(i,6)==1                        %寻找平衡节点
        node_sign(1,n)=i;                %记录该PQ节点编号并将其放置到新的导纳矩阵节点编号矩阵最后一列
    end
    if B2(i,6)==2                        %寻找PQ节点
        PQ_count=PQ_count+1;           %PQ节点个数加1
        node_sign(1,PQ_count)=i;         %记录该PQ节点编号
    end
end