《MATLAB實驗 電力系統(tǒng)暫態(tài)穩(wěn)定分析》由會員分享,可在線閱讀,更多相關《MATLAB實驗 電力系統(tǒng)暫態(tài)穩(wěn)定分析(9頁珍藏版)》請在裝配圖網(wǎng)上搜索。
1、精品文檔,僅供學習與交流,如有侵權請聯(lián)系網(wǎng)站刪除
實驗三 電力系統(tǒng)暫態(tài)穩(wěn)定分析
電力系統(tǒng)暫態(tài)穩(wěn)定計算實際上就是求解發(fā)電機轉(zhuǎn)子運動方程的初值問題,從而得出δ-t和ω-t的關系曲線。每臺發(fā)電機的轉(zhuǎn)子運動方程是兩個一階非線性的常微分方程。因此,首先介紹常微分方程的初值問題的數(shù)值解法。
一、 常微分方程的初值問題
(一)問題及求解公式的構造方法
我們討論形如式(3-1)的一階微分方程的初值問題
(3-1)
設初值問題(3-1)的解為,為了求其數(shù)值解而采取離散化方法,在求解區(qū)間[]上取一組節(jié)點
稱()為步長。在等步長的情況下,步長為
用表示在節(jié)點處解的準確值的近似值。
2、
設法構造序列所滿足的一個方程(稱為差分方程)
(3-2)
作為求解公式,這是一個遞推公式,從(,)出發(fā),采用步進方式,自左相右逐步算出在所有節(jié)點上的近似值()。
在公式(3-2)中,為求只用到前面一步的值,這種方法稱為單步法。在公式(3-2)中的由明顯表示出,稱為顯式公式。而形如(3-3)
(3-3)
的公式稱為隱式公式,因為其右端中還包括。
如果由公式求時,不止用到前一個節(jié)點的值,則稱為多步法。
由式(3-1)可得
= (3-4)
兩邊在[,]上積分,得
(3-5)
由此可以看出,如果想構造求
3、解公式,就要對右端的積分項作某種數(shù)值處理。這種求解公式的構造方法叫做數(shù)值積分法。
(二)一般的初值問題的解法
1. 歐拉法和改進歐拉法
對于初值問題(3-1),采用數(shù)值積分法,從而得到(3-5)。對于(3-5)右端的積分用矩形公式(取左端點),則得到
進而得到(3-1)的求解公式(3-2)
(=0,1,2,n-1) (3-6)
此公式稱為歐拉(Euler)格式。
如果對式(3-5)右端的積分用梯形公式
則可以得到初值問題(3-1)的梯形求解公式如式(3-7)
(=0,1,2,n-1) (3-7)
式(3-7)是個隱式公式??梢圆扇∠扔脷W拉格式求一個的初步
4、近似值,記作,稱之為預報值,然后用預報值替代式(3-7)右端的,再計算得到,稱之為校正值,這樣建立起來的預報-校正方法稱為改進歐拉格式
(3-8)
2. 龍格—庫塔方法
在單步法中,應用最廣泛的是龍格-庫塔(Runge-kutta)法,簡稱R-K法。下面直接給出一種四階的龍格-庫塔法的計算公式(3-9)
(3-9)
它也稱為標準(古典)龍格-庫塔法。
例3-1 研究下列微分方程的初值問題
解:
這是一個特殊的微分方程,其解的解析式可以給出,為
應用龍格-庫塔法,取=0.25,根據(jù)式(3-9)編寫一段程序,由零開始自左相右
5、逐步算出在所有節(jié)點上的近似值。計算結(jié)果見表3-1。計算結(jié)果表明,四階龍格-庫塔方法的精度是較高的。
表3-1
2.0
0.39995699
4.3e-5
4.0
0.23529159
2.5e-6
6.0
0.16216179
3.7e-7
8.0
0.12307683
9.2e-8
實際上,MATLAB為常微分方程提供了很好的解題指令,使得求解常微分方程變得很容易,并且能將問題及解答表現(xiàn)在圖形上。因此,我們可以不用根據(jù)式(3-9)編寫較復雜的程序,而只需應用MATLAB提供的常微分方程解題器來解決問題。下面給出用MATLAB編寫的解題程序。
首先編寫
6、描述常微分方程的ODE文件,文件名為ˊmyfunˊ,便于解題器調(diào)用它。
function dy = myfun(x,y)
dy = zeros(1,1);
dy=1/(1+x^2)-2*y^2;
再編寫利用解題器指令求解y的程序。
clear
x0=0;
for i=1:4
xm=2*i;
y0=0;
[x,y] = ode45('myfun',[x0 xm],[y0]);
format long
y(length(y))
end
plot(x,y,'-')
運行上述程序,在得到幾個點的函數(shù)值的同時,也得到函數(shù)y的曲線,如圖3-1所示。
圖3-1 根據(jù)運算結(jié)果
7、畫出y的曲線
二、 簡單電力系統(tǒng)的暫態(tài)穩(wěn)定性
(一)物理過程分析
某簡單電力系統(tǒng)如圖3-2(a)所示,正常運行時發(fā)電機經(jīng)過變壓器和雙回線路向無限大系統(tǒng)供電。發(fā)電機用電勢作為其等值電勢,則電勢與無限大系統(tǒng)間的電抗為
(3-10)
這時發(fā)電機發(fā)出的電磁功率可表示為
(3-11)
如果突然在一回輸電線路始端發(fā)生不對稱短路,如圖3-2(b)所示。故障期間發(fā)電機電勢與無限大系統(tǒng)之間的聯(lián)系電抗為
(3-12)
在故障情況下發(fā)電機輸出的電磁功率為
(3-13)
在短路故障發(fā)生之后,線路
8、繼電保護裝置將迅速斷開故障線路兩端的斷路器,如圖3-2(c)所示。此時發(fā)電機電勢與無限大系統(tǒng)間的聯(lián)系電抗為
(3-14)
發(fā)電機輸出的功率為
(3-15)
圖3-2 簡單電力系統(tǒng)及其等值電路
(a)正常運行方式及其等值電路;(b)故障情況及其等值電路;(c)故障切除后及其等值電路
如果正常時發(fā)電機向無限大系統(tǒng)輸送的有功功率為,則原動機輸出的機械功率等于。假定不計故障后幾秒種之內(nèi)調(diào)速器的作用,即認為機械功率始終保持。因此,可以得到此簡單電力系統(tǒng)正常運行、故障期間及故障切除后的功率特性曲線如圖3-3所示。
圖3-3 簡單系統(tǒng)正常運
9、行、故障期間及故障切除后的功率特性曲線
對于上述簡單電力系統(tǒng),我們可以根據(jù)等面積定則求得極限切除角。但是,實際工作需要知道在多少時間之內(nèi)切除故障線路,也就是要知道與極限切除角對應的極限切除時間。要解決這個問題,必須求解發(fā)電機的轉(zhuǎn)子運動方程。
(二)求解發(fā)電機的轉(zhuǎn)子運動方程
求解發(fā)電機轉(zhuǎn)子運動方程可以得出δ-t和ω-t的關系曲線。其中δ-t曲線一般稱為搖擺曲線。在上述簡單電力系統(tǒng)中故障期間的轉(zhuǎn)子運動方程為
(3-16)
式中,——功率角,其單位為弧度;——轉(zhuǎn)子角速度,標幺值;——轉(zhuǎn)子的同步角速度,即==314.16,其單位為弧度/秒;——發(fā)電機的慣性時間常數(shù),其
10、單位為秒;、――分別為機械和電磁功率,標幺值。
這是兩個一階的非線性常微分方程,它的起始條件是已知的,即
==0; ==1.0;==
故障切除后,由于系統(tǒng)參數(shù)改變,以致發(fā)電機功率特性發(fā)生變化,必須開始求解另一組微分方程:
(3-17)
式中變量含義同前述,其中也為標幺值。這組方程的起始條件為
其中為給定的切除時間;、為與時刻對應的和,它們可由故障期間的δ-t和ω-t的關系曲線求得(和都是不突變的)。一般來說,在計算故障發(fā)生后幾秒種的過程中,如果δ始終不超過180o,而且振蕩幅值越來越小,則系統(tǒng)是暫態(tài)穩(wěn)定的。
當發(fā)電機與無限大系統(tǒng)之間發(fā)生振蕩或失去同步時,在
11、發(fā)電機的轉(zhuǎn)子回路中,特別是阻尼繞組中將有感應電流而形成阻尼轉(zhuǎn)矩(也稱為異步轉(zhuǎn)矩)。當作微小振蕩時,阻尼功率可表達為:
== (3-18)
式中,稱為阻尼功率系數(shù);為轉(zhuǎn)子角速度的偏移量,標幺值;為轉(zhuǎn)子角速度,標幺值。阻尼功率系數(shù)除了與發(fā)電機的參數(shù)有關外,還和原始功角、的振蕩頻率有關。在一般情況下它是正數(shù)。在原始功角較小,或者定子回路中有串聯(lián)電容使定子回路總電阻相對于總電抗較大時,D可能為負數(shù)。如果考慮阻尼功率的影響,則故障后的轉(zhuǎn)子運動方程又可表達為
(3-19)
電力系統(tǒng)暫態(tài)穩(wěn)定計算包括兩類問題,一類是應用數(shù)值計算法得出故障期間的曲線
12、后,根據(jù)曲線找到與極限切除角對應的極限切除時間,此時只需要求解微分方程(3-16);另一類是已知故障切除時間,需要求出搖擺曲線來判斷系統(tǒng)的穩(wěn)定性,此時需要分段分別求解微分方程(3-16)和(3-17)。如果考慮阻尼轉(zhuǎn)矩的影響,則此時需要分段分別求解微分方程(3-16)和(3-19)。
三、 例題
例3-2 某簡單電力系統(tǒng)如圖3-4所示,取基準值=220MVA,=209KV。換算后的參數(shù)已經(jīng)標在圖中,其中一回線的電抗=0.486,=8.18秒。設電力線路某一回的始端發(fā)生兩相接地短路。假定=常數(shù)。(1)計算保持暫態(tài)穩(wěn)定而要求的極限切除角。(2)計算極限切除時間,并且作出在0.15秒切除故障時
13、的δ-t曲線。
圖3-4 某簡單電力系統(tǒng)的接線圖
解:計算系統(tǒng)正常運行方式,決定和。
由3-3(a)的正序網(wǎng)絡可得,此時系統(tǒng)的總電抗為
=0.295+0.138+0.243+0.122=0.798
發(fā)電機的暫態(tài)電勢為:
==1.41
==34.53o
(2)故障后的功率特性
又由3-3(b)的負序、零序網(wǎng)絡可得故障點的負序、零序等值電抗為
==0.222
==0.123
所以在正序網(wǎng)絡故障點上的附加電抗為:
于是故障時等值電路如圖3-3(c)所示,則
因此,故障期間發(fā)電機的最大功率為:
(3)故障切除后的功率特性
故障切除后的等值電路如圖3-3(d)所示
此時
14、最大功率為
圖3-5 例題7-12的等值電路
(a)正常運行等值電路;(b)負序和零序等值電路;(c)故障時等值電路;(d)故障切除后等值電路
(4)計算極限切除角
=0.458
(5)找出極限切除時間
根據(jù)(3-16),首先計算初值
令y(1)=,y(2)=。編寫描述故障期間轉(zhuǎn)子運動方程的ODE文件,文件名為ˊmyequˊ。
function dy = myequ(t,y)
dy = zeros(2,1);
f=50;w1=2*pi*f;
dy(1) = (y(2)-1)*w1;
dy(2) = (1/8.18)*(1.0-0.504*sin(y(1)));
再編
15、寫利用解題器指令求解y的程序。
clear
t0=0;tm=0.25;
d0=(34.53/180)*pi;w0=1;
[T,Y] = ode45('myequ',[t0 tm],[d0 w0]);
plot(T,(Y(:,1)/pi)*180,'-',0.194,62.76,'*')
text(0.194,60,'\delta_{cmax}=62.76\circ','FontSize',10)
text(0.194,56,'t_{cmax}=0.194s','FontSize',10)
圖3-6 例題7-12的δ-t曲線
圖3-6給出短路發(fā)生后0秒到0.25秒期間的δ-t
16、計算曲線,根據(jù)最大切除角()找到極限切除時間為0.194秒。由圖3-6可見,如果故障切除時間大于0.194秒,則發(fā)電機的功角將不斷地增大,最終失去暫態(tài)穩(wěn)定。在極限切除時間之前切除故障,發(fā)電機的搖擺曲線的狀況將在下面作計算、分析。
(6)不考慮阻尼轉(zhuǎn)矩的影響,當故障切除時間為0.15秒時通過計算得出δ-t曲線
首先編寫描述故障期間轉(zhuǎn)子運動方程的ODE文件,文件名為”myfun01”。
function dy = myfun01(t,y)
f=50; w1=2*pi*f;
TJ=8.18; Pt=1.0; P2m=0.504;
dy = zeros(2,1);
dy(1) =
17、(y(2)-1)*w1;
dy(2) = (1/TJ)*(Pt-P2m*sin(y(1)));
再編寫描述故障切除后轉(zhuǎn)子運動方程的ODE文件,文件名為”myfun02”。
function dy = myfun02(t,y)
f=50; w1=2*pi*f;
TJ=8.18; Pt=1.0; P3m=1.35;
dy = zeros(2,1);
dy(1) = (y(2)-1)*w1;
dy(2) = (1/TJ)*(Pt-P3m*sin(y(1)));
編寫利用解題器指令求解y的小程序。
clear
t0=0; tc=0.15; tm=2.0;
d0=(34
18、.53/180)*pi; w0=1.0;
[T1,Y1] = ode45('myfun01',[t0 tc],[d0 w0]);
dc=Y1(length(Y1),1);
wc=Y1(length(Y1),2);
[T2,Y2] = ode45('myfun02',[tc tm],[dc wc]);
plot(T1,(Y1(:,1)/pi)*180,'-',T2,(Y2(:,1)/pi)*180,'-',tc,(dc/pi)*180,'*')
text(0.28,50,'\it{t}_{c}=0.15s','FontSize',8)
text(0.28,43,'\it{\del
19、ta}_{c}=51.71\circ','FontSize',8)
xlabel('\it{t}')
ylabel('\it{\delta}')
計算結(jié)果表明,功角沿著故障切除后的功角特性曲線根據(jù)等面積定則作等幅振蕩,如圖3-7所示。實際上,由于阻尼轉(zhuǎn)矩的影響,振蕩的幅度是逐漸衰減的,功角最終運行在=47.8o。因此,發(fā)電機能夠保持暫態(tài)穩(wěn)定。
圖3-7 不考慮阻尼轉(zhuǎn)矩影響,當0.15秒切除故障時發(fā)電機的曲線
(7)考慮阻尼轉(zhuǎn)矩的影響,當故障切除時間為0.15秒時通過計算得出δ-t曲線
描述故障期間轉(zhuǎn)子運動方程的ODE文件與(6)相同,文件名也為”myfun01”。
重新編寫描述故
20、障切除后轉(zhuǎn)子運動方程的ODE文件,文件名為”myfun03”,阻尼功率系數(shù)D取為15。
function dy = myfun03(t,y)
f=50; w1=2*pi*f;
TJ=8.18; Pt=1.0; P3m=1.35;
D=15;
dy = zeros(2,1);
dy(1) = (y(2)-1)*w1;
dy(2) = (1/TJ)*(Pt-D*(y(2)-1)-P3m*sin(y(1)));
再編寫利用解題器指令求解y的小程序。
clear
t0=0; tc=0.15; tm=4;
d0=34.53*3.14/180; w0=1;
[T1,Y1
21、] = ode45('myfun01',[t0 tc],[d0 w0]);
dc=Y1(length(Y1),1);
wc=Y1(length(Y1),2);
[T2,Y2] = ode45('myfun03',[tc tm],[dc wc]);
plot(T1,(Y1(:,1)/pi)*180,'-',T2,(Y2(:,1)/pi)*180,'-',tc,(dc/pi)*180,'*')
text(0.3,50,'\it{t}_c=0.15s','FontSize',8)
text(0.28,43,'\it{\delta}_{c}=51.71\circ','FontSize',8)
xlabel('\it{t}')
ylabel('\it{\delta}')
圖3-8 不考慮阻尼轉(zhuǎn)矩影響,當0.15秒切除故障時發(fā)電機的曲線
計算結(jié)果表明,功角沿著故障切除后的功角特性曲線作減幅振蕩,如圖3-8所示。功角最終運行在=47.8o。因此,發(fā)電機能夠保持暫態(tài)穩(wěn)定。
【精品文檔】第 - 27 - 頁