（合肥工業(yè)大學(xué) 電氣與自動化工程學(xué)院，合肥 230009）

0 引言

電力系統(tǒng)正常運行的破壞多半是由短路故障引起的。發(fā)生短路時，系統(tǒng)從一種狀態(tài)巨變到另一種狀態(tài)，并伴隨產(chǎn)生復(fù)雜的暫態(tài)現(xiàn)象[1,2]。狀態(tài)方程常常被用來研究電力系統(tǒng)暫態(tài)穩(wěn)定的動態(tài)特性，由于實際電力系統(tǒng)建立的狀態(tài)方程難以求得其解析解，通常要借助數(shù)值分析法。電力系統(tǒng)暫態(tài)穩(wěn)定計算常用的數(shù)值方法包括梯形法和Runge-Kutta法，而Runge-Kutta法要求方程的解具有足夠的光滑性，其迭代步數(shù)的增加會導(dǎo)致舍入誤差的積累。近年來，應(yīng)用于工程力學(xué)領(lǐng)域的微分求積法（Differential Quadrature Method，DQM），能用較少的數(shù)值點求得高精度的數(shù)值解，但目前還沒有基于DQM求解狀態(tài)方程的研究。本文給出了基于DQM求解狀態(tài)方程的算法，以電力系統(tǒng)常見的暫態(tài)過程—發(fā)電機發(fā)生三相短路為例，建立了系統(tǒng)的暫態(tài)穩(wěn)定狀態(tài)方程，將DQM和Runge-Kutta法用于求解系統(tǒng)的暫態(tài)穩(wěn)定狀態(tài)方程，通過測試及對比分析，驗證了DQM算法適用于電力系統(tǒng)暫態(tài)穩(wěn)定性的數(shù)值計算。

1 微分求積法

微分求積法的實質(zhì)是用全域上的全部節(jié)點的函數(shù)值進行加權(quán)求和來表示函數(shù)極其導(dǎo)數(shù)在給定點處的值，因而可以將微分方程變成以節(jié)點處的函數(shù)值為未知數(shù)的一組代數(shù)方程組[3,4]。不失一般性，考慮一維函數(shù)f(x)，它在區(qū)間[a,b]上連續(xù)可微，則有：

式中L是線性微分算子；Wj(x)是插值基函數(shù)，當(dāng)采用多項式插值時，它應(yīng)是N-1次多項式；xj為N個互異節(jié)點a=x1＜x2＜…＜XN=b中第j個節(jié)點的坐標(biāo)值。

我們稱Aij為函數(shù)f(x)一階導(dǎo)數(shù)的加權(quán)系數(shù)，稱[Aij]為一階導(dǎo)數(shù)的加權(quán)系數(shù)矩陣。

微分求積法權(quán)系數(shù)的顯式表達式的確定可通過拉格朗日（Lagrange）插值的方法得到[3,4]，即Aij可由下式計算：

根據(jù)插值函數(shù)各階導(dǎo)數(shù)的遞推關(guān)系，可得到二階及二階以上導(dǎo)數(shù)的權(quán)系數(shù)顯式表達式。

2 微分求積法求解狀態(tài)方程的誤差分析

設(shè)：

則狀態(tài)方程可表為：

其中，f1,f2,…,fn為n個狀態(tài)變量。

已知初值向量：[f1(0),f2(0),…,fn(0)]，應(yīng)用微分求積法求解區(qū)間[a,b]上的數(shù)值解。在求解區(qū)間[a,b]取N+1個節(jié)點，a=x0＜x1＜…＜xN=b，將狀態(tài)方程用微分求積的形式代替并移項整理，則導(dǎo)出的微分求積方程組為：

其中，Aij為DQM一階導(dǎo)數(shù)加權(quán)系數(shù)矩陣。

式(5)為微分求積法轉(zhuǎn)化后的線性方程組，其系數(shù)矩陣T為：

則式(5)可表為：

再設(shè)DQM所求的解向量為x1，誤差解向量為x2，由此可得：

由式(7)和式(8)可得解向量x1=T-1B，誤差解向量x2=-T-1b，于是有：

式(10)是狀態(tài)方程相對誤差的上界，它由兩部分組成，其中，為矩陣T的條件數(shù)，另一部分為B的相對誤差，的值與選取節(jié)點的方式有關(guān)，因此可通過節(jié)點的選取方式來減小的值使相對誤差減小。若系統(tǒng)的輸入為零，即求解系統(tǒng)的零輸入響應(yīng)時，則式(9)是其狀態(tài)方程誤差的上界。由此可見，對于一個確定的狀態(tài)方程而言，使用DQM法來求解時，可預(yù)先判斷轉(zhuǎn)化后線性方程組的系數(shù)矩陣的條件數(shù)是否過大，若過大，則對結(jié)果的相對誤差有較大影響，還可通過增加節(jié)點的方式來減小數(shù)值解的相對誤差。

3 電力系統(tǒng)暫態(tài)穩(wěn)定分析

電力系統(tǒng)暫態(tài)穩(wěn)定是指電力系統(tǒng)受到打擾動后，各同步發(fā)電機保持同步運行并過渡到新的或恢復(fù)到原來穩(wěn)態(tài)運行的能力，其計算分析的目的是在規(guī)定的運行方式和故障形態(tài)下，對系統(tǒng)的暫態(tài)穩(wěn)定性進行校驗，研究保證電網(wǎng)安全運行穩(wěn)定的控制策略，并對繼電保護裝置提出相應(yīng)的要求[5,6,7]。電力系統(tǒng)正常運行的破壞多半是由短路引起的。電力系統(tǒng)短路后的暫態(tài)分析一般是通過建立系統(tǒng)的狀態(tài)方程，然后用數(shù)值解法進行求解，計算出故障電流及其他故障變量[8]。

現(xiàn)以電力系統(tǒng)常見的暫態(tài)穩(wěn)定過程-同步電機發(fā)生三相短路為例進行分析。由于同步電機大量電感系數(shù)隨時間周期變化，其電壓方程為時變系數(shù)的非線性方程組，難以求得其解析解[9,10,11]。當(dāng)系統(tǒng)發(fā)生對稱故障時，易采用Park變換把坐標(biāo)變換到轉(zhuǎn)子側(cè)dq0坐標(biāo)中，而時變系數(shù)矩陣就變?yōu)槌?shù)矩陣。在求出dq0坐標(biāo)系中的電流分量后，再經(jīng)過Park逆變換求出a,b,c三相故障電流。由于本文分析的是系統(tǒng)的三相短路，其零軸分量為零，所以刪去了零軸分量之間的關(guān)系式。因突然短路而出現(xiàn)的各繞組電流變化量Δid、Δiq、Δif、ΔiD、 ΔiQ為狀態(tài)變量的關(guān)系式可表為：

式(11)可簡寫為：

在勵磁電壓不可調(diào)的前提下，以：

以及短路后最初瞬間的電流變化量：

代入，就可計算。

由以上分析可知，對于同步電機三相短路電流的計算，歸結(jié)為式(11)或式(12)狀態(tài)方程的初值問題。這類問題的數(shù)值計算方法有很多，例如改進歐拉法、Runge-Kutta法、預(yù)測-校正法等[3]。本文以DQM和Runge-Kutta法來進行求解。

求解狀態(tài)方程得到各繞組電流變化量與短路前各繞組電流疊加，得到電流全量：

其中，id[0]為定子d軸電流初始值，id[0]為定子q軸電流初始值，if[0]為勵磁繞組電流初始值。對其中的id、iq進行Park逆變換，即可得到定子三相繞組電流：

將它們代入式(15)：

4 微分求積法求解狀態(tài)方程的算法

考慮n維一階狀態(tài)方程組：

用微分求積法求解區(qū)間[a,b]上的數(shù)值解。在定義區(qū)間取N個節(jié)點，a=x1＜x2＜…＜xN=b將微分方程的導(dǎo)數(shù)用微分求積形式代替，則導(dǎo)出的微分求積方程組為：

其中A為微分求積一階導(dǎo)數(shù)的系數(shù)矩陣。將式(17)展開成線性方程組的形式并移項整理可得：

由式(18)可得所求解的線性方程組的系數(shù)矩陣TEM和B為：

其中TEM為nN×nN階矩陣，B為nN×1階矩陣。由于初值條件的存在，需對線性方程組的系數(shù)矩陣做修正。去掉由每一個初始微分方程形成線性方程組的第一個方程，且將與初值有關(guān)的項移到等式的右端，則TEM變?yōu)閚(N-1)×n(N-1)階矩陣，B變?yōu)閚(N-1)×1階矩陣，即：

5 微分求積法在同步電機三相短路計算中的應(yīng)用

由以上分析可知，同步電機三相短路計算歸結(jié)為對狀態(tài)方程的求解。將式(11)轉(zhuǎn)化為標(biāo)準(zhǔn)形式的狀態(tài)方程：

式(19)為標(biāo)準(zhǔn)形式的狀態(tài)方程，可由上述DQM求解狀態(tài)方程的算法進行求解。

為便于結(jié)果對比，以文獻[8]提供的實例進行分析。同步電機有如下參數(shù)：

在勵磁電壓不變的前提下，有：

原始運行條件為額定負(fù)載運行時：u[0]=1，i[0]=1，?[0]=0.5548 ，短路時的轉(zhuǎn)子位置角：0=3.1416，可求得：

原始運行條件為空載運行時：u[0]=1，i[0]=0.0，?[0]=0.0，短路時的轉(zhuǎn)子位置角：0=3.1416，可求得：

運用Runge-Kutta法和DQM求解各個時刻的Δidq0后，與idq0[0]進行疊加，即可得不同時刻的idq0，由式(14)、式(15)可求得定子a相電流ia和電磁轉(zhuǎn)矩T。

同步電機短路前額定負(fù)載和空載條件下的a相電流和轉(zhuǎn)矩的數(shù)值計算結(jié)果如圖1、圖2所示，計算結(jié)果為標(biāo)幺值。

圖1 額定負(fù)載下短路時的定子電流ia和轉(zhuǎn)矩T

圖2 空載下短路時的定子電流ia和轉(zhuǎn)矩T

由此可看出，DQM數(shù)值計算結(jié)果和等步長的Runge-Kutta法數(shù)值計算結(jié)果都與仿真結(jié)果基本一致，說明DQM可用于電力系統(tǒng)的暫態(tài)分析。

6 結(jié)論

本文提供了基于DQM求解狀態(tài)方程的初值算法，在求解狀態(tài)方程時，易于編程實現(xiàn)，并以電力系統(tǒng)常見的暫態(tài)過程-同步電機三相短路為例，通過本文提出的基于DQM的狀態(tài)方程初值算法來進行求解其短路電流與轉(zhuǎn)矩，得到的結(jié)果與仿真結(jié)果基本一致。得出結(jié)論，DQM數(shù)值解法可用于電力系統(tǒng)的暫態(tài)分析。

[1]葉瑞麗,劉瑞葉,劉建楠,等.直驅(qū)風(fēng)電機組風(fēng)電場接入后的電力系統(tǒng)暫態(tài)穩(wěn)定計算[J].電工技術(shù)學(xué)報,2014,06:211-218.

[2]仲悟之,湯涌.電力系統(tǒng)微分代數(shù)方程模型的暫態(tài)電壓穩(wěn)定性分析[J].中國電機工程學(xué)報,2010,25:10-16.

[3]汪芳宗,廖小兵,謝雄.微分求積法的特性及其改進[J].計算力學(xué)學(xué)報,2015,06:765-771.

[4]SS.A.EFTEKHARI,A. A. JAFARI. Mixed finite element and differential quadrature method for free and forced vibration and buckling analysis of rectangular plates[J].Applied Mathematics and Mechanics(English Edition),2012,01:81-98.

[5]劉彥呈,王川,魏一.電力系統(tǒng)大擾動下暫態(tài)電壓穩(wěn)定的動態(tài)特征分析[J].電力系統(tǒng)保護與控制,2011,10:11-17.

[6]趙晉泉,鄧暉,吳小辰,等.基于廣域響應(yīng)的電力系統(tǒng)暫態(tài)穩(wěn)定控制技術(shù)評述[J].電力系統(tǒng)保護與控制,2016,05:1-9.

[7]恰勒哈爾·吾肯,加瑪力汗·庫馬什.電力系統(tǒng)暫態(tài)穩(wěn)定影響因素的分析[J].電氣技術(shù),2012,10:5-9.

[8]陳珩.同步電機運行基本理論與計算機算法[M].北京:水利電力出版社,1992.

[9]王達達,張少泉,陳曉云,等.感應(yīng)發(fā)電機接入配網(wǎng)的三相短路電流峰值評估[J].電力系統(tǒng)保護與控制,2013,18:25-31.

[10]孫宇光,黃子果,陳麗,等.十二相整流同步發(fā)電機定子匝間短路故障計算[J].電工技術(shù)學(xué)報,2014,03:57-64.

[11]朱景偉,刁亮,任寶珠,等.具有冗余特性的永磁容錯電機短路故障分析與控制[J].電工技術(shù)學(xué)報,2013,03:80-86.