編程實(shí)現(xiàn)戴維南諾頓等值牛頓拉夫遜潮流計(jì)算

1、高等電力網(wǎng)絡(luò)分析編程實(shí)現(xiàn)戴維南諾頓等值牛頓拉夫遜潮流計(jì)算目錄摘要2一、諾頓等值和戴維南等值21、諾頓等值和戴維南等值原理22、例題解析即程序結(jié)果分析2二、牛頓-拉夫遜法潮流分析21、牛頓-拉夫遜算法原理22、編程計(jì)算及結(jié)果分析2附錄 程序清單2參考文獻(xiàn)2摘要本文主要實(shí)現(xiàn)了兩種算法，一是通過對戴維南等值和諾頓等值電路的推導(dǎo)過程研究，總結(jié)所需求解電力系統(tǒng)參數(shù)，編程實(shí)現(xiàn)計(jì)算機(jī)計(jì)算網(wǎng)絡(luò)等值，使研究網(wǎng)絡(luò)規(guī)模大大減小，提高計(jì)算速度。二是利用牛頓-拉夫遜法進(jìn)行對電力系統(tǒng)進(jìn)行潮流計(jì)算，使潮流計(jì)算的速度和收斂性得以快速提高。關(guān)鍵字 戴維南等值和諾頓等值 牛頓-拉夫遜法潮流計(jì)算 Matlab一、諾頓等值和戴維南等

2、值1、諾頓等值和戴維南等值原理在電網(wǎng)分析中，有時(shí)需要研究從網(wǎng)絡(luò)的某一端口或多個(gè)端口看進(jìn)去時(shí)該網(wǎng)絡(luò)的表現(xiàn)。每個(gè)端口都是由感興趣的一對網(wǎng)絡(luò)節(jié)點(diǎn)組成的，其中一個(gè)節(jié)點(diǎn)還可以是公共參考節(jié)點(diǎn)（即地節(jié)點(diǎn)）。這時(shí)，可以把該電網(wǎng)在端口處看成一個(gè)等值的電流源或電壓源，但要求等值前后端口的電氣特性是相同的。這就是常規(guī)諾頓等值和戴維南等值的做法。應(yīng)用諾頓等值和戴維南等值對網(wǎng)絡(luò)進(jìn)行化簡，需要滿足兩個(gè)條件：（1）、被觀察的網(wǎng)絡(luò)是線性的；（2）、每個(gè)端口上的凈流入電流為零，即要求每個(gè)端口所連接的外部電路與被觀察網(wǎng)絡(luò)沒有電磁耦合，各個(gè)端口所連接的外部電路之間也沒有電氣耦合?；驹恚喝缦聢D所示，令原來的電力網(wǎng)絡(luò)有N個(gè)節(jié)點(diǎn)，地

3、節(jié)點(diǎn)作為參考節(jié)點(diǎn)不包括在內(nèi)。從中抽出m個(gè)感興趣的端口，這m個(gè)端口分別用下標(biāo)“，m”來表示，相應(yīng)端口上的節(jié)點(diǎn)對用（p，q），（k，l）等來表示。每個(gè)端口上第一個(gè)節(jié)點(diǎn)的電流以流出網(wǎng)絡(luò)為正方向，第二個(gè)節(jié)點(diǎn)的電流以流入網(wǎng)絡(luò)為正方向，二者大小相等。第一二個(gè)節(jié)點(diǎn)和第二個(gè)節(jié)點(diǎn)之間的電壓作為端口電壓的正方向。另外不是一般性，第二個(gè)節(jié)點(diǎn)還可能是參考節(jié)點(diǎn)（即地節(jié)點(diǎn)）。首先引入節(jié)點(diǎn)-端口關(guān)聯(lián)矢量和節(jié)點(diǎn)-端口關(guān)聯(lián)矩陣的概念。以端口為例，其對應(yīng)的N×1維節(jié)點(diǎn)-端口矢量為：P、q都不是參考點(diǎn)：M=0 1p -1q 0 Tq是參考點(diǎn)：M=0 1p 0 0 T把所有節(jié)點(diǎn)-端口關(guān)聯(lián)矢量按列排起來，就構(gòu)成了N×

4、;m維的節(jié)點(diǎn)-端口關(guān)聯(lián)矩陣ML=M M Mm T設(shè)系統(tǒng)原來的網(wǎng)絡(luò)方程是或者式中，為節(jié)點(diǎn)電壓列矢量；為節(jié)點(diǎn)注入電流列矢量；Y、Z分別為節(jié)點(diǎn)導(dǎo)納矩陣和節(jié)點(diǎn)阻抗矩陣。下圖（c）的多端口戴維南等值電路的m×m階等值阻抗矩陣為：戴維南等值電動勢即為原網(wǎng)絡(luò)的m個(gè)開路電壓，圖（b）的多端口諾頓等值電路的m×m階等值導(dǎo)納矩陣為諾頓等值電流源為（c）的網(wǎng)絡(luò)中各端口短路時(shí)的短路電流，根據(jù)前面規(guī)定的正方向，定義端口上的電流矢量和電壓矢量分別如下：從這些端口向原網(wǎng)絡(luò)看進(jìn)去，節(jié)點(diǎn)注入電流由兩部分組成，其一是圖（a）網(wǎng)絡(luò)內(nèi)部的節(jié)點(diǎn)注入電流，其二是與它連接的外部電路從端口注入的電流，因此，可以寫出網(wǎng)絡(luò)的

5、節(jié)點(diǎn)電壓方程如下：由此可得兩邊同乘，并考慮到，及。再考慮戴維南等值阻抗矩陣，則有多端口戴維南等值電路方程為：則諾頓等值電路方程為：2、例題解析即程序結(jié)果分析高等電力網(wǎng)絡(luò)分析P132 例5.3 如圖所示的電力系統(tǒng)，之路阻抗和節(jié)點(diǎn)注入電流都標(biāo)在圖上，試以節(jié)點(diǎn)和節(jié)點(diǎn)為一個(gè)端口，節(jié)點(diǎn)和地為一個(gè)端口，建立兩端口諾頓等值和戴維南等值。解：首先建立以地為參考節(jié)點(diǎn)的節(jié)點(diǎn)導(dǎo)納矩陣：其逆矩陣為兩個(gè)端口的關(guān)聯(lián)矢量組成了關(guān)聯(lián)矩陣。其中節(jié)點(diǎn)和 組成的端口 的關(guān)聯(lián)矢量，節(jié)點(diǎn)與地組成的端口的關(guān)聯(lián)矢量矩陣分別為： 則戴維南等值阻抗矩陣為：為求戴維南等值電動勢，首先求各節(jié)點(diǎn)電壓：求端口戴維南等值電動勢：求諾頓等值導(dǎo)納和諾頓等值

6、電流：程序流程見附錄對教材133頁例5.3驗(yàn)證，輸出結(jié)果為：請輸入節(jié)點(diǎn)數(shù):n=3請輸入支路數(shù):nl=3請輸入端口數(shù)m=2請輸入由支路參數(shù)形成的矩陣:B1=1,1,2i;2,2,2i;3,3,2i;1,2,0.2i;1,3,0.1i;2,3,0.2i;請輸入各節(jié)點(diǎn)參數(shù)形成的矩陣:B2=1,1;2,0.5,;3,-1.2;節(jié)點(diǎn)導(dǎo)納矩陣 Y= 0 -15.5000i 0 + 5.0000i 0 +10.0000i 0 + 5.0000i 0 -10.5000i 0 + 5.0000i 0 +10.0000i 0 + 5.0000i 0 -15.5000i端口節(jié)點(diǎn)發(fā)點(diǎn)m1=1端口節(jié)點(diǎn)收點(diǎn)m2=2端口節(jié)

7、點(diǎn)發(fā)點(diǎn)m1=3端口節(jié)點(diǎn)收點(diǎn)m2=0關(guān)聯(lián)矩陣M= 1 -1 0 0 0 1戴維南等值戴維南等值阻抗Zeq= 0 + 0.1164i 0 + 0.0127i 0 + 0.0127i 0 + 0.6970i戴維南等值電動勢Veq= 0 + 0.0044i 0 + 0.1440i諾頓等值諾頓等值導(dǎo)納Yeq= 0 - 8.6094i 0 + 0.1563i 0 + 0.1563i 0 - 1.4375i諾頓等值電流Ieq= 0.00230.1004分析與總結(jié)：（1）程序運(yùn)行結(jié)果與理論計(jì)算相一致，驗(yàn)證了程序在一定范圍內(nèi)的正確性，因程序的設(shè)計(jì)是在通用的基礎(chǔ)上，故在初始時(shí)刻根據(jù)所給電力系統(tǒng)的具體參數(shù)，對程序參

8、數(shù)進(jìn)行設(shè)置。（2）矩陣B1為支路參數(shù)矩陣，前N行為N個(gè)節(jié)點(diǎn)的對地阻抗，第一列為支路首端序號，第二列為支路末端序號，第三列為支路阻抗。矩陣B2為節(jié)點(diǎn)參數(shù)矩陣，第一列為節(jié)點(diǎn)序號，第二類為節(jié)點(diǎn)注入電流。二、牛頓-拉夫遜法潮流分析1、牛頓-拉夫遜算法原理電力系統(tǒng)潮流計(jì)算是電力系統(tǒng)運(yùn)行和規(guī)劃中最基本和最經(jīng)常的計(jì)算，其任務(wù)是在已知某些運(yùn)行參數(shù)的情況下，計(jì)算出系統(tǒng)中全部的運(yùn)行參數(shù)，一般來說，各個(gè)母線所供負(fù)荷的功率是已知的，各個(gè)節(jié)點(diǎn)電壓是未知的（平衡節(jié)點(diǎn)除外），可以根據(jù)網(wǎng)絡(luò)結(jié)構(gòu)形成節(jié)點(diǎn)導(dǎo)納矩陣，然后由節(jié)點(diǎn)導(dǎo)納矩陣和網(wǎng)絡(luò)拓?fù)浣Y(jié)構(gòu)列寫功率方程，由于功率方程里功率是已知的，電壓的幅值和相角是未知的，這樣潮流計(jì)算的問

9、題就轉(zhuǎn)化為求解非線性方程組的問題了。為了便于用迭代法解方程組，需要將上述功率方程改寫成功率平衡方程，并對功率平衡方程求偏導(dǎo)，得出對應(yīng)的雅可比矩陣，給未知節(jié)點(diǎn)賦電壓初值，一般為額定電壓，將初值帶入功率平衡方程，得到功率不平衡量，這樣由功率不平衡量、雅可比矩陣、節(jié)點(diǎn)電壓不平衡量（未知的）構(gòu)成了誤差方程，解誤差方程，得到節(jié)點(diǎn)電壓不平衡量，節(jié)點(diǎn)電壓加上節(jié)點(diǎn)電壓不平衡量構(gòu)成新的節(jié)點(diǎn)電壓初值，將新的初值帶入原來的功率平衡方程，并重新形成雅可比矩陣，然后計(jì)算新的電壓不平衡量，這樣不斷迭代，不斷修正，給定收斂精度最后即能收斂。牛頓拉夫遜算法的核心內(nèi)容是雅可比矩陣，采用極坐標(biāo)，節(jié)點(diǎn)電壓表示為節(jié)點(diǎn)功率將寫成式中，

10、是、兩節(jié)點(diǎn)電壓的相角差。實(shí)際上，對于每一個(gè)節(jié)點(diǎn)或每一個(gè)節(jié)點(diǎn)都可以列寫一個(gè)有功功率不平衡量方程式而對于每一個(gè)節(jié)點(diǎn)還可以再列寫一個(gè)無功功率不平衡量方程式對于上述方程式可以寫出修正方程式如下式中是階方陣，其元素為；是階矩陣，其元素為；是階矩陣，其元素為；是階方陣，其元素為?？梢缘玫窖趴杀染仃囋氐谋磉_(dá)式如下當(dāng)時(shí)當(dāng)時(shí)2、編程計(jì)算及結(jié)果分析對如圖所示系統(tǒng)進(jìn)行潮流計(jì)算分析，結(jié)果如下：導(dǎo)納矩陣 Y= 150 -50 0 -100 -50 75 -25 0 0 -25 50 -25 -100 0 -25 125雅可比矩陣 J= -150.4618 0 50.0541 0 0 0 0 149.8625 0 -5

11、0.0541 0 0 49.8625 0 -74.4929 0 24.9312 0 0 -49.8625 0 75.0946 0 -24.9312 0 0 25.0405 0 -50.2308 0 0 0 0 -25.0405 0 49.9312迭代次數(shù) K= 3每次迭代沒有達(dá)到精度要求節(jié)點(diǎn)個(gè)數(shù) N= 3 3 0最后一次的不平衡量 DW= 1.0e-007 * -0.0319 0 -0.1907 0 -0.0397 0最后一次的修正量 DV= 1.0e-009 * -0.1829 0 -0.4860 0 -0.3214 014各節(jié)點(diǎn)電壓 U: 1.0011 0.9972 1.0016 1.00

12、0014各節(jié)點(diǎn)的功率 S： 0.3000 -0.3000 0.1500 -0.1486分析與總結(jié)（1）利用牛頓-拉夫遜法對電力系統(tǒng)進(jìn)行潮流分析是比較復(fù)雜的一個(gè)過程，因此未能與前程序一樣對多有電力系統(tǒng)通用，而是在某個(gè)電力系統(tǒng)基礎(chǔ)上進(jìn)行編程，若有其他系統(tǒng)，可在此基礎(chǔ)上做局部修改。（2）矩陣B1為支路參數(shù)矩陣，第一列為支路首端號，第二列為支路末端號，第三列為支路串聯(lián)阻抗。矩陣B2為節(jié)點(diǎn)參數(shù)矩陣，第一列為注入發(fā)電功率，第二列為節(jié)點(diǎn)電壓參數(shù)，第三列為節(jié)點(diǎn)類型標(biāo)號，1為平衡節(jié)點(diǎn)，2為PQ節(jié)點(diǎn)。附錄 程序清單參考文獻(xiàn)1 張伯明,陳壽孫,嚴(yán)正.高等電力網(wǎng)絡(luò)分析2 何仰贊,溫增銀.電力系統(tǒng)分析.華中科技大學(xué)出版

13、社,2002.一、戴維南諾頓等值計(jì)算clearn=input('請輸入節(jié)點(diǎn)數(shù):n=');nl=input('請輸入支路數(shù):nl=');m=input('請輸入端口數(shù)m=');B1=input('請輸入由支路參數(shù)形成的矩陣:B1=');B2=input('請輸入各節(jié)點(diǎn)參數(shù)形成的矩陣:B2='); Y=zeros(n);Z=zeros(n);M=zeros(m,n);f=zeros(1,n);Zeq=zeros(m);V=zeros(n,1);Veq=zeros(m,1);Yeq=zeros(m);I=zeros(n,

14、1);Ieq=zeros(m,1);%-求節(jié)點(diǎn)導(dǎo)納矩陣for i=(n+1):(n+nl) p=B1(i,1);q=B1(i,2);Y(p,q)=Y(p,q)-1./B1(i,3); Y(q,p)=Y(p,q); Y(q,q)=Y(q,q)+1./B1(i,3); Y(p,p)=Y(p,p)+1./B1(i,3); endfor p=1:n Y(p,p)=Y(p,p)+1./B1(p,3);end disp('節(jié)點(diǎn)導(dǎo)納矩陣 Y=');disp(Y)%-求節(jié)點(diǎn)端口關(guān)聯(lián)矩陣for j=1:m m1=input('端口節(jié)點(diǎn)發(fā)點(diǎn)m1='); m2=input('

15、端口節(jié)點(diǎn)收點(diǎn)m2='); for k=1:n if k=m1 f(1,k)=1; else if k=m2&&m2=0 f(1,k)=-1; else f(k)=0; end end end M(j,:)=f;end disp('關(guān)聯(lián)矩陣M=');disp(M)%-求端口戴維南等值阻抗disp('戴維南等值')Z=inv(Y);Zeq=M*Z*M'disp('戴維南等值阻抗Zeq=');disp(Zeq)%-求各節(jié)點(diǎn)電壓for l=1:n I(l,1)=B2(l,2);endV=Z*I;%-求端口戴維南等值電動勢V

16、eq=M*V;disp('戴維南等值電動勢Veq=');disp(Veq)%-求端口諾頓等值disp('諾頓等值')Yeq=inv(Zeq);Ieq=Zeq'*Veq;disp('諾頓等值導(dǎo)納Yeq=');disp(Yeq);disp('諾頓等值電流Ieq=');disp(Ieq)二、牛頓-拉夫遜潮流計(jì)算clear;n=4;%節(jié)點(diǎn)數(shù)：n;nl=4;%支路數(shù)：nl;isb=4;%平衡母線節(jié)點(diǎn)號：isb;pr=0.00001;%誤差精度：pr=;B1=1 2 0.02;1 4 0.01;2 3 0.04;3 4 0.04;B2

17、=0.3 1 2;-0.3 1 2;0.15 1 2;0 1 1;disp('一、潮流分析結(jié)果：');Y=zeros(n);e=zeros(1,n);f=zeros(1,n);I=zeros(1,n);%-求導(dǎo)納矩陣for i=1:nlp=B1(i,1);q=B1(i,2);Y(p,q)=Y(p,q)-1./B1(i,3);Y(q,p)=Y(p,q); Y(q,q)=Y(q,q)+1./B1(i,3);Y(p,p)=Y(p,p)+1./B1(i,3);enddisp('導(dǎo)納矩陣 Y=');disp(Y)%-K=0;N=1; %K迭代次數(shù)；N不滿足收斂要求的節(jié)點(diǎn)數(shù)

18、;for i=1:ne(i)=real(B2(i,2);f(i)=imag(B2(i,2);S(i)=B2(i,1); endG=real(Y);B=imag(Y);P=real(S);Q=imag(S);%-求雅可比矩陣 while N=0 N=0;K=K+1;for i=1:nif i=isbC(i)=0;D(i)=0;for j1=1:nC(i)=C(i)+G(i,j1)*e(j1)-B(i,j1)*f(j1);%(Gij*ej-Bij*fj)D(i)=D(i)+G(i,j1)*f(j1)+B(i,j1)*e(j1);%(Gij*fj+Bij*ej)endP0=C(i)*e(i)+f(i

19、)*D(i);%節(jié)點(diǎn)功率P=ei(Gij*ej-Bij*fj)+fi(Gij*fj+Bij*ej)Q0=C(i)*f(i)-e(i)*D(i);%節(jié)點(diǎn)功率Q=fi(Gij*ej-Bij*fj)-ei(Gij*fj+Bij*ej)DP=P(i)-P0;DQ=Q(i)-Q0; %-求不平衡量-v=2*i-1;w=2*i;DW(v)=DP;DW(w)=DQ; %-for j1=1:nif j1=isb&&j1=iX1=-G(i,j1)*e(i)-B(i,j1)*f(i);% X1=dP/de X2=B(i,j1)*e(i)-G(i,j1)*f(i);% X2=dP/df dP/de=

20、-dQ/dfX3=X2; % X3=dQ/de dP/df=dQ/deX4=-X1; % X4=dQ/dfp=2*i-1;q=2*j1-1;J(p,q)=X1;m=p+1; J(m,q)=X3;q=q+1; J(p,q)=X2;J(m,q)=X4; elseif j1=i&&j1=isbX1=-C(i)-G(i,i)*e(i)-B(i,i)*f(i);X2=-D(i)+B(i,i)*e(i)-G(i,i)*f(i);X3=D(i)+B(i,i)*e(i)-G(i,i)*f(i); X4=-C(i)+G(i,i)*e(i)+B(i,i)*f(i);p=2*i-1;q=2*j1-1;J(p,q)=X1;m=p+1;J(m,q)=X3;q=q+1;J(p,q)=X2;J(m,q)=X4;endendendendDV=-J(-1)*DW'%disp('每次迭代雅可比矩陣 J=