%高斯——賽爾德潮流計算
%[s_nodeu_node]=gaosai（YmnsS_PQP_PVU_PVusQmaxQminexpl）
%返回值：s_node:??節(jié)點功率
%???????u_node：??節(jié)點電壓
%輸入值：Y：???導(dǎo)納矩陣
%????????m：???PQ節(jié)點數(shù)包含平衡節(jié)點
%????????n：???總結(jié)點數(shù)
%????????s：???平衡節(jié)點編號，一般定為1號編號
%????????S_PQ:?PQ節(jié)點的視在功率[1:m]其中平衡節(jié)點的視在功率暫定為0；
%????????P_PV：PV節(jié)點的有功功率[m+1:n]
%????????U_PV:?PV節(jié)點的電壓給定值只有幅值，沒有相位角，是實數(shù)[m+1:n]
%????????us:???平衡節(jié)點的電壓給定值，復(fù)數(shù)。
%????????Qmax；PV節(jié)點無功功率最大值
%????????Qmin：PV節(jié)點無功功率最小值
%????????epxl:?最大電壓變量dumax的最大值，即規(guī)定的潮流計算的精度，實數(shù)。
function?[u_node]=DYgaosai（YmnsS_PQP_PVU_PVusQmaxQminepxl）
%**************************************************************************
%以下為高斯塞爾德潮流計算程序
?fprintf（‘高斯——塞爾德潮流計算開始\n‘）
U=ones（1n）;%節(jié)點電壓的初始值；
U（s）=us;%將電壓初始值的平衡節(jié)點的初始電壓設(shè)為給定值
k=1000;%設(shè)置迭代次數(shù)
dumax=0;%最大電壓變化量，實數(shù)
flag=zeros（1n-m）;%PV節(jié)點在計算過程中轉(zhuǎn)化為PQ節(jié)點的標(biāo)志，PV轉(zhuǎn)化為PQ節(jié)點時，相應(yīng)的標(biāo)志位置1
Q_PV=ones（1n-m）;%PV節(jié)點的無功功率；
Q_PV=Q_PV.*0.2;
Ss=0;%平衡節(jié)點的視在功率
for?a=m+1:n
????U（a）=U_PV（a-m）;%第a個節(jié)點是PV節(jié)點，暫時將PV節(jié)點電壓給定值賦給節(jié)點電壓初始值其相角暫時為0
end
while（k>=1）
for?i=1:n
????u=U（i）;%保存U（i）的初值，為求du和dumax做準(zhǔn)備；
????if（i==s）%第i個節(jié)點是否為平衡節(jié)點？
???????%?i=i+1;%d第i個節(jié)點是平衡節(jié)點，跳過。
????else????%第i節(jié)點不是平衡節(jié)點，進行以下處理
????????A=0;
????????for?a=1:n
????????????if（a~=i）
?????????????A=A+Y（ia）*U（a）;
????????????end
????????end
????????a=1;
????????if（i>=m+1?&&?flag（i-m）==0）?%第i各節(jié)點是否為PV節(jié)點？
????????????Q_PV（i-m）=-imag（U（i）‘*Y（ii）*U（i）+U（i）‘*A）;
????????????if（Q_PV（i-m）<=Qmax?&&?Q_PV（i-m）>=Qmin）
????????????????Upv=（（P_PV（i-m）-Q_PV（i-m）*j）/U（i）-A）/Y（ii）;
????????????????ang=angle（Upv）;
????????????????U（i）=abs（U（i））*（cos（ang）+sin（ang）*j）;
%???????????????Q_PV（i-m）=-imag（U（i）‘*Y（ii）*U（i）+U（i）‘*A）;
????????????else
????????????????if（Q_PV（i-m）>Qmax）
????????????????????Q_PV（i-m）=Qmax;
????????????????????flag（i-m）=1;
????????????????????fprintf（‘PV節(jié)點%gQ大于Qmax轉(zhuǎn)化為PQ節(jié)點\n‘i）;
????????????????elseif（Q_PV（i-m）????????????????????Q_PV（i-m）=Qmin;
????????????????????flag（i-m）=1;
????????????????????fprintf（‘PV節(jié)點%gQ小于Qmin轉(zhuǎn)化為PQ節(jié)點\n‘i）;
????????????????end
????????????????U（i）=（（P_PV（i-m）-Q_PV（i-m）*j）/U（i）‘-A）/Y（ii）;
????????????end
????????else
????????????if（i<=m）%第i個節(jié)點是PQ節(jié)點
????????????????U（i）=（S_PQ（i）‘/U（i）‘-A）/Y（ii）;
????????????elseif（i>=m+1?&&?flag（i-m）==1）
????????????????U（i）=（（P_PV（i-m）-Q_PV（i-m）*j）/U（i）‘-A）/Y（ii）;
????????????end
????????end
????end
?du=abs（U（i）-u）;
??if（du>dumax）
??dumax=du;
??end
end
i=1;
if（dumax????fprintf（‘迭代次數(shù)%g\n‘1000-k+1）
????break;
else
????k=k-1;
????dumax=0;
end
end
if（k==0）
????fprintf（‘高斯——塞爾德潮流計算不能收斂！！！\n‘）
end
for?i=1:n?%計算平衡節(jié)點的視在功率
????????Ss=Ss+Y（si）‘*U（i）‘;
end
Ss=Ss*U（s）;%平衡節(jié)點的視在功率
????for?i=1:n%計算各個節(jié)點的視在功率
????????if（i==s）
????????????s_node（i）=Ss;
????????else
????????????if（i<=m）
????????????????s_node（i）=S_PQ（i）;
????????????elseif（i>=m+1）
????????????????s_node（i）=P_PV（i-m）+Q_PV（i-m）*j;
????????????end
????????end
????end
????s_node
????u_node=U
????fprintf（‘高斯——塞爾德潮流計算結(jié)束\n‘）

????
????????????
????????????

?屬性????????????大小?????日期????時間???名稱
-----------?---------??----------?-----??----

?????文件???????3531??2011-05-07?20:29??newton\DYgaosai.m

?????文件??????21808??2011-06-06?09:25??newton\newton.m

?????文件???????2597??2011-05-12?08:59??newton\write.m

?????文件??????16384??2011-06-06?09:25??newton\輸出數(shù)據(jù).xls

?????目錄??????????0??2011-10-18?08:58??newton

-----------?---------??----------?-----??----

????????????????44320????????????????????5