【正文】 end end % 初始化并計(jì)算 Qi和 Pi 29 for i = 1:nb1 for j = 1:nb1 if i~=j H(i,j)=bus(i,2)*bus(j,2)*(real(Y(i,j))*sin(bus(i,3)bus(j,3))imag(Y(i,j))*cos(bus(i,3)bus(j,3)))。 for i = 1:nb1 for j = 1:nb Pi(i,1)=Pi(i,1)+bus(i,2)*bus(j,2)*(real(Y(i,j))*cos(bus(i,3)bus(j,3))+imag(Y(i,j))*sin(bus(i,3)bus(j,3)))。 K L]，并初始化 Qi = zeros(nb1,1)。 L = zeros(nPQ,nPQ)。 N = zeros(nb1,nPQ)。 end end end % 利用循環(huán)計(jì)算求取 dP 和 dQ 子程序 5 作用為計(jì)算雅克比矩陣 function J = Jac_(bus,Y,nPQ) [nb,mb]=size(bus)。 % 對 dP 和 dQ賦初值 PV節(jié)點(diǎn)不需計(jì)算 dQ 平衡節(jié)點(diǎn)不參與計(jì)算 for i = 1:n1 for j = 1:n dP(i,1) = dP(i,1)bus(i,2)*bus(j,2)*(real(Y(i,j))*cos(bus(i,3)bus(j,3))+imag(Y(i,j))*sin(bus(i,3)bus(j,3)))。 % 總節(jié)點(diǎn)個(gè)數(shù) dP = bus(1:n1,4)。 Y(J,I)=Y(I,J)。 Y(J,J)=Y(J,J)+K*K*Yt。 Y(J,I)=Y(I,J)。 Y(J,J)=Y(J,J)+Yt/K/K。(J==0) % 對地支路 : K=0,J=0,R=X=0 Y(I,I)=Y(I,I)+Ym。 Y(J,I)=Y(I,J)。 Y(J,J)=Y(J,J)+Yt+Ym。 if (K==0)amp。 % 接地支路不計(jì)算 Yt end Ym=line(k,5)+j*line(k,6)。 Zt=line(k,3)+j*line(k,4)。 % 對導(dǎo)納矩陣賦初值 0 for k=1:nl I=line(k,1)。 [nl,ml]=size(line)。% 生成新舊節(jié)點(diǎn)對照表 for i=1:nl for j=1:2 for k=1:nb if line(i,j)==nodenum(k,2) line(i,j)=nodenum(k,1)。 % 對 bus 矩陣按 PQ、 PV、平衡節(jié)點(diǎn)的順序重新排序 nodenum=[[1:nb]39。PV。 PQ(nPQ,:)=bus(i,:)。 PV(nPV,:)=bus(i,:)。 SW(nSW,:)=bus(i,:)。 % nPQ 為 PQ節(jié)點(diǎn)個(gè)數(shù) for i = 1:nb, % nb 為總節(jié)點(diǎn)數(shù) type= bus(i,6)。 % nSW 為平衡節(jié)點(diǎn)個(gè)數(shù) 26 nPV = 0。 [nl,ml]=size(line)。 % 讀取文件字符串長度 eval_r(dfile(1:1file2))。 you must select a valid data file39。)。,39。 % 迭代結(jié)束后繼續(xù)在 “”寫入節(jié)點(diǎn)計(jì)算結(jié)果和線路計(jì)算結(jié)果 程序結(jié)束 子程序 1 作 用為打開數(shù)據(jù)文件 function [bus,line] = OpDF_ [dfile,pathname]=uigetfile(39。,real(S_res(i,5)),imag(S_res(i,5)))。 fprintf(myf,39。.6f + j .6f 39。,real(S_res(i,3)),imag(S_res(i,3)))。 fprintf(myf,39。% 39。,S_res(i,1))。 for i = 1:nl fprintf(myf,39。\n 線路計(jì)算結(jié)果： \n 節(jié)點(diǎn) I 節(jié)點(diǎn) J 線路功率 S(I,J) 線路功率S(J,I) 線路損耗 dS(I,J)\n39。,real(bus_res(i,4)),imag(bus_res(i,4)))。 fprintf(myf,39。.6f 39。,bus_res(i,2))。 fprintf(myf,39。% 39。)。 fprintf(myf,39。a39。39。 % 計(jì)算節(jié)點(diǎn)數(shù)據(jù)結(jié)果的子程序 S_res = S_res_(bus,line,YtYm)。 % 對節(jié)點(diǎn)恢復(fù)編號的子程序 YtYm = YtYm_(line)。(max(abs(dAng))EPS) break end % 判斷是否滿足 精度誤差，如滿足則跳出，否則返回繼續(xù)迭代 end bus = PQ_(bus,Y,nPQ,nPV)。 % 修正電壓 bus(1:nb1,3) = bus(1:nb1,3) dAng。 % 計(jì)算相角修正量 dU = UD*(dAngU(nb:nb+nPQ1,1))。dQ]。 for i = 1:nPQ UD(i,i) = bus(i,2)。 % 計(jì)算功率偏差 dP 和 dQ 的子程序 J = Jac_(bus,Y,nPQ)。 % 在當(dāng)前目錄下生成 “”文件，寫入節(jié)點(diǎn)導(dǎo)納矩 24 陣 format long EPS = 。)。 fprintf(myf,39。w39。39。 % 對節(jié)點(diǎn)重新排序的子程序 Y = y_(bus,line)。 [nl,ml]=size(line)。 程序 設(shè)計(jì) [bus,line] = OpDF_。在設(shè)計(jì)研究單位和工業(yè)開發(fā)部門， MATLAB 被廣泛的應(yīng)用于研究和解決各種具體問題。在歐美等國家的高校， MATLAB 已成為線性代數(shù)、自動控制理論、數(shù)理統(tǒng)計(jì)、數(shù)字信號處理、時(shí)間序列分析、動態(tài)系統(tǒng)仿真等高級課程的基本教學(xué)工 具。 MATLAB 具有 編程效率高 、用戶使用方便 、 擴(kuò)充能力強(qiáng) 、 語句簡單，內(nèi)涵豐富 、 高效方便的矩陣和數(shù)組運(yùn)算 、 方便的繪圖功能 等特點(diǎn)，給用戶帶來了極大的方便。 22 輸 入 原 始 數(shù) 據(jù)形 成 節(jié) 點(diǎn) 導(dǎo) 納 矩 陣設(shè) 電 壓 初 值 e( 0 )、 f( 0 )設(shè) 迭 代 次 數(shù) k = 0計(jì) 算 誤 差 向 量? P( k )、 ? Q( k )、 ? U2( k )收 斂 否求 雅 可 比 矩 陣 元 素解 修 正 方 程 ，求 解 ? e( k )、 ? f( k )修 正 節(jié) 點(diǎn) 電 壓e( k + 1 )= e( k ) ? e( k )f( k + 1 )= f( k ) ? f( k )K = K + 1K K m a x不 收 斂 停 機(jī)求 P U 節(jié) 點(diǎn) 無 功 功 率 ，求 平 衡 節(jié) 點(diǎn) 功 率求 支 路 功 率 分 布 和 損 耗停 機(jī)是否 圖 3— 2 牛頓拉夫遜法計(jì)算潮流程序框圖 23 第四章程序 設(shè)計(jì) 及 算例分析 Matlab 概述 MATLAB (Matrix Laboratory)為美國 Mathworks 公司 1983年首次推出的一套高性能的數(shù)值分析和計(jì)算軟件，其功能不斷擴(kuò)充，版本不斷升級 。 圖中 3— 2 中示出牛頓拉夫遜法計(jì)算潮流的程序框圖。當(dāng)采用直角坐標(biāo)表式時(shí)，要增加一個(gè)對應(yīng)于該節(jié)點(diǎn)的無功功率不平衡量 ? ?Q? 的關(guān)系式，同時(shí)要減少一個(gè)對應(yīng)于設(shè)定節(jié)點(diǎn)電壓約束條件 ? ?2U? 的形式。 需要指出的是，當(dāng)計(jì)算過程中發(fā)生 PU節(jié)點(diǎn)的無功功率越限時(shí)， PU節(jié)點(diǎn)要轉(zhuǎn)化成 PQ節(jié)點(diǎn)。 （ 4）具有結(jié)構(gòu)對稱性，但數(shù)值不對稱。 21 （ 3）與導(dǎo)納矩陣具有相似的結(jié) 構(gòu)，當(dāng) 0?ijY 時(shí)， ijH 、 ijI 、 ijJ 、 ijL 均為 0，因此，也是高度稀疏的矩陣。當(dāng)節(jié)點(diǎn)電壓以極坐標(biāo)表式時(shí)，該矩陣為 mn ?? )1(2階方陣；節(jié)點(diǎn)電壓以直角坐標(biāo)表式時(shí)，為 2( 1)n? 階。 牛頓 — 拉夫遜法潮流計(jì)算的步驟 求取雅可比矩陣是牛頓 — 拉夫遜法的一項(xiàng)重要的工作。其雅可比矩陣元素為 ????????????????????????????????????????????????iiiiiiijiiiijiiiiiiiijiiiiiiiijiiiiiiiijiiiiiiiijiijijijiijiijijijeSfRjiSReBfGbLfBeGaJfBeGaNeBfGbHjieBfGJNjieBfGLH22)(0)()( （ 332） 上式中 ia 和 ib 分別是節(jié)點(diǎn)注入電流 iI? 的實(shí)部和虛部，即 1 (ni ij j ij jja G e B f?? ? ??????? （ 333） 1 ()ni ij j ij jjb G f B e???? （ 334） 由式 4— 34 可求得第 1?k 次迭代的修正量 )1( ??ke 和 )1( ??kf ，從而可得到新的解 ( 1 ) ( ) ( 1 )( 1 ) ( ) ( 1 )k k kk k ke e ef f f??? ? ? ? ? ????? ? ? ? ? ??? ? ? ? ? ? （ 335） 這樣反復(fù)計(jì)算，直到收斂至要求的精度。 對于除平衡節(jié)點(diǎn)以外的所有節(jié)點(diǎn) 1 1 1 1n n n niii j i ij i ij jj j j jjjPPP f e H f N efe? ? ? ?? ? ? ?? ? ? ? ? ? ? ? ???? ? ? ? （ 328） 對于 PQ 節(jié)點(diǎn) 1 1 1 1n n n niii j i ij i ij jj j j jjjQ f e J f L efe? ? ? ?? ? ? ?? ? ? ? ? ? ? ? ???? ? ? ? （ 329） 對于 PU 節(jié)點(diǎn) 222 iii i i ii i ii iUUU f e R f S efe? ? ?? ? ? ? ? ? ? ? ??? （ 330） 如系統(tǒng)具有 n個(gè) 節(jié)點(diǎn)，其中 m 個(gè)為 PU 節(jié)點(diǎn)，則式 428~429 式可寫成如下矩陣形式 2P H N fQ J L eRSU??? ????? ????? ??? ? ??? ??????? ?????? ???? （ 331） 20 式中 P? 、 e? 和 f? 為 1)1( ??n 向量， Q? 為 1)1( ???mn 向量， 2U? 為 1?m 向量。 節(jié)點(diǎn)電壓以直角坐標(biāo)表示，即 iii jfeU ??? ，