基于VMD 和改進聚類算法的配電網(wǎng)故障選線方法

王遠川，李澤文，夏翊翔，毛紫玲，郭欣玉

（長沙理工大學電氣與信息工程學院，長沙 410114）

目前我國中低壓配電網(wǎng)普遍采用小電流接地方式，大數(shù)據(jù)統(tǒng)計顯示，一年中單相接地故障發(fā)生的次數(shù)達到配電網(wǎng)所有發(fā)生故障的80%左右[1]，故障可能會引起系統(tǒng)局部電壓過高，進而對電氣設備造成破壞[2]，嚴重時甚至會危及人員的生命安全，因此需及時找出故障線路并予以切除。

基于主動注入法[3-5]和被動檢測法的故障選線方法是目前常用的兩種方法。通常，主動注入法易受過渡電阻的影響，且成本相對較高，適用性不強；被動檢測法可分為基于穩(wěn)態(tài)信息[6-7]和基于暫態(tài)信息的故障選線法[8-10]。實際運行中，諸多因素都會對傳統(tǒng)的基于穩(wěn)態(tài)信息的故障分析法造成干擾[11]，使得該方法很不可靠[12]。相比之下，基于暫態(tài)信息的故障分析法，憑借自身諸多優(yōu)點[13-14]，被廣泛應用于電力系統(tǒng)故障選線、行波定位等電網(wǎng)保護領域。許多故障特性，例如信號的相位、幅值、極性變化等，都存在于高頻暫態(tài)信號中[15-17]，如果能被有效地分解、提取出來，就能夠為選線提供思路。例如，文獻[18]采用小波分解將故障零序電流分解后，通過比較各線路在特征頻帶中的幅值和極性來鑒別故障線路。但在特定分解尺度下，小波變換不能同時在時間和頻率上都具有很高的精度，以及小波變換是非自適應性的，基函數(shù)選定了之后就不能在分解過程中改變，其選擇可能會對分解結果產生影響。針對這些問題，經(jīng)驗模態(tài)分解EMD（empirical mode decomposition）算法被運用在故障選線領域，其能依據(jù)數(shù)據(jù)自身的時間尺度特征來將信號分解為若干模態(tài)分量，如文獻[19]使用EMD算法對各線路導納進行分解，提取出特征分量，通過比較故障線路和非故障線路的特征分量符號進行選線。但是EMD 分解存在端點效應和頻率混疊等問題。此后，又相繼誕生了集合經(jīng)驗模態(tài)分解EEMD（ensemble empirical mode decomposition）和互補集合經(jīng)驗模態(tài)分解CEEMD（complementary ensemble empirical mode decomposition）等改進型算法，但始終存在分解效果不理想的問題。

有一類利用形態(tài)學選線的方法，其對分解后的零序電流采用相關性聚類算法來鑒別故障線路，如文獻[20]使用凝聚層次聚類HAC（hierarchical agglomerative clustering）算法來計算不同線路零序電流相似性測度并聚類，但HAC 算法時間復雜度高，且具有貪心算法的缺點。又如，文獻[21]使用的基于劃分算法的K-means 聚類算法具有步驟簡單、時間復雜度、空間復雜度低的優(yōu)點，但該算法易受初始聚類中心選取的影響[22]，且運算效率不及在其基礎之上改進的一系列算法。

針對上述情況，為了高效、準確地選出故障線路，本文提出了一種基于變分模態(tài)分解VMD 和Kmeans++聚類算法的故障選線方法。首先對提取到的零序電流信號進行VMD 分解，得到有限個自適應頻帶特征的本征模態(tài)函數(shù)IMF（intrinsic mode function）；而后通過對其中的低頻分量作相關性計算，對高頻分量測算其初始極性，從而確定代表各出線零序電流信號的散點在二維平面內的橫、縱坐標；最后通過K-means++聚類算法檢測點集的分布情況，將故障線路篩選出來。

1 故障零序電流的分析計算

1.1 小電流接地系統(tǒng)解析

配電網(wǎng)發(fā)生單相接地故障示意如圖1 所示。圖1 中，rp、Lp分別為消弧線圈的等效電阻和等效電感值；Rf為接地過渡電阻；C0j為第j條線路Lj的對地電容大小；i0j為線路Lj上流過的零序電流；iL為消弧線圈上流過的感性電流；U0f為零序電壓源。

圖1 零序電流分布Fig.1 Zero-sequence current distribution

由圖1 可知，發(fā)生單相接地故障后，系統(tǒng)中各條健全線路首端流過其本身產生的容性電流，這些電流都將經(jīng)過接地點流入故障線路中。因此，零序電流在故障和健全線路中的流向剛好相反，在非故障線路中為出線端向線路末端，而在故障線路中則為線路末端向出線端，這也導致故障線路中的容性電流與健全線路中的容性電流初始極性理論上也是相反的。值得注意的是，消弧線圈中此時將產生一個與接地點處容性電流相位相反的感性電流，該電流經(jīng)電感線圈注入大地再經(jīng)接地故障點流入故障線路中。

1.2 故障零序電流的分析

建立單相接地故障時的暫態(tài)分析等值電路模型，該電路如圖2所示。

圖2 系統(tǒng)暫態(tài)分析等值電路Fig.2 Equivalent circuit of transient analysis of system

圖2 中，R0、L0分別為配電網(wǎng)的等值電阻、電感；C∑為系統(tǒng)零序電容之和；if為總故障電流；ic為暫態(tài)電容電流。

根據(jù)基爾霍夫電壓定理列寫圖2 中等效電路回路的微分方程為

式中：φ為故障時刻電壓初始相角；Um為系統(tǒng)相對地的電壓幅值；ω為工頻角頻率；t為時間。

對上述微分方程求解，可得故障點處流過的容性電流和感性電流的計算公式為

式中：、分別為容性電流的高頻振蕩與工頻交流分量；為感性電流的工頻交流分量；iLdc為感性電流的衰減直流分量；Icm、ILm分別為容性電流、感性電流的幅值大小；ωf為暫態(tài)自由振蕩分量的角頻率；β為自由振蕩分量的衰減系數(shù)；τL為回路的時間常數(shù)。

容性電流和感性電流相加可得故障點處流過的總電流if，即

由于衰減直流分量只存在于感性電流中，因此其只會隨著感性電流從接地故障點流入故障線路。根據(jù)接地故障點電流表達式可以看出，故障零序電流的特征不僅與配電線路本身一系列固有參數(shù)相關，其各組成部分的幅值也會隨著故障初相角取值大小的改變而變化。根據(jù)容性電流高頻振蕩分量的表達式可知，其波形振幅很大程度上取決于sinφ的大小；而對于感性電流中的衰減直流分量，其波形振幅大小取決于cosφ的大小。因此當故障電流直流分量數(shù)值較大時，高頻振蕩分量一般較小；反之，當高頻振蕩分量數(shù)值較大時，衰減直流分量一般較小。

例如，當φ=0°時，故障電流表達式可化簡為

從式（6）可以看出，接地故障點處零序電流含有全部的3 個分量，且衰減直流分量幅值較大，而高頻分量幅值較小。

當φ=90°時，故障電流表達式可化簡為

從式（7）可以看出，此時故障零序電流主要由高頻分量與工頻分量所組成，且高頻分量有著較大的幅值，而幾乎不含有直流分量。

2 變分模態(tài)分解基本原理

VMD算法具有自適應和完全非遞歸性[23]，其將信號的分解問題轉化為構建變分問題并求解，尋找變分問題的最優(yōu)解，保證分解后得到的序列為具有中心頻率和有限帶寬的模態(tài)分量，實現(xiàn)各頻率下本征模態(tài)函數(shù)IMF（intrinsic mode function）的有效分離，得到預設分級層數(shù)下的K個模態(tài)分量。其核心步驟由變分約束問題的構造和求解構成。

首先，構造變分約束問題，其數(shù)學表達式為

式中：f為待分解的原始信號；?t為梯度運算；δ(t)為單位脈沖；j為虛數(shù)單位；{uk} ={u1,u2,…,uK} 為分解后得到的K個本征模態(tài)函數(shù)；{ωk} ={ω1,ω2,…,ωK} 為各模態(tài)分量的中心頻率。通過引入懲罰因子α和拉格朗日算子λ，構造增廣拉格朗日表達式，即

通過使用交替乘法算子迭代、計算式（9），得到其鞍點即為變分約束問題的最優(yōu)解。各個模態(tài)分量和其中心頻率的迭代公式為

如果滿足，則將目前得到的各IMF作為輸出量輸出；若不滿足，則回到式（10）和（11）繼續(xù)進行迭代更新。

3 選線相關算法原理及判據(jù)構造

結合第1.2節(jié)對發(fā)生單相接地故障時故障零序電流組成及其相關特征的解析可知，含有衰減直流分量的感性電流只流入故障線路中，這將使故障發(fā)生初期故障線路零序電流低頻分量的幅值遠大于健全線路，從而產生波形上的差異；又由于零序電流在故障和健全線路上的流向相反，故障線路與健全線路零序電流信號高頻分量的初始極性也應當是相反的，因此無論是低頻分量還是高頻暫態(tài)分量，故障線路相比于非故障線路均存在著明顯的差異，如果能將這種差異投射到二維平面或者三維空間中，將有利于篩選出故障線路。

第1.2 節(jié)提到，故障電流衰減直流分量數(shù)值較大時，高頻振蕩分量往往較小；反之當高頻振蕩分量較大時，衰減直流分量往往較小。換言之，一方差異較為明顯時，另一方差異有可能難以辨識。因此，本文選擇綜合考慮二者的差異，將其投射到二維平面中，通過對電流信號低頻分量相關系數(shù)與高頻分量初始極性差異構造的二維平面點集進行聚類，從而準確選出故障線路。

3.1 VMD 分解相關參數(shù)取值的確定

對采樣信號做VMD 分解之前，需要提前設定好分解層數(shù)K、懲罰因子α以及預設精度ε[24-25]。

故障零序電流信號經(jīng)過VMD 分解后將得到K個子模態(tài)分量，每個子模態(tài)都有各自的中心頻率。VMD算法與EMD等其他分解方法的一個很大的不同點在于，其需要人為預先設定K的取值。

零序電流信號中的高頻振蕩分量的振蕩頻率浮動區(qū)間大致為1 000～3 000 Hz，與工頻分量和直流分量相差巨大；而工頻分量與直流分量之間頻率相差不大。考慮到實際采集的信號往往含有噪聲，因此本文所有的仿真信號都將注入25 dB的噪聲作為干擾。

本文最終采取的分解方式為：先使用FIR 濾波器濾去頻率過高的噪聲信號（在仿真實驗時取5 500 Hz 作為濾波界線），而后截取故障發(fā)生時刻后2 個信號周期（40 ms）的零序電流信號先做2 層分解，得到的模態(tài)分量分別記為IMF1和IMF2。高頻振蕩分量與過濾后剩下的噪聲信號將全部分解到IMF2中，而衰減直流、工頻交流分量和低次諧波則全部存在于IMF1中。將得到的IMF1分量再次做2層分解，得到的模態(tài)分量分別記為IMF1-1與IMF1-2，這樣衰減直流分量和工頻分量被全部分解到IMF1-1中，低次諧波分量則被分解到了IMF1-2中。由此各模態(tài)分量根據(jù)頻率由低到高（IMF1-1、IMF1-2、IMF2）的主要成分如表1所示。

表1 各模態(tài)分量的主要成分Tab.1 Main components of each modal component

懲罰因子α的取值會對模態(tài)分量的帶寬產生影響，研究表明α值的大小與模態(tài)分量的帶寬成反比例關系。α的值越小，信號帶寬越大，各模態(tài)分量的函數(shù)會有很大一部分頻率重疊部分，即模態(tài)混疊；α的值越大，信號帶寬越窄，模態(tài)分量向中心頻率收斂速度越快，故障波形中的重要信息存在丟失的可能性。根據(jù)以往的VMD分解經(jīng)驗，α的取值范圍通常在2 000～4 000。本文經(jīng)過多次實驗，結合電網(wǎng)實際運行中的工況，選擇α的取值為4 000，預設精度ε的值取為1×10-6。圖3 是某條件下的故障線路信號波形，其分解結果如圖4所示。

圖3 某條件下的故障線路信號波形Fig.3 Fault line signal waveform under certain conditions

圖4 信號分解結果Fig.4 Signal decomposition results

3.2 二維平面散點圖的構造

選擇采用數(shù)學形態(tài)學中皮爾遜相關系數(shù)來反映不同出線上零序電流低頻分量IMF1-1的相關性，其公式為

式中：ai、bi分別為a、b兩個信號的第i個采樣點的值；n為采樣點個數(shù)；ρ表示兩者的相關性，-1≤ρ≤1。ρ越趨近于-1，兩波形之間相關性越低；ρ越趨近1，兩波形之間的相關性越高。假設配電網(wǎng)共有m回出線，對各出線首個信號周期（20 ms）內的IMF1-1分量進行相關性計算，得到它們之間的相關系數(shù){ρij} ={ρ11,ρ12,…,ρmm} ，ρij表示出線Li與Lj低頻信號分量IMF1-1之間的相關系數(shù)。將每一個出線看作二維平面的一個散點，各散點橫坐標為該點所代表的出線IMF1-1分量與其他出線（包括其自身）的IMF1-1分量相關系數(shù)之和。例如，代表線路Li的點的橫坐標xi表達式為

各散點縱坐標構造方法如下。首先對各線路IMF2分量做一階差分計算，即

式中：diff為一階差分計算符號；IMFi-2表示出線Li的高頻分量；n為采樣點個數(shù)，當n=1時，式（15）為信號的一階差分初值。

然后使用符號函數(shù)計算其一階差分初值的極性，確定各線路高頻分量的極性，以表示其縱坐標。例如，代表線路Li的點的縱坐標yi表達式為

式中，sgn為符號函數(shù)，當信號初始極性為正時，sgn的值為1；當信號極性為負時，其值為-1。

故障特征較為明顯時，故障線路零序電流低頻IMF1-1分量波形與健全線路存在明顯差異，因此它們之間的相關系數(shù)也將較小，而健全線路之間的波形差異則相對較小，相關系數(shù)較大并趨近于1。根據(jù)式（14），做相加運算后，故障線路的散點橫坐標值會較小，而健全線路的散點橫坐標值則較大，且數(shù)值比較接近；又由于高頻IMF2分量初始極性存在差異，故障線路與健全線路散點的縱坐標值將為相反數(shù)。即使某一頻段下故障特征不明顯導致單一坐標無法區(qū)分故障線路，在二維平面中代表故障線路的散點依然是一個離群點。

3.3 K-means++聚類算法

基于劃分聚類的K-means++算法的原理簡單來說就是對一堆散點進行聚類，首先確定將這些散點聚成幾類，然后挑選幾個點作為初始中心點，并根據(jù)一定規(guī)則迭代重置聚類中心點，直到達到類中的點都足夠近、類間的點都足夠遠的目標效果。傳統(tǒng)的K-means 聚類算法采取隨機選取的方式確定初始聚類中心，但最終的聚類結果可能會受到不同初始聚類中心選取的影響，從而對選線結果產生影響。K-means++的改進之處在于，其在選取初始聚類中心時可以確保初始聚類中心之間的距離盡可能地遠，從而更高效、準確地進行聚類分析。

在本文中故障線路和非故障線路各自的散點分別形成一個簇（記為故障線路簇和健康線路簇），因此需要將所有散點劃分為2類。由于故障線路所屬的散點理論上屬于離群點，K-means++聚類算法在首輪分類時即可將代表故障線路的散點單獨劃分出來，從而大大提高了運算速度和選線準確率。

該算法通過對各出線零序電流信號IMF1-1分量相關系數(shù)和IMF2分量初始極性差異構造的點集P={P1,P2,…,Pm} 進行聚類時的具體步驟為如下。

步驟1從輸入的數(shù)據(jù)點集P中隨機選一個點作為故障線路簇Q1和健康線路簇Q2的首個初始聚類中心c1。

步驟2分別計算點集中其他的點到該聚類中心的歐氏距離，根據(jù)距離盡可能遠地原則選擇一個點作為另一個初始聚類中心c2。用xci、yci表示聚類中心ci的橫、縱坐標，第i個點到其的歐氏距離dis表示為

步驟3根據(jù)點集中所有點與這兩個聚類中心的距離最短原則將它們劃分到相應的線路簇。

步驟4重新計算兩個簇各自的聚類中心，即

式中：z表示單個樣本向量；|Qi|表示簇的大小。

步驟5將點集P中全部點根據(jù)與聚類中心歐氏距離最小原則進行劃分。

步驟6計算標準測度函數(shù)，當滿足收斂條件時，函數(shù)停止迭代，算法終止；如不滿足則重復步驟3～步驟5。

3.4 選線流程

選線流程如圖5所示，具體步驟如下。

圖5 故障選線流程Fig.5 Flaw chart of fault line selection

步驟1利用檢測裝置持續(xù)監(jiān)測系統(tǒng)母線零序電壓U0是否超過額定電壓UN的115%，符合條件即啟動故障選線裝置，提取各出線在故障發(fā)生時刻后2個信號周期（40 ms）內的零序電流信號波形。

步驟2將提取到的各出線上的零序電流信號用VMD 算法按照前文所述方法先過濾后再進行3層分解，得到預設個數(shù)的本征模態(tài)函數(shù)。

步驟3根據(jù)式（13）和（16），通過對故障發(fā)生時刻后首個信號周期（20 ms）內的低頻分量做相關性計算，對高頻分量測算初始極性，進而確定各點橫、縱坐標，作出二維平面散點分布圖。

步驟4使用K-means++算法進行聚類分析，篩選出故障線路。

4 仿真分析

4.1 仿真模型

本節(jié)將在仿真的基礎上分析驗證本文所提算法的實用性。通過仿真軟件Pscad 建立配電網(wǎng)模型，獲得各出線單相接地條件下的零序電流信號波形；將零序電流信號導入Matlab工作區(qū)中進行VMD算法分解，獲得故障信號不同頻段下的模態(tài)分量；使用本文所提出的方法進行模擬選線，得到最終結果數(shù)據(jù)。

在Pscad中搭建一個6條出線6個負荷的10 kV配電網(wǎng)模型，其結構如圖6所示。系統(tǒng)中，線路L1、L3、L4、L6為架空線路，除L1長度為15 km之外，其余長度均為12 km；線路L2為電纜線路，長度8 km；線路L5為纜線混合線路，電纜長4 km，架空線長12 km。各出線相關參數(shù)如表2所示，各負荷功率均設為2 MW，消弧線圈電感值采用過補償10%進行設置。

表2 線路參數(shù)Tab.2 Parameters of lines

圖6 仿真電路模型Fig.6 Simulation circuit model

4.2 不同初始條件下的仿真

仿真實驗選定的3 個變量分別是故障距離、故障初相角和接地過渡電阻。首先選取一個典型案例進行詳細分析，隨后在選定兩個變量而不斷對剩余的單一變量做出改變的情況下進行仿真實驗，以驗證本文所提出選線方法的適用性。

4.2.1 典型案例分析

假設距離出線L2始端4 km處發(fā)生A相單相接地短路故障，故障時刻為0.04 s，故障初相角設置為90°，接地過渡電阻大小設置為2 kΩ，通過故障錄波裝置提取出各出線上故障發(fā)生時刻后2 個信號周期（40 ms）的零序電流波形如圖7所示。

圖7 零序電流波形Fig.7 Zero-sequence current waveforms

將提取到的電流進行濾波處理后，按照前文所述方式進行3層VMD分解，提取故障條件下各出線零序電流分解結果中的IMF1-1分量和IMF2分量分別如圖8和圖9所示。

圖8 各出線IMF1-1波形Fig.8 IMF1-1 waveform of each line

圖9 各線路IMF2波形Fig.9 IMF2 waveform of each line

觀察圖8和圖9中各出線上零序電流分解后得到的低頻分量IMF1-1和高頻分量IMF2可知，在故障初相角為90°的情況下，故障線路與非故障線路低頻分量之間的波形差異不大，說明故障線路L2上流過的零序電流低頻分量中幾乎不含衰減直流分量；除故障線路L2的高頻分量初始極性為正外，其余各線路中零序電流高頻分量的初始極性全部為負。

算例的分解結果與本文第1 節(jié)中關于零序電流的成分及特性分析相吻合。

計算各出線IMF1-1分量之間的相關系數(shù)，所得結果如表3所示。

表3 各線路IMF1-1之間的相關系數(shù)Tab.3 Correlation coefficient between IMF1-1 of each line

計算各出線IMF2分量的初始極性，所得結果如表4所示。

表4 各線路IMF2的初始極性Tab.4 Initial polarity of IMF2 of each line

計算代表各出線的散點的橫、縱坐標，并做出所有散點在二維平面分布，結果如圖10所示。

圖10 二維平面的散點分布Fig.10 Scatter distribution on two-dimensional plane

從圖中可以看到各散點的橫坐標值均較為接近，代表僅通過低頻分量之間的波形相關性并不足以使故障線路的散點成為離群點，必須通過縱坐標代表的高頻分量初始極性之間的差異構造離群點。由于初相角φ=90°，高頻分量初始極性差異較為明顯，保證了故障線路與健全線路散點的縱坐標值為相反數(shù)，因此故障線路L2的點在平面中仍然成為了一個離群點。

使用K-means++聚類算法對圖10 中的點集進行聚類分析，得到輸出聚類標簽為[2 1 2 2 2 2]，表示出線L2與其他線路不同，單獨劃分為了一類，因此可以得出結論為出線L2發(fā)生了單相接地故障。

4.2.2 適用性分析

選取故障初相角、故障距離以及過渡電阻為3個變量，在選定兩個變量而不斷對剩余的單一變量做出改變的情況下進行仿真實驗。受限于篇幅不再詳細分析，仿真結果如表5所示。

表5 仿真結果Tab.5 Results of simulation

從表5 中可以看出，不論故障初始條件如何改變，包括過渡電阻超過2 kΩ 等對選線不利的工況下，本文提出的算法依然能夠準確地選擇出故障線路。除表5中列舉的情況，還進行了大量的仿真實驗，包括改變噪聲強度等，本文提出的算法均能正確地選出故障線路，表明該選線方法不受條件改變的影響，具有較強的適用性。

4.2.3 方法優(yōu)越性比較

為驗證本文所提出的方法的優(yōu)越性，現(xiàn)假設在距離出線L5始端4 km 處發(fā)生單相接地故障，故障初相角設為60°，過渡電阻依次設置為100 Ω、500 Ω、1 000 Ω。分別用本文的方法和進行EMD分解再使用K-means++聚類（方法2）以及進行VMD分解再使用K-means聚類（方法3）進行比較，結果如表6所示。

表6 不同方案的選線結果Tab.6 Line selection results under different schemes

通過對比3 種方法的結果可以看出，隨著過渡電阻的增大，高頻分量幅值的減小，本文的算法均正確地選擇出了故障線路，而方法2 有錯誤記錄。方法3 雖然也能正確選線，但是相較本文算法，Kmeans 算法在執(zhí)行聚類任務時運行更加耗時（約為40 ms），約為K-means++算法（約為10 ms）的4 倍，因此本文設計的選線方法優(yōu)于基于EMD分解或Kmeans聚類的相關性分析方法。

5 結 論

通過對諧振接地系統(tǒng)單相接地短路故障時故障零序電流組成特性的分析，本文提出了一種基于VMD 分解和改進K-means++相關性聚類的故障選線方法，實現(xiàn)了精確故障選線的目的。主要結論如下。

（1）本文將VMD分解和K-means++聚類算法相結合，利用在構造的散點圖中，故障線路信號所屬的散點屬于離群點的特點，使用由K-means 聚類算法改進而來的K-means++聚類算法進行分析，相比于K-means算法，提高了選線速率。

（2）本文將故障零序電流信號不同模態(tài)分量之間的差異投射到二維平面中，從兩個維度進行聚類分析，極大提高了選線準確率和保護裕度。

（3）加入噪聲進行仿真實驗，結果表明本文方法對噪聲具有一定的魯棒性。