周博英，唐可翾，沈沉

（清華大學 電機工程與應用電子技術(shù)系，北京 100084）

含高比例可再生能源發(fā)電和高比例電力電子設備（“雙高”）正成為電力系統(tǒng)發(fā)展的重要趨勢和關(guān)鍵特征。“雙高”電力系統(tǒng)的動態(tài)特性中，隨機性和非線性特性進一步凸顯，其中隨機性越來越成為影響電力系統(tǒng)安全和穩(wěn)定運行的重要因素。傳統(tǒng)電力系統(tǒng)穩(wěn)定性往往使用確定性的時域仿真法對系統(tǒng)的微分方程進行數(shù)值分析，并基于概率論采用蒙特卡洛等分析方法對系統(tǒng)狀態(tài)進行重復模擬及統(tǒng)計分析。這種基于給定的系統(tǒng)模型和數(shù)據(jù)的方法很難全面刻畫系統(tǒng)中的隨機因素。

隨機響應分析側(cè)重于研究系統(tǒng)受動態(tài)時變隨機干擾的情況，通過注入隨機噪聲并對系統(tǒng)隨機響應統(tǒng)計特征進行分析，可在電磁暫態(tài)仿真中實現(xiàn)電力系統(tǒng)隨機穩(wěn)定性動態(tài)分析。電力系統(tǒng)隨機因素大多滿足維納過程（布朗運動），可通過高斯型噪聲進行建?！，F(xiàn)有的研究很少關(guān)注如何在電磁暫態(tài)仿真中隨仿真步長變化生成隨機噪聲，且缺少對噪聲注入后隨機響應仿真結(jié)果正確性的驗證。本文擬在電磁暫態(tài)仿真中實現(xiàn)對高斯型隨機噪聲的正確仿真，為電力系統(tǒng)基于隨機響應的動態(tài)分析提供支撐。

本文使用的電磁暫態(tài)仿真平臺為CloudPSS 平臺［1］，由清華大學和清華四川能源互聯(lián)網(wǎng)研究院共同研發(fā)，是一套面向新型電力系統(tǒng)、多時間尺度建模的高性能電磁暫態(tài)仿真云平臺。

1 電力系統(tǒng)高斯白噪聲建模

在電力系統(tǒng)研究中，通常假設隨機過程為高斯分布的白噪聲。對于受到隨機噪聲干擾的系統(tǒng)，如果注入的隨機噪聲干擾滿足維納過程，則在動力學分析時可以通過高斯白噪聲模型進行模擬［2］。白噪聲是一種最簡單的隨機過程，是由一系列不相關(guān)的隨機變量組成的理想化隨機過程。

電力系統(tǒng)的隨機因素大多滿足維納過程，連續(xù)隨機擾動因素大部分可以用高斯白噪聲隨機擾動進行描述［3］，例如系統(tǒng)諧波的處理、負荷的隨機波動、風電功率的隨機變化、機械功率的隨機扭振等，均可視為高斯白噪聲擾動［4］。

高斯白噪聲的標準矢量隨機微分方程如式（1）所示［5］：

其中，X(t)為n維矢量隨機變量，W(t)表示n維維納過程，W(t)的形式導數(shù)定義為W(t)/dt=η(t)，為高斯過程。本文高斯白噪聲統(tǒng)一用符號η表示。

線性系統(tǒng)高斯白噪聲隨機擾動可視作系統(tǒng)輸入，向系統(tǒng)的狀態(tài)方程和輸出方程進行注入。注入擾動后線性系統(tǒng)的隨機微分方程如式（2）所示：

其中，x為狀態(tài)變量，y為輸出變量，η為服從標準正態(tài)分布的高斯白噪聲，a，b，c，d為對應參數(shù)。

電力系統(tǒng)通常由高維非線性微分代數(shù)方程組描述，注入噪聲后的微分代數(shù)方程如式（3）所示：

發(fā)電側(cè)的高斯型隨機擾動主要指原動機隨機變化導致的發(fā)電機輸入機械功率的隨機波動［6］。以單機無窮大系統(tǒng)為例［7］，加入隨機擾動項后的系統(tǒng)方程如式（4）所示：

用電側(cè)負荷也具有較強的不確定性和波動性。在研究隨機擾動對電力系統(tǒng)電壓穩(wěn)定的影響時，需要研究電力系統(tǒng)無功負荷的隨機波動。在無功負荷中添加高斯白噪聲后，其定義如式（5）所示［8］：

此項無功負荷的波動會反映到系統(tǒng)電壓幅值和相角的方程中，導致系統(tǒng)電壓出現(xiàn)隨機波動。

電力系統(tǒng)在日常運行過程中，母線電壓、相角和線路功率等信號均存在類似噪聲的微小幅值波動。這類信號定義為類噪聲信號。若將負載變化看做完全隨機的過程，則類噪聲信號可看做白噪聲信號激勵電力系統(tǒng)產(chǎn)生的響應。通過對類噪聲信號進行隨機系統(tǒng)建模分析，有望估計出系統(tǒng)的動態(tài)模型［9］。

2 含高斯白噪聲微分方程一般求解

為探究電磁暫態(tài)仿真過程中高斯白噪聲的注入方法及其正確性，首先向只含一個狀態(tài)變量的一階RL電路中注入高斯白噪聲隨機擾動。含噪聲的系統(tǒng)微分方程如式（6）所示：

在電磁暫態(tài)仿真中，動態(tài)元件經(jīng)過諾頓等效可視為諾頓等效電導與諾頓等效歷史電流源并聯(lián)形式［10］，電路方程如式（7）所示：

其中，諾頓等效電導Geff=0.5ΔtL?1，諾頓等效歷史電流源

根據(jù)隨機微分方程求解法則，當積分步長Δt較小時，(Wt?Wt?Δt)可以用半階無窮小的高斯過程來代替。(Wt?Wt?Δt)為獨立變量，且與當前時刻電路狀態(tài)無關(guān)。得到遞推公式如式（8）所示：

其中，σN(0，1)為添加的隨機擾動項。

為了最小程度地改變仿真拓撲，隨機擾動項可以視作電流源iweff并聯(lián)于原電感元件，即可以將含有隨機擾動的電流源拆分成2 個元件，得到含高斯白噪聲的RL電路諾頓等效電路圖如圖1所示。

圖1 含高斯白噪聲的RL電路諾頓等效電路圖Fig.1 Norton equivalent circuit diagram of RL circuit with Gaussian white noise

該電路的約束方程如式（9）所示：

原本的電感支路iL0(t)記錄了原電感元件暫態(tài)電流，而iweff(t)記錄了每一時刻高斯白噪聲注入后隨機擾動等效電流源對系統(tǒng)的電流貢獻。按照該等效電路圖的思路，在CloudPSS平臺上對含高斯白噪聲的RL電路進行建模與仿真，仿真電路拓撲結(jié)構(gòu)如圖2所示。

圖2 CloudPSS平臺仿真含高斯白噪聲的RL電路Fig.2 CloudPSS simulation of RL circuit with Gaussian white noise

根據(jù)隨機分析的數(shù)學知識，一階RL簡單電路的隨機微分方程解析解如式（10）所示：

一階隨機微分方程的狀態(tài)變量中注入高斯白噪聲，iL(t)解析解也為高斯過程，服從正態(tài)分布，即

計算該高斯過程的方差及自協(xié)方差函數(shù)理論值，得到一階RL電路系統(tǒng)統(tǒng)計特征解析解如式（12）所示：

將仿真結(jié)果的方差和自協(xié)方差函數(shù)與理論值進行對比，結(jié)果如圖3 所示，驗證了該仿真方法在對含高斯白噪聲的RL電路進行仿真時的正確性。

圖3 RL電路統(tǒng)計特征仿真值與理論值對比Fig.3 Comparison between simulation value and theoretical value of RL circuit statistical characteristics

3 高斯白噪聲注入微分代數(shù)方程所引發(fā)的問題

沿用上節(jié)中的處理方法，在單機單負荷系統(tǒng)的負荷側(cè)有功功率中注入高斯白噪聲，即在這個電路中將高斯白噪聲注入代數(shù)變量中進行仿真。

在圖4 所采用的單機單負荷電路中，E=30 kV，R+jX=(0.1+1.6×10?6j)Ω，負荷數(shù)據(jù)為Pset=65 MW，Qset=20 MVar。在電力系統(tǒng)由微分代數(shù)方程描述的情況下，同樣用式（8）的半階無窮小方法對高斯白噪聲進行處理，對代數(shù)變量中注入隨機擾動后的系統(tǒng)隨機響應進行仿真。

圖4 含高斯白噪聲的單機單負荷電路拓撲圖Fig.4 Single machine single load circuit topology with Gaussian white noise

仿真預期的正確結(jié)果為變量的方差不隨仿真步長的變化而變化，但在CloudPSS仿真平臺上測試后發(fā)現(xiàn)（仿真結(jié)果如圖5 所示），無論是否進行半階無窮小處理，系統(tǒng)無功功率的方差均會隨著仿真步長的變化而變化。

圖5 無功功率方差隨仿真步長變化圖Fig.5 Reactive power variance changes with simulation step

仿真步長變化導致噪聲響應的性質(zhì)發(fā)生變化是不符合物理規(guī)律的。這一現(xiàn)象說明，第2 節(jié)中針對隨機微分方程的求解方法并不適用于隨機微分代數(shù)方程的求解。實際上，當系統(tǒng)由微分代數(shù)方程描述時，無論是否采用式（8）的半階無窮小處理方法，均會出現(xiàn)噪聲響應統(tǒng)計特征隨步長變化而變化的問題。因此，必須研究新的噪聲處理方法以求解微分方程與代數(shù)方程共同約束的系統(tǒng)。

4 高斯有色噪聲的構(gòu)造

為了解決含高斯白噪聲的微分代數(shù)方程求解問題，嘗試將注入的高斯白噪聲替換為高斯有色噪聲。事實上，現(xiàn)實中白噪聲過程并不能很好地描述自然界的各種擾動［11］，理想的白噪聲在物理上是很難實現(xiàn)的，并且實際工程中測量的數(shù)據(jù)所包含的噪聲往往是有色噪聲。使用PMU（Phase Measurement Unit）對實際電力系統(tǒng)進行數(shù)據(jù)采集時，會在理想高斯白噪聲的基礎上混雜一部分高斯有色噪聲［12］，這部分高斯有色噪聲是PMU 防混疊濾波器產(chǎn)生的［13］。同時，高斯白噪聲經(jīng)過低通濾波處理也將形成高度相關(guān)的高斯有色噪聲［14］。

高斯白噪聲的功率譜密度呈全頻域均勻分布，其自相關(guān)函數(shù)為δ沖激函數(shù)，方差為無窮大。而高斯有色噪聲的功率譜密度和自相關(guān)函數(shù)則如圖6所示［15］，其方差為有限值。不同于高斯白噪聲，有色噪聲序列中不同時刻的數(shù)據(jù)具有相關(guān)性。有色噪聲可視為白噪聲序列經(jīng)過線性環(huán)節(jié)的輸出，是理想白噪聲的一種近似。

圖6 高斯白噪聲與有色噪聲的功率譜密度及自相關(guān)函數(shù)Fig.6 Curves of power spectral density and autocorrelation function of Gaussian white noise and colored noise

考慮到電磁暫態(tài)仿真計算過程的單向遞推特性，高斯有色噪聲的構(gòu)造采取由高斯白噪聲插值的方法。圖7 展示了如何利用插值方法由高斯白噪聲構(gòu)造高斯有色噪聲。(dT)?1為高斯白噪聲注入的頻率，為保證噪聲的性質(zhì)相同，dT是不隨仿真步長變化而變化的，即在不同的仿真步長下，高斯白噪聲注入頻率相同。dt為實際仿真步長，則(dT/dt?1)即為構(gòu)造有色噪聲時每個dT時段中插值點的數(shù)量。隨著仿真步長的減小，插入數(shù)據(jù)的密度增加。

圖7 插值構(gòu)造高斯有色噪聲示意圖Fig.7 Diagram of Gaussian colored noise construction by interpolation

由這種方法構(gòu)造出的有色噪聲序列相關(guān)性與白噪聲略有差異，不同dT時段的數(shù)據(jù)之間沒有相關(guān)性，但每段dT中插入的數(shù)據(jù)點之間具有相關(guān)性。

為了方便在電磁暫態(tài)仿真中實現(xiàn)，選取階梯插值與分段線性插值2 種插值方法。經(jīng)插值后產(chǎn)生的高斯有色噪聲序列如圖8所示。

圖8 不同插值方法產(chǎn)生的高斯有色噪聲Fig.8 Gauss colored noise produced by different interpolation methods

圖9 不同插值點數(shù)的噪聲序列自相關(guān)函數(shù)Fig.9 Autocorrelation function of noise sequence with different interpolation points

圖10 不同插值點數(shù)的噪聲序列方差Fig.10 Variance of noise sequence with different interpolation points

圖11 不同仿真步長下輸入噪聲的自相關(guān)函數(shù)Fig.11 Autocorrelation function of input noise sequence with different simulation steps

圖12 不同插值點數(shù)的輸入噪聲方差Fig.12 Variance of input noise sequence with different interpolation points

圖13 不同插值點數(shù)的系統(tǒng)無功功率響應方差Fig.13 Variance of systemreactive power response with different interpolation points

圖8中，高斯白噪聲的注入頻率為5 × 105Hz，2 個相鄰的白噪聲數(shù)據(jù)中插值199 個數(shù)據(jù)點，即生成的高斯有色噪聲序列頻率為1×108Hz。

對不同插值方法產(chǎn)生的序列方差進行理論分析。由于階梯插值所插入的數(shù)據(jù)點均與原白噪聲相同，因此所產(chǎn)生的有色噪聲均值和方差應當均與原白噪聲相同。而對于線性插值產(chǎn)生的有色噪聲，首先考慮隨機變量X與Y相互獨立且均服從正態(tài)分布N（0，1）時，隨機變量aX+(1?a)Y的方差為=a2+(1?a)2。類似的，在服從正態(tài)分布N（0，1）的高斯白噪聲序列中的兩點(xi，yi)與(xi+1，yi+1)之間進行線性插值時，兩端的數(shù)據(jù)點可視為相互獨立的隨機變量，中間的插值數(shù)據(jù)點[kyi+(1?k)yi+1]，k∈(0，1)服從正態(tài)分布，其均值為0，方差為=k2+(1?k)2。因此，當插值點數(shù)為m時，每段dT中的全部數(shù)據(jù)點為[j(m+1)?1yi+(m+1?j)(m+1)?1yi+1]，j=0，1，2，3，…，m。經(jīng)過必要的近似處理（如忽略樣本方差與整體方差的差別及忽略首末端數(shù)據(jù)的影響等）之后，可以計算出線性插值有色噪聲序列的方差如式（13）所示：

計算結(jié)果表明，對于線性插值后的有色噪聲，當插值點數(shù)較大時，有色噪聲序列方差比原高斯白噪聲序列方差有所減小，并隨插值點數(shù)的增加趨近于固定值的2/3。

為驗證噪聲插值后統(tǒng)計特征的變化規(guī)律，在Matlab中對不同仿真步長（插值點數(shù)）下的有色噪聲統(tǒng)計特征進行數(shù)值驗證，結(jié)果如圖9-10所示。

在數(shù)值實驗中，為避免插值前噪聲序列的偶然隨機性對方差結(jié)果產(chǎn)生持續(xù)的影響，在不同的插值點數(shù)下需要重新生成正態(tài)分布的隨機噪聲序列。

Matlab 實驗結(jié)果表明，插值前后噪聲序列的統(tǒng)計特征變化符合預期，這為進一步在電磁暫態(tài)仿真中將有色噪聲作為系統(tǒng)輸入提供了正確性保證。

對于階梯插值，同一噪聲序列插值所得到的新序列與原白噪聲序列的自相關(guān)函數(shù)趨勢及方差大小均基本相同，且不隨插值點數(shù)量的變化而變化。

對于線性插值，在插值點數(shù)較大（≥20）時，插值后序列的統(tǒng)計特征基本不隨插值點數(shù)的變化而變化，且方差大致穩(wěn)定在原序列方差的2/3；在插值點數(shù)較小（<20）時，插值點數(shù)的變化對于插值后序列的方差及自相關(guān)函數(shù)的峰值有著較大的影響，插值后序列的方差可由式（13）進行定量計算。在實際仿真中使用線性插值方法得到有色噪聲序列時，若想要保持插值后序列的統(tǒng)計特征穩(wěn)定（不隨仿真步長變化而變化），則仿真頻率應比噪聲頻率高20 倍以上；若想要保持插值后序列的方差大小與原白噪聲序列基本相同，則可以基于式（13）的定量計算結(jié)果通過其他方式補償。

5 高斯有色噪聲的注入

使用2 種插值方法分別構(gòu)造有色噪聲，并在CloudPSS 平臺上分別作為隨機干擾注入到圖4 所示的單機單負荷系統(tǒng)的負荷有功功率中。電路參數(shù)與上文相同，隨機數(shù)據(jù)的注入頻率為1 × 104Hz。采用不同仿真步長（插值點數(shù)）進行仿真，每次仿真均重新生成高斯白噪聲序列，插值后輸入系統(tǒng)的有色噪聲序列和系統(tǒng)響應（負荷無功功率）的統(tǒng)計特征如圖11-13所示。

由圖11-13 可見，仿真結(jié)果的理論計算結(jié)果和數(shù)值實驗結(jié)果基本相符。

對于階梯插值，插值點數(shù)（仿真步長）不同時，插值后噪聲序列的自相關(guān)函數(shù)形態(tài)及方差大小基本不變，系統(tǒng)無功功率響應的方差在不同的插值點數(shù)下也基本保持不變。對于線性插值，在插值點數(shù)較大時，所得到的有色噪聲方差趨近于原高斯白噪聲方差的2/3，噪聲的性質(zhì)以及注入噪聲后系統(tǒng)響應的統(tǒng)計特征在不同的仿真步長下均基本保持不變；在插值點數(shù)較小時，仿真步長的變化對噪聲序列及系統(tǒng)響應的統(tǒng)計特征則有較大影響。仿真結(jié)果驗證了通過插值方法生成高斯有色噪聲并注入電磁暫態(tài)仿真的可行性與合理性。

6 結(jié)論

電力系統(tǒng)中的隨機擾動多為高斯型噪聲。在電磁暫態(tài)仿真中注入高斯型噪聲隨機擾動，需要保證系統(tǒng)的隨機響應統(tǒng)計特征不隨仿真步長變化而變化。在RL基礎電路中可以將高斯白噪聲隨機擾動等效為并聯(lián)半階無窮小電流源的形式加入系統(tǒng)，噪聲注入系統(tǒng)后的隨機響應統(tǒng)計特征與理論值相符。但該方法無法應用于求解含高斯白噪聲的微分-代數(shù)系統(tǒng)。采用階梯插值以及插值點數(shù)較大的線性插值方法所產(chǎn)生的有色噪聲統(tǒng)計特征在一定程度上與插值前的白噪聲相似，通過將注入的白噪聲替換為經(jīng)插值產(chǎn)生的有色噪聲，能夠規(guī)避隨機微分代數(shù)方程求解中所遇到的問題。在電磁暫態(tài)仿真中，通過向系統(tǒng)的微分代數(shù)方程組中注入本文方法所產(chǎn)生的隨機有色噪聲，可以實現(xiàn)對系統(tǒng)隨機響應統(tǒng)計特征的正確分析。