彭慧敏, 李 峰, 袁虎玲, 鮑顏紅, 李 威
(1. 國電南瑞科技股份有限公司, 江蘇省南京市 211106; 2. 南瑞集團(tuán)(國網(wǎng)電力科學(xué)研究院)有限公司, 江蘇省南京市 211106)
大電網(wǎng)穩(wěn)定計算過程中,發(fā)電機(jī)勵磁、電力系統(tǒng)穩(wěn)定器、靜止同步補(bǔ)償器(STATCOM)等無功功率—電壓控制模型,以及原動機(jī)和調(diào)速器等有功控制模型,由于控制系統(tǒng)結(jié)構(gòu)或控制參數(shù)的不合理,均可能導(dǎo)致大電網(wǎng)機(jī)電暫態(tài)數(shù)值仿真振蕩。同時,頻率和電壓大范圍變化下的直流和感應(yīng)電動機(jī)負(fù)荷模型,引起直流傳輸功率、無功消耗和感應(yīng)電動機(jī)負(fù)荷的等值注入電流急劇振蕩變化,也將導(dǎo)致大電網(wǎng)數(shù)值仿真振蕩。嚴(yán)重時甚至引起穩(wěn)定計算過程中數(shù)值穩(wěn)定性問題。
隨著大規(guī)模新能源并網(wǎng)、分布式發(fā)電和大量電力電子控制設(shè)備的應(yīng)用,大電網(wǎng)穩(wěn)定計算中的收斂性、振蕩等問題將更加趨于復(fù)雜,這給電網(wǎng)調(diào)度運(yùn)行和方式穩(wěn)定計算帶來了巨大的挑戰(zhàn)[1]。
電力系統(tǒng)機(jī)電暫態(tài)仿真求解的微分—代數(shù)方程組,通常采用顯式或隱式積分方法將微分方程差分,再與網(wǎng)絡(luò)代數(shù)方程交替或聯(lián)立求解。目前,國內(nèi)常用的電力系統(tǒng)分析綜合程序(PSASP)[2]、中國版BPA電力系統(tǒng)分析程序(PSD-BPA)[3]以及電力系統(tǒng)安全穩(wěn)定量化分析與優(yōu)化控制軟件(FASTEST)[4],微分方程采用隱式梯形積分法,并交替求解微分—代數(shù)方程組[5-7]。文獻(xiàn)[8]提出了一種新的動態(tài)仿真交替迭代方法,文獻(xiàn)[9]進(jìn)一步在含分頻風(fēng)電的電力系統(tǒng)中實現(xiàn)了機(jī)電暫態(tài)交替迭代求解,但交替求解的交接誤差問題仍不容忽視。
采用微分—代數(shù)方程交替求解時,當(dāng)系統(tǒng)電壓偏低,使得狀態(tài)量、代數(shù)量急劇變化時,迭代次數(shù)顯著增多,嚴(yán)重影響計算速度。為提高計算收斂性,目前仿真程序通常采用潮流穩(wěn)態(tài)結(jié)果初始化時域仿真狀態(tài)量,將直流輸電系統(tǒng)、光伏發(fā)電系統(tǒng)、STATCOM、靜止無功補(bǔ)償器(SVC)等電力電子設(shè)備以及負(fù)荷節(jié)點(diǎn)的初始有功功率和無功功率,折算成相應(yīng)大小的節(jié)點(diǎn)導(dǎo)納并入系統(tǒng)導(dǎo)納陣中,各仿真時步僅采用因節(jié)點(diǎn)電壓變化而引起的注入功率的變化量修正求解網(wǎng)絡(luò)代數(shù)方程。同時,在電壓偏低時,采用負(fù)荷模型轉(zhuǎn)化策略,將恒功率和恒電流負(fù)荷轉(zhuǎn)化為恒阻抗型進(jìn)行近似模擬。但是以上工程化實用的措施,仍無法解決大規(guī)模電力電子化設(shè)備引入系統(tǒng)后的仿真收斂性問題,更無法定位影響仿真收斂性的關(guān)鍵元件和原因。
本文首先由機(jī)電暫態(tài)仿真微分代數(shù)方程交替求解公式推導(dǎo)并分析引起仿真收斂性問題的基本要素,進(jìn)而選擇狀態(tài)量x中對應(yīng)發(fā)電機(jī)節(jié)點(diǎn)的勵磁電壓、機(jī)械功率、電磁功率和動態(tài)元件的等值注入電流,基于交替求解各積分時步迭代中上述變量的方差與均值之比,提出了發(fā)電機(jī)和動態(tài)負(fù)荷的仿真數(shù)值收斂性指標(biāo),自動識別影響仿真收斂性的關(guān)鍵動態(tài)元件,并從迭代初值改進(jìn)、積分步長自動調(diào)整、迭代合理逼近和迭代解合理選取方面,給出了提高仿真收斂性的方法。
微分代數(shù)方程組主要有聯(lián)立求解和交替求解兩種方法。聯(lián)立迭代求解時,由于系數(shù)矩陣為變系數(shù)的矩陣,計算量大,因此,電力系統(tǒng)大型暫態(tài)仿真軟件一般采用交替求解。交替求解的迭代格式如下:
(1)
式中:x表示狀態(tài)變量;V表示節(jié)點(diǎn)電壓;A為系統(tǒng)等值導(dǎo)納矩陣;I為節(jié)點(diǎn)等值注入電流;h為仿真步長;n和k分別為仿真積分時步數(shù)和每一積分時步的交替迭代次數(shù)。
(2)
(3)
(4)
(5)
由第1節(jié)交接誤差分析可知,正是由于暫態(tài)仿真交替求解前后兩積分步電壓變化較大,或者發(fā)電機(jī)組存在不平衡功率,使得顯式積分法預(yù)估的初值與真值存在明顯差異,導(dǎo)致暫態(tài)仿真數(shù)值收斂性變差,交替求解迭代次數(shù)明顯增加。
目前,大型機(jī)電暫態(tài)仿真程序常用迭代中最大注入電流偏差、加速功率偏差和兩次迭代間功角偏差綜合判別本輪迭代是否滿足收斂條件。但缺乏整個交替求解迭代過程中的數(shù)值收斂性自動識別指標(biāo)。事實上,本輪迭代偏差(注入電流或加速功率或功角)最大的節(jié)點(diǎn),并不一定后續(xù)各輪迭代偏差最大。同時,即使確定了迭代偏差最大節(jié)點(diǎn),仍無法明確收斂性差的原因。
暫態(tài)仿真迭代求解中等值注入電流反映了對應(yīng)動態(tài)元件微分方程組的收斂性,而發(fā)電機(jī)節(jié)點(diǎn)的勵磁電壓Efd、機(jī)械功率Pm、電磁功率Pe則分別反映了該機(jī)組勵磁系統(tǒng)、調(diào)速系統(tǒng)和轉(zhuǎn)子運(yùn)動的動態(tài)迭代收斂性。本文選擇狀態(tài)量x中分別對應(yīng)發(fā)電機(jī)的勵磁系統(tǒng)、調(diào)速系統(tǒng)和發(fā)電機(jī)輸出的勵磁電壓Efd、機(jī)械功率Pm、電磁功率Pe,以及所有對應(yīng)動態(tài)元件的等值注入電流輸出Ie,基于交替求解中各積分步迭代中上述變量的方差與均值之比,通過迭代求解中各控制系統(tǒng)的輸出收斂均值的離散程度,實現(xiàn)自動識別暫態(tài)仿真數(shù)值收斂性,以及制約收斂的主要動態(tài)元件和控制系統(tǒng)。
定義y=[Efd,Pm,Pe,Ie]T,交替迭代第i0次(暫態(tài)仿真交替迭代各積分時步正常收斂不超過i0次)至最大交替迭代求解次數(shù)N后y中第j個分量的均值和方差分別為:
(6)
(7)
Dmax(yn+1)的大小反映了該積分時步暫態(tài)仿真數(shù)值收斂的難易程度。
若各關(guān)鍵節(jié)點(diǎn)集對應(yīng)指標(biāo)均為Ie最大,且發(fā)電機(jī)節(jié)點(diǎn)指標(biāo)均為發(fā)電機(jī)輸出Pe偏大(相對于Efd和Pm),則網(wǎng)絡(luò)方程表現(xiàn)的交替誤差(電壓預(yù)估初值偏離真值)對數(shù)值收斂性影響明顯。
若關(guān)鍵節(jié)點(diǎn)集僅限為發(fā)電機(jī)節(jié)點(diǎn),則依據(jù)發(fā)電機(jī)的Efd,Pm,Pe對應(yīng)指標(biāo)大小確定影響因素。
依據(jù)Efd和Pm指標(biāo)大小,可以確定是動態(tài)無功支撐問題,還是有功輸出問題。
通常關(guān)鍵節(jié)點(diǎn)集中近區(qū)節(jié)點(diǎn)恒功率和恒電流占比大,或存在大規(guī)模直流或動態(tài)負(fù)荷,在電壓和頻率快速波動時,表現(xiàn)為數(shù)值收斂性差。
實際電力系統(tǒng)機(jī)電暫態(tài)仿真除網(wǎng)絡(luò)故障或操作時步外,一般每一積分時步迭代2至3次,最大交替迭代求解次數(shù)一般設(shè)為25次,即一般取i0=3,N=25可以快速識別仿真數(shù)值收斂性。
傳統(tǒng)的隱式積分方法中步長h是固定的,h值小,則相同的電壓修正量變化,引起的狀態(tài)量x的修正量變化小,仿真數(shù)值收斂性高。另一方面,較小的數(shù)值仿真步長h,雖然收斂性更好且一般仿真結(jié)果更精確,但將會顯著增加計算量和耗時。
適當(dāng)?shù)亟档头抡娣e分步長確實可以減小各次迭代交接誤差,減小各積分時步的迭代次數(shù),提高數(shù)值仿真收斂性[10-14]。以增加計算量為代價,導(dǎo)致相同的仿真觀察時段內(nèi)累計的迭代次數(shù)增加的步長降低,一般不可取。同時,何時啟動變步長,如何選擇變步長,何時恢復(fù)原積分步長,仍值得探討。
建議僅在識別數(shù)值收斂性差時(如網(wǎng)絡(luò)故障或操作階段),采用減小步長的變步長措施,即依據(jù)積分時步的暫態(tài)仿真數(shù)值收斂性指標(biāo)Dmax(yn+1)和迭代次數(shù)m,當(dāng)Dmax(yn+1)>Dmax,set且m大于設(shè)定次數(shù)mmax(如設(shè)Dmax,set=0.000 05,mmax=16)時,則自動減小步長。否則,數(shù)值收斂有所改善后自動恢復(fù)原積分步長。具體變步長過程如下。
若第n+1積分時步滿足Dmax(yn+1)>Dmax,set且m>mmax,則第n+2積分時步攝動減小步長為Δh(如0.001),繼續(xù)仿真n+2積分時步的迭代次數(shù)為M。當(dāng)M 通常Dmax,set和nmax設(shè)置不宜太小,否則增加計算量和耗時。實際工程中需兼顧兩者,考慮減小步長。 采用以上狀態(tài)量和節(jié)點(diǎn)電壓初值預(yù)估,一般每一積分時步交替迭代3至5次基本能滿足暫態(tài)仿真收斂性要求,但是在受擾期間或受擾之后出現(xiàn)電壓和頻率變化率較大,或者出現(xiàn)電壓和頻率嚴(yán)重偏離正常運(yùn)行值,使得直流等快速調(diào)節(jié)型電力電子型動態(tài)負(fù)荷、恒功率負(fù)荷的等值節(jié)點(diǎn)注入電流對節(jié)點(diǎn)電壓的靈敏度急劇增大,較小的節(jié)點(diǎn)電壓變化將引起較大的等值注入電流修正,導(dǎo)致交替迭代收斂性差。 為此,本文充分利用仿真程序已有計算信息,采用具有三階代數(shù)精度的Simpson求積計算式提供轉(zhuǎn)子運(yùn)動方程中發(fā)電機(jī)功角δ和頻差ω狀態(tài)量的初值,即發(fā)電機(jī)組j當(dāng)前時步n+1的功角δj,n+1和頻差ωj,n+1可以利用n-1時步的值復(fù)化求積求得,具體如式(8)和式(9)所示。 ωj,n+1= (8) (9) 式中:ΔPj,n為發(fā)電機(jī)組j第n步時的功率偏差;yj,n+1=ωj,n-1+4ωj,n;Tj為機(jī)組時間常數(shù);Dj為阻尼系數(shù)。 采用合理的初值預(yù)估,給出更加近似真值的預(yù)估值,可以明顯減小迭代誤差,提高數(shù)值仿真收斂性。不影響仿真結(jié)果的正確性。 實際機(jī)電暫態(tài)數(shù)值仿真中,數(shù)值仿真收斂問題更多出現(xiàn)在頻率和電壓大范圍變化下的動態(tài)負(fù)荷(直流和感應(yīng)電動機(jī)負(fù)荷模型)急劇振蕩變化,如故障瞬間和故障清除時刻對應(yīng)的網(wǎng)絡(luò)操作和變化,引起等值節(jié)點(diǎn)注入電流對節(jié)點(diǎn)電壓的靈敏度急劇增大,較小的節(jié)點(diǎn)電壓變化將引起較大的等值注入電流修正,導(dǎo)致交替迭代收斂性差,甚至出現(xiàn)數(shù)值收斂問題。 對于此類問題,一般可以通過暫態(tài)仿真數(shù)值收斂性識別方法,比較準(zhǔn)確地定位。此時,初值預(yù)估一方面困難,另一方面誤差較大,最好的解決方案是采用合理的迭代逼近策略提高數(shù)值收斂性,具體說明如下。 1)直流等快速調(diào)節(jié)型電力電子型動態(tài)負(fù)荷 對于直流等電力電子型動態(tài)負(fù)荷,為準(zhǔn)確模擬電力電子動態(tài)負(fù)荷的快速調(diào)節(jié)特性,機(jī)電暫態(tài)仿真時通常將其等值為可變動態(tài)負(fù)荷,而交流側(cè)接入點(diǎn)等值為可變電壓源,分別采用不同的迭代步長交替求解[16]。交流積分步內(nèi)相同的等效節(jié)點(diǎn)電壓變化,可變動態(tài)負(fù)荷側(cè)積分步長越小,電力電子的快速調(diào)節(jié)模擬越精確,引入的不同速率下的交替求解誤差越小,明顯提高了數(shù)值仿真收斂性。 2)恒功率負(fù)荷 目前,仿真程序中采用了靜特性負(fù)荷自動轉(zhuǎn)換為恒阻抗策略,根據(jù)設(shè)定的低電壓門檻值,自動轉(zhuǎn)換。本文建議識別因大型恒功率負(fù)荷節(jié)點(diǎn)存在引起數(shù)值收斂性問題后,為提高暫態(tài)仿真的收斂性,可依據(jù)設(shè)定的多輪低電壓門檻值啟動轉(zhuǎn)換,待電壓恢復(fù)后再恢復(fù)原靜態(tài)負(fù)荷特性。正常時采用較低的電壓門檻值(如0.5(標(biāo)幺值))啟動轉(zhuǎn)換,當(dāng)識別的各關(guān)鍵節(jié)點(diǎn)集對應(yīng)指標(biāo)均為Pe最大,且Dmax(yn+1)>Dmax,set,m>mmax,此時可以采用較高的電壓門檻值(如0.7(標(biāo)幺值))啟動轉(zhuǎn)換,以提高數(shù)值仿真收斂性。 工程中采用如下策略判別是否可能存在迭代達(dá)到Nmax并選取迭代解,然后給出仿真告警提示。 統(tǒng)計各迭代次數(shù)內(nèi)的最大網(wǎng)絡(luò)注入電流偏差和加速功率偏差之和Adpi,當(dāng)?shù)霈F(xiàn)Adpi增大發(fā)散,或者前后兩次迭代Adpi減小甚微且迭代次數(shù)達(dá)到Nmax/2,則不再進(jìn)行交替迭代求解,采用當(dāng)前迭代解值作為最接近真實解值,進(jìn)行下一積分時步迭代。 以上策略僅僅是工程實用的措施,通過近似選擇最小二乘解,可以避免迭代發(fā)散和Adpi趨于穩(wěn)態(tài)值的無效迭代。 合理的迭代逼近和合理的迭代解選取,均會影響仿真結(jié)果。如低電壓下靜特性負(fù)荷自動轉(zhuǎn)換策略,實質(zhì)上認(rèn)為電壓嚴(yán)重偏低時,繼續(xù)采用恒功率或恒電流模型仿真可能不合理,采用轉(zhuǎn)換恒阻抗方式(即等值的負(fù)荷電流存在最大值)繼續(xù)仿真模擬。 合理的迭代解選取,則改變迭代收斂可能達(dá)到最大迭代次數(shù)階段的仿真結(jié)果,以一種避免選擇迭代振蕩中的某次累計偏差更大的解的方式,在識別累計偏差開始振蕩時結(jié)束迭代,類似如該積分步放松迭代收斂條件,但對整個仿真階段結(jié)果影響較小,工程實用中是可取的。 將暫態(tài)仿真數(shù)值收斂性自動識別和提高調(diào)整嵌入到仿真計算流程中,可以實現(xiàn)大電網(wǎng)機(jī)電暫態(tài)仿真數(shù)值收斂性自動識別和提高調(diào)整(流程見附錄A)。 以IEEE 14系統(tǒng)數(shù)據(jù)為例,設(shè)置擾動場景為:t=0時支路9(首、末節(jié)點(diǎn)號分別為9和7)首端三相故障,t=0.1 s時切除支路9。t=1 s時,負(fù)荷節(jié)點(diǎn)2,3,4分別按每秒增加初始節(jié)點(diǎn)負(fù)荷的100%速率增負(fù)荷,負(fù)荷節(jié)點(diǎn)12,13,14分別按每秒減少初始節(jié)點(diǎn)負(fù)荷的100%速率減負(fù)荷,t=4 s時增減負(fù)荷均結(jié)束。仿真觀察時間20 s,暫態(tài)仿真各時刻的交替迭代次數(shù)如圖1所示,發(fā)電機(jī)頻率曲線(對應(yīng)節(jié)點(diǎn)8)和母線電壓曲線如圖2所示,其中母線電壓為標(biāo)幺值。 圖1 暫態(tài)仿真各積分時刻的交替迭代次數(shù)Fig.1 Number of alternate iterations at each integration step in transient simulation 圖2 暫態(tài)母線電壓和發(fā)電機(jī)頻率曲線Fig.2 Curves of transient bus voltage and generator frequency 由圖1和圖2可知,隨著模擬負(fù)荷的調(diào)整,發(fā)電機(jī)頻率降低,交替迭代次數(shù)明顯增加。當(dāng)頻率開始恢復(fù)時,交替迭代次數(shù)開始減少。隨著8 s后頻率第二次降低并恢復(fù),交替迭代次數(shù)再一次增加并趨于穩(wěn)定。 如圖3所示,采用本文提出的暫態(tài)仿真數(shù)值收斂性指標(biāo),自動識別本例數(shù)值收斂性指標(biāo)最大值Dmax(yn+1)對應(yīng)的狀態(tài)量為節(jié)點(diǎn)6的發(fā)電機(jī)注入電流,而引起迭代中注入電流方差較大的主要原因是電磁功率波動較大。圖3也給出了5臺分別連接于對應(yīng)節(jié)點(diǎn)(圖中數(shù)字標(biāo)示)的發(fā)電機(jī)電磁功率收斂性指標(biāo)曲線。 從圖3并對照圖1和圖2可知,采用暫態(tài)仿真數(shù)值收斂性指標(biāo)可以準(zhǔn)確定位影響仿真收斂性的地點(diǎn)(對應(yīng)節(jié)點(diǎn))、影響仿真收斂性的原因(對應(yīng)的狀態(tài)量),而且收斂性指標(biāo)與迭代次數(shù)以及頻率或電壓的變化波動正相關(guān),即某時段內(nèi)收斂性指標(biāo)越大,則交替迭代次數(shù)越多,對應(yīng)的頻率或電壓的變化越大。某時段內(nèi)收斂性指標(biāo)趨于穩(wěn)態(tài)非零值(或略有小幅波動),則交替迭代次數(shù)趨于穩(wěn)態(tài)值(或略有小幅波動),對應(yīng)的頻率或電壓趨于穩(wěn)態(tài)(或略有小幅波動)。 圖3 各積分時刻的暫態(tài)仿真數(shù)值收斂性指標(biāo)Fig.3 Numerical convergence index of transient simulation at each integration step 圖4 關(guān)鍵節(jié)點(diǎn)注入電流、電磁功率、勵磁電壓指標(biāo)Fig.4 Injection current, electromagnetic power and excitation voltage indices of key nodes 由圖4可知,關(guān)鍵節(jié)點(diǎn)集對應(yīng)指標(biāo)均為Ie最大,且發(fā)電機(jī)節(jié)點(diǎn)指標(biāo)均為Pe偏大,識別的關(guān)鍵節(jié)點(diǎn)均為青海省西寧、日月山、海西區(qū)的發(fā)電機(jī)節(jié)點(diǎn)。河武線故障潮流轉(zhuǎn)移至疆電東送二通道,導(dǎo)致青海750 kV線路母線電壓快速大幅下降,動態(tài)無功支撐不足,導(dǎo)致網(wǎng)絡(luò)方程求解數(shù)值收斂性差。 仍以上節(jié)西北電網(wǎng)某時段在線方式數(shù)據(jù)為例,圖5給出了采用本文方法前后故障清除前后0.3 s時段內(nèi)暫態(tài)仿真各步迭代次數(shù)比較。其中,在0~0.1 s故障期間,由于綜合采用了合理自動減小步長和合理選取的迭代解,整個迭代累積迭代次數(shù)明顯減小(若仿真步長選擇更小的值,單次迭代收斂提高了,但累計的迭代次數(shù)可能明顯增加)。 圖5 仿真各步迭代次數(shù)曲線Fig.5 Number of iterations at each simulation step 由圖5可知,采用本文方法在第一積分時步判別收斂性差啟動后,通過合理選取迭代解、改變發(fā)電機(jī)狀態(tài)變量預(yù)估方法以及短暫調(diào)整系統(tǒng)仿真步長,可以明顯改善后續(xù)故障期間各時步的數(shù)值收斂性,將交替迭代次數(shù)減小至10次以內(nèi)。整個仿真觀察時段內(nèi)減少迭代66次。 圖6給出了采用本文方法前后,發(fā)電機(jī)搖擺曲線差異最大的3臺機(jī)組功角差曲線。最大功角差異不超過0.25°,且隨著仿真收斂性改善,仿真收斂性提高措施自動退出,經(jīng)一段時間迭代后,兩者差異逐漸減小,并趨近于0。 圖6 采用本文方法前后3臺機(jī)組功角差曲線Fig.6 Power angle difference curves of three units before and after adopting the proposed method 本文分析了交替求解引起數(shù)值仿真收斂性問題的要素,基于各積分步中發(fā)電機(jī)狀態(tài)量和動態(tài)負(fù)荷的等值功率的迭代變化曲線方差與均值之比,提出了發(fā)電機(jī)和動態(tài)負(fù)荷的仿真數(shù)值收斂性指標(biāo)以及識別和提高仿真數(shù)值收斂性的方法,并將其嵌入現(xiàn)有大型仿真軟件中,自動識別影響仿真收斂的關(guān)鍵動態(tài)元件,依據(jù)識別的原因,自動選擇迭代初值改進(jìn)、積分步長自動調(diào)整等提高仿真收斂性措施和切換退出。仿真分析結(jié)果表明,本文提出的指標(biāo)和方法能準(zhǔn)確有效識別暫態(tài)仿真數(shù)值收斂性,通過適當(dāng)?shù)淖詣诱{(diào)整有助于收斂性改善。后續(xù)將從工程實用化角度進(jìn)一步完善。 附錄見本刊網(wǎng)絡(luò)版(http://www.aeps-info.com/aeps/ch/index.aspx)。3.2 狀態(tài)量和節(jié)點(diǎn)電壓初值預(yù)估
3.3 合理的迭代逼近
3.4 合理的迭代解選取
4 算例分析
4.1 暫態(tài)仿真數(shù)值收斂性自動識別實例
4.2 暫態(tài)仿真數(shù)值收斂性提高實例
5 結(jié)語