劉永清，蔡金錠

（福州大學(xué) 電氣工程與自動化學(xué)院，福建 福州 350108）

目前，輸電線路的距離保護(hù)主要利用工頻極化相量和補償相量的相位差別來構(gòu)成保護(hù)原理。相量的計算通常采用全波傅立葉積分的算法來實現(xiàn)，可以濾除直流以及工頻的倍頻諧波分量，具有相對較好的效果。然而，當(dāng)故障發(fā)生后，故障電氣量中將包含大量衰減的直流和振蕩分量，這些衰減的暫態(tài)量是在故障暫態(tài)過程中產(chǎn)生的，嚴(yán)重影響相量的有效值以及相位角的計算精度，進(jìn)而導(dǎo)致保護(hù)的暫態(tài)超越。

為了防止暫態(tài)超越，很多文獻(xiàn)提出采用各種濾波方法來消除暫態(tài)分量的影響，然而衰減暫態(tài)分量的頻譜是連續(xù)的，不可能用濾波的方法將之徹底濾除［1］。這樣，不得不以犧牲速動性和靈敏性為代價，或者通過時間延遲躲避暫態(tài)分量，或者通過減少保護(hù)范圍保證不發(fā)生超越。當(dāng)在保護(hù)范圍末端附近故障時，保護(hù)動作時間超過40ms，嚴(yán)重影響系統(tǒng)的安全穩(wěn)定性［2-3］。因此，需要研究一種能夠在短時間窗內(nèi)（比如半個周期內(nèi)）準(zhǔn)確、快速地從衰減的暫態(tài)分量中提取工頻相量的計算方法，兼顧可靠性和速動性。

Prony算法能夠從含有大量振蕩衰減分量的故障信號中提取工頻分量，從本質(zhì)上說，Prony算法是通過2次線性估計將非線性的參數(shù)（頻率、相位、幅值和衰減常數(shù)）計算出來［4］，其思想是：將測量信號當(dāng)作某個線性系統(tǒng)的輸出，其振蕩頻率和衰減常數(shù)的信息與系統(tǒng)的輸入無關(guān)，只與系統(tǒng)本身的參數(shù)（線性系統(tǒng)的極點）有關(guān)，因此，可以用零輸入響應(yīng)作為方程，對其極點所在的多項式系數(shù)（z變換的分母多項式）進(jìn)行一次估計，然后再利用采樣值對幅值和相位信息進(jìn)行估計［5-7］。

自1795年由Baron de Prony提出的Prony算法以來，該算法在電力系統(tǒng)中得到了廣泛的應(yīng)用，但幾乎沒有應(yīng)用于距離保護(hù)，因為保護(hù)裝置很難實現(xiàn)2次參數(shù)估計的計算［8-10］（最小二乘法求解超定方程）。筆者對傳統(tǒng)的Prony算法進(jìn)行改造，精確估算工頻相量的幅值和相位。改造后的Prony算法只需要一次參數(shù)估計，其余的計算均為短窗的FIR濾波計算方法，大大減少了計算量，能在8ms內(nèi)準(zhǔn)確估算出工頻相量［11-12］。PSCAD仿真結(jié)果對該算法進(jìn)行了測試，測試結(jié)果表明，保護(hù)在8～10ms準(zhǔn)確動作而不受暫態(tài)分量的影響。

1 Prony算法基本原理

1.1 問題的提出

系統(tǒng)故障后的電氣量可以用多個幅值、頻率、相位和衰減的信號疊加表示，即衰減直流分量（P個衰減常數(shù)）、衰減高頻分量（Q個頻率和衰減常數(shù)）和工頻分量：

式（1）為模擬的信號；A0，A1k，A2k分別為工頻、衰減直流、衰減高頻分量的幅值；σk為衰減直流和衰減高頻分量的衰減因子；ω0和ωk分別為工頻和高頻分量的角頻率；ψ0和ψk分別為工頻和高頻分量的初始相位。

很顯然，式（1）中未知的參數(shù)是2×P+4×Q+3個，從理論上說，只要用多于參數(shù)個數(shù)的采樣點就可以通過解方程確定這些參數(shù)（最小二乘法參數(shù)估計）。但問題是，這是非線性方程組。Prony算法提供了一種利用線性方法求解式（1）非線性方程組的方法。

1.2 信號模型的歐拉公式表示

式中 N=p+2q+2；σk為各頻率分量的衰減因子；ωk為各頻率分量的角頻率（正、負(fù)頻率分別用不同的ωk來表示）；對于衰減的直流分量，Ck就是其幅值，而對于衰減的振蕩分量，Ck則是一個包含幅值和相位的復(fù)數(shù)：

1.3 信號的z變換

將方程式（2）進(jìn)行采樣，采樣間隔為Ts，采樣點方程為

對式（4）中的分式通分，可以發(fā)現(xiàn)信號的z變換是有理多項式的分式：

式中 a1～aN是分母的有理多項式系數(shù)。很顯然，分母多項式的系數(shù)是由極點值決定的，即僅與極點有關(guān)；分子也是有理多項式，但不必過多考慮（分子多項式系數(shù)是由Ck和分母zk共同決定）。

1.4 頻率和衰減因子的估計

信號的z變換（式5）可看作是某個線性系統(tǒng)的響應(yīng)。對于一個線性系統(tǒng)，其振蕩頻率和衰減因子是系統(tǒng)本身固有的，與輸入無關(guān)，因此，可以用零輸入方程來確定分母的多項式系數(shù)（從數(shù)學(xué)角度看，頻率和衰減因子也可以看作是差分方程的通解參數(shù)，零輸入響應(yīng)方程自然也是齊次差分方程）：

將式（6）反變換為時域形式可得到方程：

為了確定未知的參數(shù)，采用M個采樣值方程，且M＞N，即

用最小二乘法求解方程式（8）。

根據(jù)有理多項式系數(shù)｛a1，a2，…，aN｝可以確定頻率和衰減因子，因為zk=（σk+jωk）Ts為式（5）分母有理多項式等于零的根（極點）。

1.5 幅值和相位的估計

在求出式（6）中未知參數(shù)后，可根據(jù)式（3）用采樣值構(gòu)成M（M＞N）組方程組，求解Ck，若對應(yīng)的是復(fù)數(shù)，則Ck也是復(fù)數(shù)，即包含幅值和相位；若是實數(shù)，則Ck也是實數(shù)，只包含幅值信息。采樣值方程組為

這樣，通過求解方程（9）可計算出幅值和相位。

2 Prony算法的改造

通過上述分析，不難發(fā)現(xiàn)，Prony算法有2次最小二乘的參數(shù)估算、一次求解高次多項式方程，所有這些都不適合實時的距離保護(hù)的計算，因此，需對其進(jìn)行改造才能應(yīng)用于保護(hù)算法中。

首先，信號的z變換中肯定包含工頻極點，因此，無需再對工頻極點進(jìn)行計算；其次，可利用將式（5）消除其他高頻極點的方法，去掉衰減振蕩分量，而不需要再進(jìn)行求解高次多項式來確定高頻振蕩分量的頻率和衰減因子，同時，不需要再求解方程式（9）來確定各頻率分量的幅值和相位，可以直接獲得工頻量的相量值。

2.1 消除工頻極點

Prony算法中，極點的計算是通過求解齊次差分方程式（7）的有理多項式系數(shù)得到。所有的振蕩分量、衰減分量的頻率和衰減因子都包含在有理多項式等于零的極點中，當(dāng)然也包含工頻極點。因此，為了簡化計算，減少計算的階數(shù)以及降低由于計算誤差帶來的巨大影響，可將已知的工頻極點消除掉。消除工頻極點最簡單的方法就是將工頻極點從分母中移除，即將式（5）乘以一個工頻的零點：

2.2 消除衰減自由分量

經(jīng)過上述處理后，信號變?yōu)?/p>

式中 系數(shù)｛b1，b2，...，bN-2｝為濾除工頻極點后信號的z變換分母有理多項式系數(shù)。由式（11）不難發(fā)現(xiàn)，將原始信號經(jīng)過該式處理后，實際上是濾除了工頻分量，因為分母中的工頻極點已被消掉。

經(jīng)濾除工頻分量后，式（11）中只包含衰減的自由分量。同理，如果求出式（11）的系數(shù)｛b1，b2，...，bN-2｝，將原始信號X（z）乘上以該系數(shù)組成的有理多項式后，同樣也可以濾除衰減的自由分量，因為衰減自由分量的極點被消除后，其響應(yīng)自然不存在。

濾除工頻極點后，信號的修正模型可以表示為采樣值方程組的形式：

式中 系數(shù)｛b1，b2，...，bN-2｝為濾除工頻極點后信號的z變換分母有理多項式系數(shù)。其中包含除了工頻以外的衰減直流分量和高頻分量的極點，因此，將原信號X（z）乘以其非工頻分量的極點多項式，即可濾除衰減的非周期分量和衰減的高頻分量，濾除公式為

G（z）即為以衰減非周期分量的極點為零點的濾波器，系數(shù)｛b1，b2，...，bN-2｝為消除工頻極點后剩余極點的有理多項式系數(shù)。

2.3 工頻相量的計算

濾除非周期分量后，信號中只包含工頻分量，即

將式（14）轉(zhuǎn)化為時域形式：

其矩陣形式為

因此，工頻相量為

其中，Npc為每周期的采樣點數(shù)。

2.4 改造后的Prony算法

假設(shè)采樣的離散信號為｛x（k）｝，采樣頻率為每周期采樣Npc個點，經(jīng)過改造后的Prony算法流程：

1）濾除信號中的工頻分量，濾除公式為

2）利用式（12），計算衰減的自由分量的極點系數(shù)｛b1，b2，...，bN-2｝。N為Prony算法計算的階數(shù)，其值越大，擬合越精確，計算量當(dāng)然就越大。

3）利用上述極點系數(shù)濾除信號中衰減的自由分量，濾除公式為

4）利用式（17）計算信號中的工頻相量。

2.5 基于Prony算法的相間距離保護(hù)方案

距離保護(hù)算法需要構(gòu)造在線實時的保護(hù)方案，并應(yīng)綜合考慮計算精度和計算效率的協(xié)調(diào)問題。據(jù)此構(gòu)造的相間距離保護(hù)方案流程如圖1所示。

1）利用采樣中斷測量故障相間電壓和電流的瞬時值uφφ（k）和iφφ（k）；

當(dāng)ΔU＞0時，判斷為區(qū)外故障；當(dāng)ΔU＜0時，判斷為區(qū)內(nèi)故障。

圖1 相間距離保護(hù)方案流程Figure 1 Flow chart of interphase distance protection

3 仿真驗證

3.1 仿真模型建立

為了驗證該文提出的Prony算法在輸電線路相間距離保護(hù)中的有效性和準(zhǔn)確性，建立ATP仿真模型，如圖2所示，圖中各元件參數(shù)：500kV三相電源EM和EN的初始相位分別為0°，-20°；2個系統(tǒng)阻抗的正序電感值分別為39.8，35mH，忽略電阻和電容值；線路總長為100km，故障點F與M和N端的距離分別為80，20km；輸電線路采用分布參數(shù)模型，R0=0.184 7Ω／km，L0=3.601mH／km，C0=7.52nF／km，R1=0.017 4Ω／km，L1=0.967mH／km，C1=12.03nF／km。采樣頻率取1kHz（每周期20點），計算時間窗的寬度為半個周期（10點），階數(shù)N選擇5。

圖2 距離保護(hù)仿真模型Figure 2 Simulation model for distance protection

3.2 波形處理過程的驗證

故障后測量的電壓、電流波形如圖3所示，經(jīng)過工頻濾波器后的波形如圖4所示，圖4中的電流只剩下衰減的自由分量。僅對暫態(tài)量進(jìn)行極點的計算，不僅簡化了計算量，而且使得Prony算法的階數(shù)N的選取變得很容易，因為它不會影響到工頻極點的計算，只會影響計算的精確度（選擇模型階數(shù)N=5）。

利用計算出的暫態(tài)量極點系數(shù)，濾除衰減的自由分量，由于每一步計算的系數(shù)都不相同，因此無法具體分析其頻譜特性。濾除衰減自由分量后的波形與原始波形的比較如圖5所示，處理后的電壓、電流相量（幅值）隨時間變化曲線如圖6所示。

圖3 故障后的電壓和電流波形Figure 3 Fault voltage and current waveforms

圖4 濾除工頻量后的暫態(tài)量波形Figure 4 Transient components after 50Hz signal filtered

圖5 改進(jìn)Prony算法處理后的波形Figure 5 The waveforms after processing by improved Prony algorithm

圖6 改進(jìn)的Prony算法和傳統(tǒng)傅里葉算法處理后的電壓、電流波形比較Figure 6 Comparison of voltage and current by using improved Prony method and traditional Fourier

3.3 改造后的Prony算法的精度分析

假設(shè)在線路上發(fā)生A，B兩相接地短路故障，分別選取故障發(fā)生后M端的A相電壓和電流，應(yīng)用改進(jìn)的Prony算法對信號中的工頻分量進(jìn)行提取，分析比較其精度。在兩相接地短路故障趨于穩(wěn)定后，利用全波傅氏積分算法計算M端電壓、電流信號的有效值，作為Prony仿真結(jié)果的檢測標(biāo)準(zhǔn)。計算結(jié)果：U=336.4kV，θu=127.35°，I=6.89kA，θi=-179.97°。

3.3.1 不同時間窗的誤差

首先，根據(jù)式（3）用采樣值構(gòu)成M（M＞N）組方程組，然后，利用改造后的Prony算法解出式（9）距離保護(hù)算法的幅值和相位，計算結(jié)果與實際值的相對誤差如表1所示。通過分析可見，電壓相量的擬合誤差將隨著時間窗的縮短而增大，但半個周期（10ms）的誤差不超過3%，因此能滿足距離保護(hù)要求。對電流相量的擬合得到了相同的結(jié)果，篇幅所限，在此不再列出。

表1 不同時間窗時電壓工頻分量的Prony擬合結(jié)果Table 1 The Prony fitting results of voltage at different time window

3.2.2 不同故障位置的仿真情況

當(dāng)采用10ms時間窗時，在MN線路上不同位置發(fā)生兩相接地短路故障時的擬合結(jié)果如表2所示（表中“距離”百分值表示在線路總長的百分?jǐn)?shù)位置發(fā)生故障），可見在10ms時間窗下，MN線路不同位置發(fā)生短路故障時，利用改造后的Prony算法擬合的電壓幅值，其擬合誤差最大不超過1%，相位的擬合結(jié)果誤差最大不超過0.3%。因此，該文提出的算法的計算精度不受故障點位置不同的影響。

表2 不同故障位置時電壓工頻分量的Prony擬合結(jié)果Table 2 The Prony fitting results of voltage at different fault location

3.4 不同故障位置的保護(hù)動作情況

針對不同位置發(fā)生相間短路故障時驗證保護(hù)的動作情況，并與基于全波傅氏積分算法的結(jié)果進(jìn)行仿真對比。設(shè)線路全長為100km，整定距離為85km。

當(dāng)距離變電站10，80km處發(fā)生正方向區(qū)內(nèi)故障時，Prony算法與全波傅氏積分算法在故障發(fā)生后0.55ms時刻計算的比相值結(jié)果如圖7，8所示，圖中縱坐標(biāo)為相位比較后的余弦值，即將相位比較轉(zhuǎn)化為幅值比較并歸一化的值，原判據(jù)為相位比較：

轉(zhuǎn)化為幅值比較：

即當(dāng)幅值比較的值一旦小于零就進(jìn)入動作區(qū)，當(dāng)其等于-1時，為最大靈敏區(qū)域。

由圖7，8不難發(fā)現(xiàn)，距離保護(hù)安裝處10km的位置發(fā)生故障時，采用筆者提出的改造后的Prony算法的距離保護(hù)在故障（故障時刻為0.055s）后第8個點就進(jìn)入動作區(qū)（8ms）。而傳統(tǒng)的基于全波傅氏積分算法的距離保護(hù)則需要12個點才進(jìn)入動作區(qū)（12ms），而且該文提出的距離保護(hù)自8ms后就達(dá)到最大的靈敏度（相位比較為180°，即-1），而全波傅氏積分算法則需要一個周期才進(jìn)入最大靈敏度。

當(dāng)距離保護(hù)安裝處80km的位置發(fā)生故障時（接近整定點，整定距離為85km），筆者提出的距離保護(hù)算法在第10點（10ms）進(jìn)入動作區(qū)，而全波傅氏積分算法則在第17點（17ms）進(jìn)入動作區(qū)。通過圖5，6可以清楚地發(fā)現(xiàn)該文提出的距離保護(hù)算法的優(yōu)越性，不僅動作時間縮短，而且一旦進(jìn)入動作區(qū)即快速收斂到最靈敏的區(qū)域；從圖中還可看出，Prony算法的衰減速度快，達(dá)到穩(wěn)定的時間比較短。

圖7 10km故障處的比相算法比較Figure 7 The comparison of phase comparison algorithm at fault location of 10km from substation

圖8 80km故障處的比相算法比較Figure 8 The comparison of phase comparison algorithm at fault location of 80km from substation

筆者提出的改進(jìn)Prony算法是一種實時的算法，即每一個采樣間隔計算一個回合。在一回合計算中，共包括如下計算量：濾除工頻濾波器，3點相乘相加；參數(shù)估計計算（階數(shù)為5階），10個點矩陣計算；消除衰減分量的計算，5個點相乘相加；相量的計算，2個點相乘相加；總時間窗長度為20點，即若按照20點／周期采樣頻率，仿真結(jié)果表明，一般情況下，第8點（8ms）就進(jìn)入動作區(qū)。以上時間分析來自于基于Matlab下的仿真結(jié)果，并沒有將本算法實現(xiàn)到具體的裝置中。

現(xiàn)針對各種類型的故障，如不同故障位置的正向區(qū)內(nèi)、正向區(qū)外以及反向出口故障等進(jìn)行仿真，其動作情況如表3所示，可以看出，在不同位置發(fā)生相間短路時，筆者提出的相間距離保護(hù)能夠正確區(qū)分區(qū)內(nèi)、區(qū)外故障，準(zhǔn)確判斷不同故障位置的保護(hù)動作情況；所提出的距離保護(hù)判據(jù)在區(qū)內(nèi)任何位置故障時能夠在10ms左右準(zhǔn)確判斷，而全波傅氏積分算法則需要延遲14～17ms才能進(jìn)入動作區(qū)。

表3 不同故障處的比相算法比較Table 3 The comparison of phase comparison algorithm at different fault location

4 結(jié)語

利用Prony算法可在半個周期內(nèi)提取有效工頻信號，故障提取時間比全波傅氏積分算法縮短了半個周期，在進(jìn)行電力系統(tǒng)故障分析時具有足夠的可靠性和準(zhǔn)確性。筆者提出的基于改造后的Prony算法的距離保護(hù)與傳統(tǒng)的基于全波傅氏積分的保護(hù)算法相比，動作延遲時間較短，具有很高的穩(wěn)定性。

［1］張保會，尹項根.電力系統(tǒng)繼電保護(hù)（第二版）［M］.北京：中國電力出版社，2010.

［2］Zhang Q，Zhang Y，Song W，et al.Transmission line fault location for single-phase-to-earth fault on non-direct-ground neutral system［J］.IEEE Transations on Power Delivery，1998，13（4）：1 086-1 092.

［3］劉震.基于DSP的微機距離保護(hù)的研究與開發(fā)［D］.濟南：山東大學(xué)，2009.

［4］張新慧.基于Prony算法的小電流接地故障暫態(tài)選線技術(shù)［D］.濟南：山東大學(xué)，2008.

［5］Amono M，Watanabe M，Banjo M.Self-testing and self-tuning of power system stabilizers using prony analysis［J］.IEEE Transations on Power Delivery，1998，18（5）：38-47.

［6］常曉穎.基于Prony算法的電力系統(tǒng)故障信號快速擬合［D］.濟南：山東大學(xué)，2007.

［7］Trudnowski D J，Johnson J M，Hauer J F.Making Prony analysis more accurate using multiple signals［J］.IEEE Transations on Power System，1999，14（1）：226-231

［8］丁藍(lán)，薛安成，李津，等.基于窗口滑動改進(jìn)Prony算法的電力系統(tǒng)低頻振蕩識別［J］.電力系統(tǒng)自動化，2010，34（22）：24-28.DING Lan，XUE An-cheng，LI Jin，et al.A movingwindow Prony algorithm for power system low frequen-cy oscillation identification［J］.Automation of Electric PowerSystems，2010，34（22）：24-28.

［9］李大虎，曹一家.基于模糊濾波和Prony算法的低頻振蕩模式在線辨識［J］.電力系統(tǒng)自動化，2007，31（1）：14-19.LI Da-hu，CAO Yi-jia.An online identification method forpower system low frequency oscillation mode based on fuzzy filtering and Prony algorithm［J］.Automation of Electric Power Systems，2007，31（1）：14-19.

［10］吳旭升，馬偉明，王公寶，等.基于小波變換和Prony算法的同步電機參數(shù)辨識［J］.電力系統(tǒng)自動化，2003，27（19）：38-42.WU Xu-sheng，MA Wei-ming，WANG Gong-bao，et al.Parameter identification of synchronous machine based on wavelet transform and Prony algorithm［J］.Automation of Electric Power Systems，2003，27（19）：38-42.

［11］Tawfik M M，Morcos M M.On the use of Prony to locate faults in loop systems by utilizing modal parameters of fault current［J］.IEEE Transations on Power Delivery，2005，20（1）：532-534.

［12］Tawfik M M，Morcos M M.ANN-based techniques for estimating fault location on transmission lines using Prony method［J］.IEEE Transations on Power Delivery，2001，16（2）：219-224.