33节点潮流计算代码

采用前推回代的方法解决IEEE33节点潮流计算问题，代码内包含大量注释，方便学习，代码验证过，计算正确，都是学生共同进步。

clc
clear
Bus=[100;210060;39040;412080;56030;66020;7200100;8200100;96020;106020;114530;126035;136035;1412080;156010;166020;176020;189040;199040;209040;219040;229040;239040;24420200;25420200;266025;276025;286020;2912070;30200600;3115070;32210100;336040;];
Branch=[1120.09220.0407;2230.49300.2511;3340.36600.1864;4450.38110.1941;5560.81900.7070;6670.18720.6188;7780.71440.2351;8891.03000.7400;99101.04400.7400;1010110.19660.065;1111120.37440.1238;1212131.46801.1550;1313140.54160.7129;1414150.59100.5260;1515160.74630.5450;1616171.28901.7210;1717180.73200.5740;182190.16400.1565;1919201.50421.3554;2020210.40950.4784;2121220.70890.9373;223230.45120.3083;2323240.89800.7091;2424250.89600.7011;256260.20300.1034;2626270.28420.1447;2727281.05900.9337;2828290.80420.7006;2929300.50750.2585;3030310.97440.9630;3131320.31050.3619;3232330.34100.5302;];
[busnumrow]=size（Bus）;
[branchnumrow]=size（Branch）;
soubus=Branch（:2）;
mobus=Branch（:3）;
Vbus=12.66*ones（busnum1）;
Vbus1=Vbus;
Ploss=zeros（busnum1）;
Qloss=zeros（busnum1）;
e=1;
k=0;
Branch1=Branch;
n=1;
%% 精髓
%%%%%%%%%%%%%%%支路重新排序，各个分支线同时进行计算
while ~isempty（Branch1）%%%%T1为排好的支路矩阵。
    m=1;
    [srow]=size（Branch1）;
        while s>0
            t=find（Branch1（:2）== Branch1（s3））;%判断是否是叶子节点
                if isempty（t）%空，则是叶子节点
                    T1（n:）= Branch1（s:）;
                    %将叶子节点放入T1中，
                    %且从节点系统末端向首端进行；
                    %第二次由于Branch1排序问题从首端向末端，无差别，因为同事进行
                    n=n+1;
                else
                    T2（m:）= Branch1（s:）;%非叶子节点
                     m=m+1;
                end;
               s=s-1;
            end;
        Branch1=T2;
        T2=[];
end;
%% 
%%%%%%%%%%%%%%%%%%%%%%%
while e>1.0e-05%%%%收敛条件是精度
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    %前推功率
    %从末端向首端推
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    P=zeros（busnum1）;%%存放后续功率情况
    Q=zeros（busnum1）;
        for s=1:branchnum
                i=T1（s2）;
                j=T1（s3）;
                R=T1（s4）;
                X=T1（s5）;
                %按尾节点读取Bus数据
                Pload=Bus（j2）;
                Qload=Bus（j3）;
                II=（（Pload+P（j））^2+（Qload+Q（j））^2）/（Vbus（j）^2*1000）;
                Ploss（ij）=II*R;
                Qloss（ij）=II*X;
                P（ij）=Pload+Ploss（ij）+P（j）;
                Q（ij）=Qload+Qloss（ij）+Q（j）;
                P（i）=P（i）+P（ij）;
                Q（i）=Q（i）+Q（ij）;
         end;
     %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
    %后推电压
    %从首端向末端推
    %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
         for s=branchnum:-1:1
                i=T1（s2）;
                j=T1（s3）;%存储推导电压的顺序
                R=T1（s4）;
                X=T1（s5）;
                Vbus（j）=（Vbus（i）-（P（ij）*R+Q（ij）*X）/（Vbus（i）*1000））^2

 属性            大小     日期    时间   名称
----------- ---------  ---------- -----  ----

     文件      35328  2018-07-03 16:07  IEEE33配电网潮流计算\html\PF33.doc

     文件      35328  2018-06-24 15:20  IEEE33配电网潮流计算\html\powerflow33_.doc

     文件       3400  2018-07-03 16:06  IEEE33配电网潮流计算\PF33.m

     文件      12397  2018-06-20 20:16  IEEE33配电网潮流计算\完成分层的T1.xlsx

     文件      13189  2018-06-20 10:26  IEEE33配电网潮流计算\工作簿1.xlsx

     文件       6323  2018-06-20 10:20  IEEE33配电网潮流计算\新建 Microsoft Excel 工作表.xlsx

     文件      27689  2018-06-24 14:01  IEEE33配电网潮流计算\新建 Microsoft Visio 绘图.vsdx

     目录          0  2018-07-03 16:07  IEEE33配电网潮流计算\html

     目录          0  2018-07-08 22:38  IEEE33配电网潮流计算

----------- ---------  ---------- -----  ----

               133654                    9