%形成鄰接矩陣
data=textread（‘IEEE33.dat‘）;?
a=zeros（3333）;
for?j=1:32;
??a（data（j2）+1data（j3）+1）=1;
end
%%電壓賦初值
v=ones（331）;
v（1）=1.05;
delt=zeros（331）;
V=v.*cos（delt）+i*v.*sin（delt）;
%%節點負荷功率
SL=zeros（331）;
for?k=2:33
???SL（k）=（data（k-16）+i*data（k-17））/10;
end
SL（1）=0;
%%線路阻抗
Z=zeros（331）;
for?m=2:33
????Z（m）=（data（m-14）+i*data（m-15））/100;
end
Z（1）=0;
S=zeros（331）;
%迭代
for?dd=1:100
?????V=v.*cos（delt）+i*v.*sin（delt）;
?????U=V;
?????u=abs（U）
????for?j=1:33
????????if?a（j:）==0
???????????S（j）=SL（j）+conj（SL（j））*SL（j）/（conj（V（j））*V（j））*Z（j）;
????????end
????end
????
?