clc;clear;
%%?使用IEEE33節點配電網前推回代潮流計算程序?
N=33;
Branchdata?=?[
????1???2???0.0922??0.0470;
????2???3???0.4930??0.2511;
????3???4???0.3660??0.1864;
????4???5???0.3811??0.1941;
????5???6???0.8190??0.7070;
????6???7???0.1872??0.6188;
????7???8???0.7114??0.2351;
????8???9???1.0300??0.7400;
????9???10??1.0440??0.7400;
????10??11??0.1966??0.0650;
????11??12??0.3744??0.1238;
????12??13??1.4680??1.1550;
????13??14??0.5416??0.7129;
????14??15??0.5910??0.5260;
????15??16??0.7463??0.5450;
????16??17??1.2890??1.7210;
????17??18??0.7320??0.5740;
????2???19??0.1640??0.1565;
????19??20??1.5042??1.3554;
????20??21??0.4095??0.4784;
????21??22??0.7089??0.9373;
????3???23??0.4512??0.3083;
????23??24??0.8980??0.7091;
????24??25??0.8960??0.7011;
????6???26??0.2030??0.1034;
????26??27??0.2842??0.1447;
????27??28??1.0590??0.9337;
????28??29??0.8042??0.7006;
????29??30??0.5075??0.2585;
????30??31??0.9744??0.9630;
????31??32??0.3105??0.3619;
????32??33??0.3410??0.5302;
];?????????????????????????????????????%?支路，阻抗
?
Nodedata?=?[?1?????0???????0;
????2???100.00??60.00;
????3???90.00???40.00;
????4???120.00??80.00;
????5???60.00???30.00;
????6???60.00???20.00;
????7???200.00??100.00;
????8???200.00??100.00;
????9???60.00???20.00;
????10??60.00???20.00;
????11??45.00???30.00;
????12??60.00???35.00;
????13??60.00???35.00;
????14??120.00??80.00;
????15??60.00???10.00;
????16??60.00???20.00;
????17??60.00???20.00;
????18??90.00???40.00;
????19??90.00???40.00;
????20??90.00???40.00;
????21??90.00???40.00;
????22??90.00???40.00;
????23??90.00???50.00;
????24??420.00??200.00;
????25??420.00??200.00;
????26??60.00???25.00;
????27??60.00???25.00;
????28??60.00???20.00;
????29??120.00??70.00;
????30??200.00??600.00;
????31??150.00??70.00;
????32??210.00??100.00;
????33??60.00???40.00;
];??????????????????????????????????????%?節點，負荷
?
UB?=?12.66;?????????????????????????????%?電壓基準?kV
SB?=?10;????????????????????????????????%?功率基準?MVA
ZB?=?UB^2/SB;???????????????????????????%?阻抗基準?ohm
?
Branch?（:[45]）?=?Branchdata（:[34]）?/?ZB;?????%?阻抗標幺化
Branch?（:[23]）?=?Branchdata（:[12]）;
Branch?（:1）?=?[1:1:N-1]‘;

Bus（:[23]）?=?Nodedata（:[23]）?/?SB?/?1000;%?功率標幺化
Bus（:1）?=[1:1:N]‘;

%%?本段程序判斷所選線路是否收斂:節點電壓，支路電流功率損耗
[busnumdump]?=?size（Bus）;
Vbus?=?ones（busnum1）;?????%取額定電壓10KV????
Vbus（11）?=?1;????????????%x設定線路首端電壓為額定電壓??
Vbus2?=?Vbus;??????????????%迭代后的電壓矩陣
[branchnumdump]?=?size（Branch）;

k=0;?????????????%迭代次數
Ploss=zeros（branchnum1）;%存支路的有功損耗??
Qloss=zeros（branchnum1）;%支路無功損耗

I=zeros（branchnum1）;
F=zeros（busnum1）;
%%?本段程序將支路重新排序???s1為排好序的支路矩陣
TempBranch=Branch;
n=1;
s2=[];?????????%這段重新排列支路應該沒有問題
while?~isempty（TempBranch）
????[sdump]=size（TempBranch）;%s為支路數??
????m=1;
????while?s>0
???????i=find（TempBranch（:2）==TempBranch（s3））;%末端節點是否為其他支路首端節點
????????if?isempty（i）
????????????s1（n:）=TempBranch（s:）;%如果i是空集則該節點為葉節點
????????????n=n+1;
????????else?s2（m:）=TempBran