14节点辐射潮流计算程序

clear all 
clc
N=16;%输入系统节点数
%输入系统数据（支路、节点）矩阵A
%第一列是支路号
%第二列是支路开关状态；用K=1或0表示分段开关闭合或联络开关断开
%第三列支路首端节点号
%第四列支路末端节点号
%第五列支路电阻（标幺值）
%第六列支路电抗（标幺值）
%第七列支路末端节点负荷有功（标幺值）
%第八列支路末端节点负荷无功（标幺值）
%第九列支路末端节点连接的电容器容量（标幺值）
A=[ 1   1   1   2   0.075   0.10   0.02   0.016     0
    2   1   1   3   0.11    0.11   0.04   0.027     0
    3   1   1   4   0.11    0.11   0.01   0.009     0
    4   1   2   5   0.09    0.18   0.02   0.008     0.0025
    5   1   2   6   0.08    0.11   0.03   0.015     0.002
    6   1   5   7   0.04    0.04   0.015  0.012     0
    7   1   3   8   0.08    0.11   0.05   0.03      0.012
    8   1   8   9   0.08    0.11   0.045  0.02      0.037
    9   1   8  10   0.11    0.11   0.006  0.001     0.006
   10   1   3  11   0.11    0.11   0.01   0.009     0.018
   11   1   4  12   0.09    0.12   0.01   0.007     0
   12   1   4  13   0.08    0.11   0.01   0.009     0 
   13   1  13  14   0.04    0.04   0.021  0.01      0.018
   14   0   6  10   0.04    0.04   0.006  0.001     0.006
   15   0  11  12   0.04    0.04   0.001  0.007     0
   16   0   7  14   0.09    0.12   0.021  0.01      0.018     ];
[mAnA]=size（A）;
a=0;b=0;
for i=1:mA
   K=A（i2）;
    if K==1
        a=a+1;
        C（a1）=A（i1）;
    else if K==0
            b=b+1;
            D（b1）=A（i1）;
        end 
    end
end      
C%C为闭合开关投入的支路
D%D为联络开关打开的支路
[mCnC]=size（C）
B=zeros（mCnA）;
for i=1:mC
    B（i:）=A（i:）;
end
B
Vbus=ones（N1）;%给定电压初值
Vbus（11）=1.0;
Vbus2=Vbus;%迭代后的电压矩阵
[mBnB]=size（B）
k=0;%迭代次数
Ploss=zeros（mB1）;%存支路的有功损耗
Qloss=zeros（mB1）;%支路无功损耗
I=zeros（mB1）;%支路电流
F=zeros（N1）;
%%本段程序将支路重新排序   s1为排好序的支路矩阵
TempB=B;
n=1;
s2=[];    
while ~isempty（TempB）
    [snB]=size（TempB）;%s为支路数  
    m=1;
    while s>0
       i=find（TempB（:3）==TempB（s4））;%末端节点是否为其他支路首端节点
        if isempty（i）
            s1（n:）=TempB（s:）;%如果i是空集则该节点为叶节点
            n=n+1;
        else s2（m:）=TempB（s:）;%如果i不是空集则该节点为非叶节点
            m=m+1;
        end
        s=s-1;
   end
    TempB=s2;             %将s2赋值给TempB.重复上述判断，直到TempB为空集为止
    s2=[];
end
s1
 %前推进行支路功率计算
     while k<10  %（一般取3~5次迭代即可收敛）
    Pij1=zeros（N1）; %该支路首端节点及与其相连的其他后续支路的功率情况
    Qij1=zeros（N1）;
      for s=1:mB        
        ii=s1（s3）;      %取排好序的支路首节点   
        jj=s1（s4）;      %取排好序的支路尾节点
        Pload=s1（s7）;
        Qload=s1（s8）-s1（s9）; %节点有功和无功负荷
        R