基于IRAM的大型電力系統(tǒng)特征值在線并行計(jì)算策略

王永昆，鮑顏紅，2，方勇杰，鄭 偉，徐 偉

（1.國(guó)網(wǎng)電力科學(xué)研究院/南京南瑞集團(tuán)公司，江蘇 南京 210003；2.華北電力大學(xué) 電氣與電子工程學(xué)院，北京 102206；3.中國(guó)電力科學(xué)研究院，北京 100192）

0 引言

區(qū)域電網(wǎng)互聯(lián)，特別是特高壓交直流電網(wǎng)的建設(shè)以及大規(guī)模風(fēng)電的接入使得電網(wǎng)的低頻振蕩問(wèn)題非常突出，進(jìn)行電力系統(tǒng)小干擾穩(wěn)定在線分析對(duì)預(yù)防互聯(lián)電網(wǎng)弱阻尼引起的低頻振蕩和連鎖故障具有重要的現(xiàn)實(shí)意義[1-2]。

自20世紀(jì)80年代以來(lái)，國(guó)內(nèi)外學(xué)者對(duì)大型電力系統(tǒng)特征值計(jì)算的方法進(jìn)行了大量的研究[3-5]，目前小干擾穩(wěn)定特征值離線計(jì)算方法已經(jīng)很成熟。但是基于離線數(shù)據(jù)的分析計(jì)算并不能反映系統(tǒng)時(shí)變的特性，為滿足在線小干擾分析[6-7]的要求，提出了并行計(jì)算和分布式計(jì)算的方法。文獻(xiàn)[8]介紹了BI迭代法、BI迭代法/逆迭代混合算法的并行化處理，該算法是求解一個(gè)特征值的并行實(shí)現(xiàn)方法，對(duì)原有串行代碼改動(dòng)較大，實(shí)際的加速比并不高。文獻(xiàn)[9]提出了按區(qū)域分網(wǎng)的分布式并行算法，該算法將每個(gè)區(qū)域設(shè)為一個(gè)子網(wǎng)，子網(wǎng)僅計(jì)算本區(qū)域數(shù)據(jù)，通過(guò)主控機(jī)的協(xié)調(diào)和交互來(lái)完成全網(wǎng)的小干擾穩(wěn)定分析，但是該算法也需要對(duì)原有串行算法進(jìn)行大量改動(dòng)，并且網(wǎng)絡(luò)劃分越不均，則各節(jié)點(diǎn)機(jī)相互等待時(shí)間越長(zhǎng)，計(jì)算效率也不高。文獻(xiàn)[10]將AESOPS算法與帶邊界分區(qū)的互聯(lián)電網(wǎng)切分方法相結(jié)合，提出了分布式的特征值算法，但該算法在迭代過(guò)程中需要多次交換數(shù)據(jù)，網(wǎng)絡(luò)通信的延遲是分布式計(jì)算主要的性能瓶頸。文獻(xiàn)[11]基于多核CPU并行計(jì)算技術(shù)提出了多進(jìn)程的特征值并行搜索算法，該算法將特征值搜索點(diǎn)根據(jù)CPU核的數(shù)量和設(shè)定的任務(wù)分配策略分配給各進(jìn)程，各CPU核的計(jì)算與串行計(jì)算完全相同，不但保持了串行程序的完整性，而且保證了計(jì)算效率。

本文在文獻(xiàn)[11]的理論基礎(chǔ)上，提出了一種基于位移求逆隱式重啟動(dòng)Arnoldi算法IRAM（Implicitly Restarted Arnoldi Method）的特征值在線并行計(jì)算策略，利用前一時(shí)刻斷面分析結(jié)果對(duì)特征值位移點(diǎn)（即搜索點(diǎn)）進(jìn)行了優(yōu)化選取，并且介紹了一種新的任務(wù)并行分配策略，可進(jìn)一步提高計(jì)算速度。

1 位移求逆IRAM

IRAM在收斂性能和計(jì)算速度上較優(yōu)越[12-13]，但一次位移求逆IRAM只能計(jì)算給定位移點(diǎn)附近一定數(shù)量的特征值，為求解搜索區(qū)域內(nèi)所有關(guān)鍵特征值，需要執(zhí)行多次位移求逆IRAM，以使搜索圓盤(pán)盡量覆蓋整個(gè)區(qū)域[14]（復(fù)平面上位移點(diǎn)附近的一個(gè)圓內(nèi)的特征值已被完全搜索）。目前位移點(diǎn)的選擇方法有中國(guó)電力科學(xué)研究院的PSD-SSAP①中國(guó)電力科學(xué)研究院系統(tǒng)所.PSD-SSAP電力系統(tǒng)小干擾穩(wěn)定性分析程序用戶手冊(cè).2009.采用的靜態(tài)設(shè)置方法和國(guó)網(wǎng)電力科學(xué)研究院的FASTEST②國(guó)網(wǎng)電力科學(xué)研究院穩(wěn)定所.FASTEST小擾動(dòng)穩(wěn)定分析技術(shù)手冊(cè).2007.及上海交通大學(xué)的小干擾穩(wěn)定分析包SSAP[15]采用的動(dòng)態(tài)設(shè)置方法。前者是對(duì)整個(gè)頻段等間隔地設(shè)置位移點(diǎn)，它的缺點(diǎn)是需憑經(jīng)驗(yàn)或反復(fù)嘗試來(lái)選擇位移點(diǎn)的數(shù)量，若位移點(diǎn)數(shù)量過(guò)少，將導(dǎo)致嚴(yán)重漏根現(xiàn)象，反之則會(huì)由于特征值的重復(fù)搜索而浪費(fèi)時(shí)間。后者是在頻段兩側(cè)取位移點(diǎn)，若任意2個(gè)圓盤(pán)不相交，則在中間點(diǎn)再取1個(gè)位移點(diǎn)，直至所有圓盤(pán)均相交，但該方法不易在并行計(jì)算中實(shí)現(xiàn)。

2 特征值在線并行計(jì)算策略

定義預(yù)估圓盤(pán)為復(fù)平面上以位移點(diǎn)為圓心、半徑r=ω/k的圓，如圖1所示。其中，ω為位移點(diǎn)的虛部，k為給定阻尼比范圍對(duì)應(yīng)的恒阻尼線斜率。

其中，ξ為給定阻尼比范圍邊界。

圖1 預(yù)估圓盤(pán)示意圖Fig.1 Schematic diagram of predict circle

系統(tǒng)結(jié)構(gòu)和運(yùn)行參數(shù)的變化對(duì)特征值的阻尼影響較大，對(duì)特征值的頻率影響較小[16]。假定系統(tǒng)運(yùn)行方式?jīng)]有發(fā)生大的改變，則可根據(jù)前一時(shí)刻斷面的特征值計(jì)算結(jié)果，得到指定位移點(diǎn)的預(yù)估圓盤(pán)內(nèi)的特征值數(shù)量，進(jìn)而為當(dāng)前時(shí)刻斷面計(jì)算時(shí)位移點(diǎn)及特征值數(shù)量的選取提供依據(jù)。

根據(jù)電力系統(tǒng)的特點(diǎn)，通常在小干擾穩(wěn)定分析中區(qū)內(nèi)振蕩模式較多，而區(qū)間振蕩模式較少[10]，因此中高頻段位移點(diǎn)的預(yù)估圓盤(pán)半徑與實(shí)際圓盤(pán)半徑相差不大。對(duì)于小型電網(wǎng)，特征值分布稀疏，預(yù)估圓盤(pán)與實(shí)際圓盤(pán)誤差較大，可能會(huì)導(dǎo)致重復(fù)搜索或遺漏關(guān)鍵根，效率降低。

2.1 主位移點(diǎn)的選擇策略

位移點(diǎn)的選擇和數(shù)量對(duì)并行算法的性能起著重要作用。為減少特征值的重復(fù)搜索，選擇的位移點(diǎn)應(yīng)滿足預(yù)估圓盤(pán)相交面積較少；為盡可能多地搜索到特征值，選擇的位移點(diǎn)應(yīng)滿足預(yù)估圓盤(pán)相交交點(diǎn)距復(fù)平面虛軸較遠(yuǎn)。已知第1個(gè)位移點(diǎn)位置及預(yù)估圓盤(pán)半徑，可將第2個(gè)位移點(diǎn)位置的選擇，轉(zhuǎn)化為求解兩圓盤(pán)相交面積極小且交點(diǎn)阻尼比絕對(duì)值極大的多目標(biāo)優(yōu)化問(wèn)題。

2.1.1 主位移點(diǎn)選擇的數(shù)學(xué)模型

如圖2所示，搜索區(qū)域內(nèi)兩相鄰位移點(diǎn)的搜索圓盤(pán)相交，其圓心（即位移點(diǎn)）分別為B和C，半徑分別為r1和r2，圓心距d=BC，交點(diǎn)A到兩圓心所在直線的距離h=AD，兩圓相交部分（即陰影部分）面積為S。

圖2 主位移點(diǎn)選擇示意圖Fig.2 Schematic diagram of main searching point selection

由圖2可建立如下等式：

在復(fù)平面中又有：

由式（1）—（3）消去 d 和 r2，可定義：

其中，S（h）為兩圓相交部分面積，是關(guān)于h的一元表達(dá)式，由式（1）—（3）消去 d 和 r2得來(lái)。

計(jì)算目標(biāo)是使 f1（h）盡量小，而 f2（h）盡量大。

權(quán)系數(shù)均取 0.5，且對(duì) f1（h）、f2（h）進(jìn)行無(wú)量綱化處理，根據(jù)線性加權(quán)和法建立評(píng)價(jià)函數(shù)如下：

求解 max f（h）最優(yōu)值 h，進(jìn)而求出 d，即可得到第2個(gè)位移點(diǎn)。

2.1.2 主位移點(diǎn)的選擇策略

根據(jù)第2.1.1節(jié)位移點(diǎn)選擇方法計(jì)算指定頻率及阻尼比范圍內(nèi)的各位移點(diǎn)和預(yù)估圓盤(pán)半徑，然后基于前一時(shí)刻斷面分析結(jié)果計(jì)算各預(yù)估圓盤(pán)內(nèi)特征值數(shù)量。

當(dāng)特征值數(shù)量較少時(shí)預(yù)估圓盤(pán)與實(shí)際圓盤(pán)會(huì)存在較大誤差，實(shí)際圓盤(pán)極有可能出現(xiàn)嚴(yán)重重疊，因此為了減少重復(fù)計(jì)算量，需要對(duì)位移點(diǎn)進(jìn)行修正，即合并預(yù)估圓盤(pán)內(nèi)特征值數(shù)量較少的位移點(diǎn)，取中點(diǎn)替代被合并的位移點(diǎn)，取預(yù)估特征值數(shù)量之和作為修正位移點(diǎn)的特征值數(shù)量。對(duì)于高頻段（0.7～2.0 Hz），特征值分布密集且搜索圓盤(pán)較大，建議合并2個(gè)位移點(diǎn)；而對(duì)于低頻段（0.1～0.7 Hz），特征值分布稀疏且搜索圓盤(pán)較小，可以考慮多合并1個(gè)位移點(diǎn)。由于預(yù)估圓盤(pán)與實(shí)際圓盤(pán)總會(huì)存在誤差，為防止漏根，建議采取保守計(jì)算，即對(duì)于主位移點(diǎn)實(shí)際計(jì)算的特征值數(shù)量n稍多于預(yù)估特征值數(shù)量，例如取n=+2。

2.2 輔助位移點(diǎn)的選擇策略

由于小干擾穩(wěn)定較關(guān)注的頻段范圍內(nèi)中高段附近振蕩模式較多，在線計(jì)算時(shí)預(yù)估圓盤(pán)半徑與實(shí)際圓盤(pán)半徑相差不大，可根據(jù)預(yù)估圓盤(pán)交點(diǎn)阻尼比與給定阻尼比來(lái)選擇輔助位移點(diǎn)。如圖3所示，輔助位移點(diǎn)取線段AF靠近F的黃金分割點(diǎn)E處作為輔助位移點(diǎn)，預(yù)估圓盤(pán)半徑為AE，對(duì)于輔助位移點(diǎn)的實(shí)際計(jì)算特征數(shù)取+3。

對(duì)于低頻段，由于特征值分布較少且每次計(jì)算時(shí)都進(jìn)行了保守計(jì)算，搜索圓盤(pán)可以完全覆蓋搜索區(qū)域，因此不必增設(shè)輔助位移點(diǎn)。

圖3 輔助位移點(diǎn)選擇示意圖Fig.3 Schematic diagram of assistant searching point selection

2.3 任務(wù)分配策略

在當(dāng)前的在線動(dòng)態(tài)安全分析系統(tǒng)中，通常利用分布式并行計(jì)算技術(shù)，采用同構(gòu)的計(jì)算節(jié)點(diǎn)組成計(jì)算集群，在集群計(jì)算平臺(tái)的管理節(jié)點(diǎn)上形成多個(gè)計(jì)算方案，分配到計(jì)算節(jié)點(diǎn)上進(jìn)行并行計(jì)算，然后在管理節(jié)點(diǎn)上匯總計(jì)算結(jié)果。

文獻(xiàn)[9]所提不均分配策略由于受計(jì)算任務(wù)隨機(jī)分配的影響，CPU利用率也具有很大的隨機(jī)性，本文將對(duì)其進(jìn)行改進(jìn)。

假設(shè)位移點(diǎn)數(shù)量為N，CPU核的數(shù)量為M，當(dāng)N?M時(shí)，參與并行計(jì)算的CPU核應(yīng)盡量多；當(dāng)N與M相差不大或N＜M時(shí)，考慮其他在線分析任務(wù)的需求，參與并行計(jì)算的CPU核的數(shù)量以小于N為宜。

IRAM的計(jì)算時(shí)間與位移點(diǎn)的位置和計(jì)算特征值的數(shù)量有關(guān)。可認(rèn)為每個(gè)位移點(diǎn)計(jì)算量的大小主要取決于每次計(jì)算特征值的數(shù)量，數(shù)量越多計(jì)算量越大；對(duì)于特征值數(shù)量相同的位移點(diǎn)，計(jì)算量大小則取決于位移點(diǎn)所處位置，所處位置頻率越低計(jì)算量越大。基于以上原則可估計(jì)出各位移點(diǎn)的計(jì)算量相對(duì)大小。

本文按預(yù)估計(jì)算量的大小分配任務(wù)，具體為：基于預(yù)估計(jì)算量將各任務(wù)壓入棧中，計(jì)算量大的位移點(diǎn)處于棧頂，計(jì)算量小的位移點(diǎn)處于棧底，首先根據(jù)CPU核的數(shù)量分配計(jì)算量大的位移點(diǎn)，剩余的位移點(diǎn)作為預(yù)留任務(wù)；當(dāng)某個(gè)CPU核的計(jì)算任務(wù)完成后，即到預(yù)留任務(wù)中按出棧順序自動(dòng)提取1個(gè)位移點(diǎn)進(jìn)行計(jì)算，直至所有任務(wù)計(jì)算完成。該分配策略可近似實(shí)現(xiàn)任務(wù)計(jì)算時(shí)間均分，最大限度地減少CPU核的閑置時(shí)間。

2.4 漏根判斷

由特征值計(jì)算結(jié)果很容易得出各位移點(diǎn)的實(shí)際圓盤(pán)半徑。若兩相鄰主位移點(diǎn)的圓盤(pán)不相交，則存在嚴(yán)重漏根現(xiàn)象，取兩相鄰主位移點(diǎn)的中點(diǎn)為新的位移點(diǎn)；若相交但交點(diǎn)不在臨近輔助位移點(diǎn)圓盤(pán)內(nèi)（高頻段）或交點(diǎn)阻尼比小于給定阻尼比（低頻段），則交點(diǎn)附近存在漏根，取圓盤(pán)交點(diǎn)為新的位移點(diǎn)。

2.5 計(jì)算步驟

特征值在線并行計(jì)算的整體計(jì)算流程如圖4所示。首次特征分布可采用現(xiàn)有的靜態(tài)設(shè)置位移點(diǎn)選擇方法和任務(wù)分配策略（平均分配或不平均分配策略）來(lái)得到或者采用離線分析結(jié)果。

圖4 特征值在線并行計(jì)算整體步驟Fig.4 Overall flowchart of online parallel eigenvalue calculation

3 仿真分析

測(cè)試算例采用南方電網(wǎng)2012年5月7日的在線數(shù)據(jù)：前一斷面時(shí)刻t1和當(dāng)前斷面時(shí)刻t2分別為23:49:34和23:54:35；搜索頻率范圍為0.1～2.0 Hz；阻尼比為 -0.1～0.1。

對(duì)t1時(shí)刻斷面進(jìn)行特征值分析，結(jié)果如圖5所示，圖中橫、縱軸分別為阻尼和振蕩角頻率，ε＝ξ。

圖5 t1時(shí)刻特征值分布圖Fig.5 Eigenvalue distribution at t1

根據(jù)第2節(jié)的位移點(diǎn)選擇策略選擇的位移點(diǎn)如表1所示，表中一對(duì)實(shí)部和虛部（均保留2位有效數(shù)字）組成1個(gè)位移點(diǎn)，n為對(duì)應(yīng)位移點(diǎn)處欲計(jì)算的特征值數(shù)量。

根據(jù)表1中數(shù)據(jù)對(duì)當(dāng)前（t2時(shí)刻）斷面進(jìn)行位移求逆IRAM特征值搜索計(jì)算，結(jié)果如圖6所示，為便于作圖，將虛軸和各圓盤(pán)整體進(jìn)行了縱向壓縮。

表1 位移點(diǎn)選擇情況Tab.1 Selection of searching points

圖6 t2時(shí)刻特征值分析結(jié)果Fig.6 Results of eigenvalue analysis at t2

從圖6可以看出，圓盤(pán)幾乎覆蓋整個(gè)小干擾穩(wěn)定所關(guān)注的頻段及阻尼比范圍，而且減少了特征值的重復(fù)搜索。僅在較高阻尼比且較高頻段內(nèi)會(huì)遺漏個(gè)別特征值（阻尼比分別為0.0957和0.0983），在實(shí)際中這些特征值是非關(guān)鍵的，而且如果搜索這些非關(guān)鍵特征值會(huì)耗費(fèi)更多時(shí)間。

表1中位移點(diǎn)的序號(hào)是按預(yù)估計(jì)算量由大到小的順序排列的，位移點(diǎn)的實(shí)際搜索時(shí)間如圖7所示。從圖中可以看出，預(yù)估計(jì)算量大的位移點(diǎn)計(jì)算耗時(shí)較大，預(yù)估計(jì)算量小的位移點(diǎn)計(jì)算耗時(shí)較小。由于特征值分布的不均導(dǎo)致某些位移點(diǎn)的預(yù)估計(jì)算量誤差較大，但對(duì)本文提出的任務(wù)分配策略影響不大。

圖7 各位移點(diǎn)的計(jì)算耗時(shí)Fig.7 Time consumption of different searching points

在相同計(jì)算環(huán)境下，通過(guò)計(jì)算機(jī)模擬雙CPU核下3種任務(wù)并行分配策略。忽略數(shù)據(jù)處理及數(shù)據(jù)輸出時(shí)間，各策略的計(jì)算時(shí)間取決于計(jì)算最慢的CPU核，如表2所示。從表2可以看出，本文所提分配策略耗時(shí)最短。任務(wù)分配的具體情況對(duì)平均分配策略和不均分配策略的計(jì)算時(shí)間影響很大，但這2種策略均缺乏有效的任務(wù)分配依據(jù)，只能隨機(jī)分配任務(wù)，而本文所提分配策略則能較好地解決任務(wù)分配不合理的問(wèn)題。

表2 不同策略計(jì)算時(shí)間對(duì)比Tab.2 Comparison of computing time among different strategies

同時(shí)基于本文所提分配策略也對(duì)四川電網(wǎng)和西北電網(wǎng)進(jìn)行了在線仿真分析，結(jié)果均驗(yàn)證了本文提出策略的正確性和有效性。

4 結(jié)論

本文基于位移求逆IRAM提出了一種電力系統(tǒng)小干擾特征值在線并行計(jì)算策略，該策略保持現(xiàn)有位移求逆IRAM特征值計(jì)算程序不變，充分運(yùn)用在線分析計(jì)算結(jié)果，高效、準(zhǔn)確地完成小干擾穩(wěn)定所關(guān)注的計(jì)算區(qū)域內(nèi)關(guān)鍵特征值的計(jì)算，程序?qū)崿F(xiàn)簡(jiǎn)單可靠，可為其他特征值算法的在線并行實(shí)現(xiàn)提供參考。