實測鐵磁諧振時間序列的非線性動力學分析

黃艷玲　司馬文霞　楊　鳴　楊　慶　袁　濤

(輸配電裝備及系統安全與新技術國家重點實驗室(重慶大學)　重慶　400044)

?

實測鐵磁諧振時間序列的非線性動力學分析

黃艷玲司馬文霞楊鳴楊慶袁濤

(輸配電裝備及系統安全與新技術國家重點實驗室(重慶大學)重慶400044)

摘要在僅有實測單一電壓時間序列的情況下，為準確識別出電力系統鐵磁諧振類型，運用坐標延遲方法重構與原系統拓撲意義上等價的相空間。采用相平面、龐加萊截面和關聯維數3種基于重構相空間的非線性動力學分析工具表達電壓時間序列的運動特征，由此識別其所屬的鐵磁諧振類型。對3例典型的中性點不接地系統電磁式電壓互感器鐵磁諧振實測時間序列進行了動力學特征分析。其中兩例時間序列的相平面軌跡、龐加萊截面均呈現周期運動特性，關聯維數估計值分別為1.015 0±0.003 9、1.006 1±0.000 6，由此判別它們作周期運動；另一例時間序列的相平面軌跡、龐加萊截面均未呈現周期或準周期運動特性，關聯維數估計值為2.300 2±0.061 2，由此判別其作混沌運動。實測電壓時間序列受到多種因素的影響，綜合以上3種分析方法所表達的特征，對所發生的鐵磁諧振類型作出的判斷更準確和更有說服力。

關鍵詞：鐵磁諧振時間序列非線動力學相空間重構關聯維數電磁式電壓互感器

Nonlinear Dynamic Analysis of the Measured Ferroresonance Time Series

HuangYanlingSimaWenxiaYangMingYangQingYuanTao

(State Key Laboratory of Power Transmission Equipment and System Safety and New Technology Chongqing UniversityChongqing400044China)

AbstractWith only a single voltage time series recorded in fact，for accurately identifying the ferroresonance types，the delay coordinates method was applied to obtain the reconstructed phase space that is topologically equivalent to the phase space of the original system.Nonlinear dynamics methods based on the reconstructed phase space，e.g.phase plane，Poincaré section and correlation dimension，were used to express the dynamic characteristics of the time series and identify the type of ferroresonance based on them.The nonlinear dynamic characteristics analysis was conducted to 3 typical cases of measured time series of ferroresonance occurred in inductive voltage transformer in isolated neutral system.The phase plane trajectories and Poincaré sections of two time series cases of them show periodic motion characteristics，their correlation dimension estimation values respectively are 1.015 0±0.003 9，1.006 1±0.000 6；therefore their motion modes were identified as fundamental mode and subharmonic mode.The phase plane trajectory and Poincaré section of another time series case do not show periodic or quasiperiodic motion characteristics，its correlation dimension estimation value is 2.300 2±0.061 2；therefore its motion mode was identified as chaotic mode.The measured voltage time series is affected by many factors，the identification result of ferroresonance mode occurred is more accurate and convictive through synthesizing the features characterized by the above 3 analysis methods.

Keywords：Ferroresonance，time series，nonlinear dynamics，phase space construction，correlation dimension，inductive voltage transformer

0引言

鐵磁諧振是一種復雜的電現象，它可能發生在由電力變壓器或電磁式電壓互感器(PT)及帶鐵心的電抗器等非線性電感元件與電容元件構成的回路中[1]。非線性電感元件的飽和勵磁特性使鐵磁諧振行為呈現非線性動力學特性，其穩態響應模式高度依賴于系統參數和初始條件，初始條件稍有變化就可能使響應模式發生改變[2]。穩態響應模式可能是周期、準周期或混沌模式。在合適的條件組合下，很多設備操作行為及單相接地故障消除瞬間[3]均可能激發鐵磁諧振，由此產生的過電壓或過電流嚴重危害著變電站設備的絕緣，甚至造成變壓器或PT的燒毀，繼而引發大面積的停電事故。因此，當鐵磁諧振發生時，很有必要識別出諧振類型，為采取有效的過電壓抑制措施和電網事故分析提供指導。這對保障電網的安全運行具有重要意義。

常規的線性數學方法不適于分析非線性的鐵磁諧振行為，Z.Emin和其他學者將基于相空間的非線性動力學分析工具相平面[1，2，4，5，10]、龐加萊截面[1，2，4-10]和Lyapunov指數[4，8，10]等運用到鐵磁諧振行為的特征分析和識別中，再結合功率譜密度分析[1，6，9]對識別結果進行驗證。在此基礎上，為進一步獲取系統發生諧振的參數取值范圍以及諧振類型，基于分岔理論的分岔圖[2，5，8，9]被用于預測系統隨某一參量變化可能出現的諧振類型。

盡管基于相空間的非線性動力學分析工具對鐵磁諧振行為的特征分析和識別取得了很好的效果，但當前所研究的鐵磁諧振數據主要是通過建立諧振回路的微分方程組進行數值求解[1，4，6-10]或是應用ATP-EMTP軟件搭建回路仿真[5]所得，與實際情況有一定的差距。因建立的數學仿真模型是系統結構的簡化表達，鐵心損耗及系統對地電容等系統參數難以準確獲取，尤其對電感元件的非線性特性模擬的精度有待提高。只有對正在發生的鐵磁諧振類型進行識別，才能為制定合適的過電壓抑制策略提供幫助，但實際觀測和獲取的常常只是鐵磁諧振系統單一狀態變量(如電壓的時間序列)，而以上基于相空間的研究方法[1，2，4-10]至少需要兩個狀態變量的時間序列數據，因此需要借助其他方法先構造相空間。

因一個變量隨時間的變化隱含了整個系統的運動規律，根據F.Takens提出的嵌入定理[11]，從一維時間序列可重構一個與原動力系統在拓撲意義上等價的相空間。因此，可對單一變量的時間序列運用相空間重構方法構造相空間，通過分析重構相空間軌線的分布或結構得到原系統的動力學特征。

本文在僅有某一中性點不接地系統PT鐵磁諧振實測電壓時間序列的情況下，運用相空間重構方法獲得系統相空間。在此基礎上，運用相平面、龐加萊截面和關聯維數3種基于相空間重構的非線性動力學分析工具表達系統的動力學特征，以識別鐵磁諧振類型。

1非線性動力學分析方法

1.1相空間重構

系統相空間是由狀態變量xi(i=1,2,…,n)張成的空間Rn， 相空間的一個點對應系統的一個狀態。N.H.Parkard等[12]于1980年提出了基于單一變量時間序列的兩種相空間重構方法：導數重構法和坐標延遲重構法。因實測數據混有噪聲，對其求取的微分誤差較大，因此實際應用中主要使用坐標延遲重構法進行一維時間序列的相空間重構。

1.1.1坐標延遲重構法

設有一單變量的時間序列u={ui|i=1,2,…,N}， 對于某一給定的嵌入維數m和延遲時間τ，其相空間向量序列為

V={Vk|uk,uk+τ,…,uk+(m-1)τ}

(1)

式中，k=1,2,…,n；n=N-(m-1)τ。

由此，重構相空間的軌線分布或結構便反映了系統的運動特征。將重構相空間的某兩個向量的取值序列分別作為橫、縱坐標繪制在平面中就構成重構相平面。系統作周期運動時，其相平面軌跡是一條不斷重復的閉合曲線[2]；而作混沌運動時，相平面軌跡永不重復[1]。

1.1.2坐標延遲重構參數選取

坐標延遲重構法參數即延遲時間τ和嵌入維數m的選取很關鍵。當前對它們的選取有兩種觀點[13]，一種觀點認為二者互不相關，對τ和m的選取可獨立進行，分別有不同的確定方法。最佳τ值的求取方法主要有自相關函數法[14]和互信息法[15]。最佳m值的估計方法有計算吸引子不變量法(或G-P法)[16]、奇異值分解法[17]、偽最近鄰法[18]和偽最近鄰Cao氏改進法[19]。另一種觀點認為τ和m相互關聯，即τ和m的關系與重構相空間的嵌入窗τw(τw=(m-1)τ)密切相關。D.Kugiurmtzis[20]于1996年首次提出嵌入窗思想，基于嵌入窗思想的典型算法有C-C法[21]和自動嵌入式相空間重構法[22]。

本文在確定相空間重構參數時采用第一種觀點，即τ和m的值采用不同的方法確定。采用互信息法獲取最佳τ值，因互信息刻畫的是兩個變量間的普遍(或總體)相關性，適用于非線性系統。本文采用的計算過程無需事先確定最優m值，但通過計算一組m值下的關聯維數可間接獲取最優m值(即G-P法確定最優m值，實現方法見1.3.1節)。

1.2龐加萊截面

龐加萊截面是相平面軌跡的一種簡化表達，是對相平面軌跡以一定的頻率采樣獲取的相點圖。

若對相平面軌跡采樣的頻率為系統外加激勵函數的頻率(如電源頻率)，則對作周期運動的系統，在龐加萊截面上表現為一個不動點或是多個孤立點[1]；對作準周期運動的系統，龐加萊截面是一條分布著無限個不同點的封閉曲線[6，7]；而當系統作混沌運動時，龐加萊截面呈現為局限于平面某一特定區域的一些隨機分布的相點集合[1]。

1.3關聯維數

判斷一個時間序列是否具有混沌特性最為常用的動態不變量有Lyapunov指數[4，8，23]和吸引子維數[7，23]。由于實驗獲取的時間序列數據ui受到多種因素(如采樣精度、采樣長度、噪聲和誤差等)的影響，僅憑其最大Lyapunov指數λmax的正負就對系統的運動特性作出判斷，很可能不準確[24]。自1983年P.Grassberger等[25]提出針對一維時間序列數據的關聯維數計算方法(簡稱G-P算法)后，關聯維數作為吸引子維數的一種度量得到了廣泛應用。本文采用關聯維數作為判別系統混沌特性的動態不變量。

1.3.1關聯維數G-P算法

給定一尺度r和嵌入維數m，定義相空間距離小于r的點對數與所有點對數的比例為關聯積分Cm(r)，即

(2)

式中，n為相空間的點數；‖·‖表示求兩個重構向量的歐氏范數；θ(x)為Heaviside階躍函數。

(3)

對于無限長序列(n→∞)和足夠小的r值，若Cm(r)與r呈冪律關系，即Cm(r)∝rD2， 則定義關聯維數D2為

(4)

實際所獲取的數據量有限，無法滿足極限條件n→∞， 而且尺度r的取值須滿足式(5)。

min{‖Vi-Vj‖}

(5)

式中，i≠j，i,j=1,2,…,n。

因此，計算關聯維數采取的常規步驟為：

1)在式(5)的范圍內取一尺度序列ri(i=1,2,…,l)，計算與ri對應的關聯積分序列Cm(ri)。

2)求出序列Cm(ri)和ri的對數序列lnCm(r)i和lnri，以序列lnri為橫坐標、序列lnCm(r)i為縱坐標繪制雙對數曲線。

3)將雙對數曲線線性區(或無標度區)的斜率作為D2的估計值。

因雙對數曲線或多或少呈現S形，很難通過目測或經驗確定其線性區域。文獻[26]指出一個可信的無標度區必須是清晰可辨的，表現在lnCm(r)對lnr曲線的局部斜率dsm對lnr曲線存在隨m增大而趨于飽和的平坦區域。

為使估算所得D2盡量準確，本文計算關聯維數采取的步驟與常規的不同之處為：

1)計算一組m(m=2~10)的對數序列lnCm(r)i。

2)分別求出m=2～10下lnCm(r)對lnr曲線的局部斜率dsm。

3)若dsm對lnr曲線存在隨m值增大而趨于飽和的平坦區域，則將各m值下dsm對lnr曲線的平坦區域作為無標度區，并將此區域對應的數據點對(lnCm(r)，lnr)進行直線擬合得到序列D2_m，由此獲取D2估計值。同時，D2_m趨于飽和時對應的m值可作為該時間序列重構的最優m值，即G-P法確定最優m值。

局部斜率dsm的計算式為

(6)

式中，i=1,2,…,l-1。

對于隨機過程的時間序列，其D2隨m的增大而增大，不會趨于飽和；而確定性系統的D2隨m的增大而逐漸趨于飽和。不同的D2飽和值對應于不同的系統運動狀態[27]：D2=1，系統作周期運動；D2=2，系統作準周期或擬周期運動；D2>2或D2不為整數，系統作混沌運動。

1.3.2關聯維數G-P優化算法

對于長度為n的數據集，在所有的點對中直接搜尋距離小于r的點對所需的計算量為O(n2/2)。如果需要在多個嵌入維下估算D2，則計算耗時還要增加。由于在估算D2時，只有相空間的鄰近點對對Cm(r)值有貢獻。因此，如果事先大致確定每個相點的近鄰分布位置，可大大縮短搜索鄰近點對的耗時。基于以上思想，文獻[28-30]先后提出采用盒子網格的方法幫助尋找距離小于r的鄰近點對。文獻[30]提出的基于直方圖的盒子輔助(box-assisted)近鄰搜索算法，解決了將相空間的點放入盒子而不浪費存儲空間的技術難題，而且易于實現。本文應用這一算法計算Cm(r)。

2典型鐵磁諧振時間序列的非線性動力學分析

圖1　重慶市先鋒變電站一次主接線簡圖Fig.1　Primary wiring simplified diagram of Xianfeng substation in Chongqing city

本文所研究的鐵磁諧振時間序列是重慶市先鋒變電站的實測數據。該變電站的一次主接線簡圖如圖1所示，圖中‘/’表示隔離開關或接地刀開關。該站有兩臺3繞組變壓器，繞組接線方式為YN/yn0/d11，電壓等級為110 kV/35 kV/10 kV，35 kV側中性點不接地。為監視母線絕緣，各電壓等級的Ⅰ、Ⅱ段母線均安裝有PT。為實現對該站3個電壓等級的過電壓監測，分別在10 kV開關柜母線，35 kV戶外開關場母線，2#主變壓器套管110 kV側安裝了電壓傳感器[31]。自2006年11月16日開始實施對10 kV系統Ⅰ、Ⅱ段母線和2#主變壓器110 kV系統的過電壓監測，自2008年6月23日開始實施對35 kV系統Ⅰ、Ⅱ段母線的過電壓監測。10 kV Ⅱ段母線和35 kV Ⅰ段母線上均監測到鐵磁諧振過電壓，本文選取具有代表性的鐵磁諧振電壓序列數據進行特征分析，為進一步研究其分類識別方法奠定基礎。

本文對原始采集數據進行了下采樣處理和濾波去噪。①下采樣。原始數據的采樣頻率為200 kHz，對于鐵磁諧振的穩態過程而言過高，因其所含成分的最高頻率不大于1 kHz。為使被分析的數據展現電壓變量更多的變化規律，又不使計算量過大，對原始數據采取整數倍抽取的方式獲得待分析的時間序列。抽取倍數的確定原則為：最低抽取頻率設為10 kHz，抽取頻率必須為50 Hz的整數倍，抽取后的相平面軌跡與未抽取前的一致。因此，本文采用實驗法確定合適的抽取倍數：以抽取倍數的取值范圍2、4、5、8、10、16、20依次進行抽取，并對比抽取前、后的相平面軌跡，若發現有明顯差異，則將前一次采用的抽取倍數作為最終的信號抽取倍數。②濾波去噪。因實測數據不可避免地含有噪音，對下采樣處理所得的時間序列采用sym8小波5層啟發式閾值估計法去噪。

1)時間序列1

時間序列1是2007年2月14日02∶05∶41于10 kV Ⅱ段母線處采集的A相原始數據進行5倍抽取和濾波處理所得的數據，數據長度為10 000。

圖2和圖3分別為時間序列1的波形圖和重構相平面圖(N=10 000，τ=26)。圖2呈現明顯的周期性，過電壓幅值達2.802 6 (pu)。圖3所示的相平面軌跡是一條近似橢圓形的封閉曲線，表明時間序列1作周期運動。圖4為對應于圖3的龐加萊截面，采樣頻率為電源頻率。本文中如不作特別說明，對相平面軌跡采樣的頻率均與電源頻率相同。圖4中的相點幾乎重合，可近似為一個孤立點，表明時間序列1作基頻周期運動，頻率與電源頻率相同。

圖2　時間序列1的波形圖Fig.2　Waveform plot of time series 1

圖3　時間序列1的重構相平面圖Fig.3　Reconstructed phase plane plot of time series 1

圖4　時間序列1的龐加萊截面Fig.4　Poincaré section of time series 1

圖5為時間序列1的lnCm(r)對lnr雙對數曲線圖，圖6為對應于圖5的dsm對lnr曲線圖。圖5和圖6中顏色由淺到深的線條分別表示m=2、3、…、10下的相應曲線，后文類似圖形中各線條的含義與此相同。尺度r的取值序列為：ri=0.06×1.311 6i-1，i=1,2,…,16。 由圖5可知，各條lnCm(r)對lnr曲線的下半部分線性度較好，且在lnr取值區間[-2.270 9，-1.185 9]，m=3～10的lnCm(r)對lnr曲線基本平行。在圖6中可觀察到各條dsm對lnr曲線均存在清晰的平坦區域，且在lnr取值區間[-2.270 9，-1.475 2]，m=3～10的dsm對lnr曲線基本重合。

圖5　時間序列1的lnCm(r)對lnr曲線Fig.5　Curves lnCm(r) versus lnr of time series 1

圖6　時間序列1的dsm對lnr曲線Fig.6　Curves dsm versus lnr of time series 1

m=2～10時，D2_m分別為1.024 3、1.013 3、1.017 5、1.018 9、1.018 8、1.017 0、1.017 4、1.013 6、1.012 2。m≥3時，D2_m趨于穩定，由此估算所得D2=1.015 0±0.003 9。該序列重構所需的最小嵌入維數為3。

時間序列1的D2相對于整數1的誤差低于2%，非常小，因此其D2理論值應為整數1，作周期運動。這與從相平面軌跡和龐加萊截面得出的結論一致。

對與時間序列1同一時間同一位置獲取的B相、C相數據也進行了動力學特征分析，其運動模式與時間序列1的相同，D2分別為1.007 6±0.008 7，1.028 8±0.011 5。

根據識別結果可知，2007年2月14日 02∶05∶41于10 kV Ⅱ段母線監測的鐵磁諧振類型是基頻周期諧振，其過電壓幅值為2.802 6(pu)，一般是不危險的，但若持續時間長，可能導致PT絕緣破壞，甚至燒毀或避雷器爆炸。

2)時間序列2

時間序列2是2009年7月16 日21∶52∶26于35 kV Ⅰ段母線處采集的A相原始數據進行10倍抽取和濾波處理所得的數據，數據長度為10 000。

圖7和圖8分別為時間序列2的波形圖和重構相平面圖(N=10 000，τ=44)。從圖7可知，t2內的電壓取值與t1內的基本相同，t2之后的電壓值也有重復的跡象，可初步斷定時間序列2具有周期性。t1和t2內有4 400個采樣點，t1=t2=220 ms(采樣頻率為20 kHz)。過電壓幅值為-1.876 2 (pu)。圖8所示的相平面軌跡是一條由多條軌道組成的封閉曲線，不同的軌道之間有交點。這些特征表明時間序列2作周期運動，其周期比較長。

圖7　時間序列2的波形圖Fig.7　Waveform plot of time series 2

圖8　時間序列2的重構相平面圖Fig.8　Reconstructed phase plane plot of time series 2

圖9　時間序列2的龐加萊截面Fig.9　Poincaré section of time series 2

與圖8對應的龐加萊截面如圖9所示。經仔細分析發現，最早出現在龐加萊截面上的11個相點呈離散分布，但分布在龐加萊截面次對角線兩端的相點相距較近。此時僅憑圖9難以判斷時間序列2的運動性質。此后獲取的相點依次分布在原11個相點的附近，甚至與原有相點重合。如果時間序列長度為4 400(11個工頻周期)的整數倍，則這11個相點分布區域包含的相點數相同。結合相平面軌跡的分析結論，可斷定這11個相點區域分別表示11個孤立的離散點(用不同的標記表示)。這表明時間序列2作分頻周期運動，其頻率為工頻的1/11。受測量誤差的影響，11個孤立離散區域處相點的分散性不同。

對時間序列1和2的分析表明，測量誤差的存在導致龐加萊截面上應重疊的相點未完全重疊，給分析帶來困難。因此，對周期比較長的時間序列，運用龐加萊截面時需非常謹慎，時間序列數據要足夠長，至少要包含2個運動周期。

圖10為時間序列2的lnCm(r)對lnr雙對數曲線圖。尺度r的取值序列為：ri=0.06×1.311 6i-5，i=1,2,…,16。 從圖10可知，m=2～5時lnCm(r)對lnr曲線的形狀差異較大，m≥6時lnCm(r)對lnr曲線的形狀逐漸一致，且m=9與m=10的基本重合。與圖10相對應的dsm對lnr曲線圖如圖11所示。從圖11可觀察到，m≥6時dsm對lnr曲線均存在明顯的平坦區域，且在lnr取值區間[-3.084 7，-2.270 9]，dsm對lnr曲線平坦且基本重合。

圖10　時間序列2的lnCm(r)對lnr曲線Fig.10　Curves lnCm(r) versus lnr of time series 2

圖11　時間序列2的dsm對lnr曲線Fig.11　Curves dsm versus lnr of time series 2

m=6~10時，D2_m分別為1.006 3、1.006 7、1.006 7、1.005 7、1.005 6，波動極小，由此估算可得D2=1.006 1±0.000 6。該序列重構所需的最小嵌入維數為6。

時間序列2的D2相對于整數1的誤差低于0.7%，非常小，因此其D2理論值應為整數1，作周期運動。這與從相平面軌跡和龐加萊截面得出的結論一致。

對與時間序列2同一時間同一位置獲取的B相、C相數據也進行了動力學特征分析，其運動模式與時間序列2的相同，D2分別為1.001 3±0.001 0，1.018 4±0.005 3。

根據識別結果可知，2009年7月16日21∶52∶26于35 kV Ⅰ段母線監測的鐵磁諧振類型是分頻周期諧振，過電壓幅值為-1.876 2 (pu)，電壓本身是不危險的，但分頻引起的PT過飽和電流，若長時間持續，易導致PT高壓熔絲熔斷，或引起PT因嚴重過熱而燒毀，甚至爆炸。

3)時間序列3

時間序列3是2009年4月18日01∶29∶58于35 kV Ⅰ段母線處采集的A相原始數據進行10倍抽取和濾波處理所得的數據，數據長度為10 000。

圖12和圖13分別為時間序列3的波形圖和重構相平面圖(N=10 000，τ=57)。時間序列3的各波峰之間呈現此消彼長的現象，有一定的重復性，但又不完全重復。過電壓幅值為-2.127 5 (pu)。圖13呈現的相平面軌跡不是一條封閉曲線，但相點的運動軌道被限制在一個固定的區域，即其運動軌跡并不是無限制擴散和完全雜亂無章。這些特征表明時間序列3不是作周期運動，也不是作隨機運動。

圖12　時間序列3的波形圖Fig.12　Waveform plot of time series 3

圖13　時間序列3的重構相平面圖Fig.13　Reconstructed phase plane plot of time series 3

與圖13對應的龐加萊截面如圖14所示。圖14的相點呈不均勻的離散分布，未形成一條近似閉曲線，這表明時間序列3不是作準周期運動。因圖14中的相點數較少，未呈現自相似的分層結構，因此不能判斷時間序列3是否作混沌運動。

圖14　時間序列3的龐加萊截面Fig.14　Poincaré section of time series 3

時間序列3的lnCm(r)對lnr雙對數曲線圖如圖15所示，對應于圖15的dsm對lnr曲線圖如圖16所示。尺度r的取值序列為：ri=0.06×1.3116i-1，i=1,2,…,16。 從圖15可知，在lnr取值區間[-0.643 4，0.170 3]，雙對數曲線隨m的增大而逐漸平行，即斜率值逐漸達到飽和。相應地，在圖16中可清晰的觀察到dsm對lnr曲線存在隨m的增大而逐漸飽和的平坦區域，且m=8～10時dsm對lnr曲線的平坦區域基本重合。

圖15　時間序列3的lnCm(r)對lnr曲線Fig.15　Curves lnCm(r) versus lnr of time series 3

圖16　時間序列3的dsm對lnr曲線Fig.16　Curves dsm versus lnr of time series 3

m=2～10時，D2_m分別為1.678 4、1.837 0、2.118 5、2.341 2、2.361 4、2.264 4、2.238 9、2.302 4、2.312 1。m≥5時，D2_m趨于穩定，由此估算可得D2=2.300 2±0.061 2。該序列重構所需的最小嵌入維數為5。

時間序列3的D2>2，不為整數，且小數部分較大(>0.230 0)，表明其具有混沌運動特性。從龐加萊截面得出的其不是作準周期運動的判斷，可驗證這一結論。

對與時間序列3同一時間同一位置獲取的B相、C相數據也進行了動力學特征分析，其運動模式與時間序列3的相同，D2分別為2.223 6±0.043 6，2.408 3±0.082 4。

根據識別結果可知，2009年4月18日01∶29∶58于35 kV Ⅰ段母線監測的鐵磁諧振類型是混沌諧振，其過電壓幅值僅稍微超過2 (pu)，電壓本身也是不危險的。對時間序列3的頻譜分析得知，其含24 Hz頻率分量的幅值達1.161 4 (pu)，即含較大幅值的分次諧波。因此，若長時間持續，有與分頻諧振類似的危害。

3結論

1)為準確識別實測鐵磁諧振一維電壓時間序列所屬的諧振類型，本文運用相空間重構方法構造系統相空間。在此基礎上，采用相平面、龐加萊截面和關聯維數3種非線性動力學分析方法，以不同的形式表達其動力學特征，由此對其所屬的鐵磁諧振類型作出了準確判斷。

2)對3例典型的中性點不接地系統PT鐵磁諧振實測時間序列進行了動力學特征分析。其中兩例時間序列的相平面軌跡、龐加萊截面均呈現周期運動特性，關聯維數估計值分別為1.015 0±0.003 9、1.006 1±0.000 6，由此判斷它們作周期運動。另一例時間序列的相平面軌跡、龐加萊截面均未呈現周期或準周期運動特性，關聯維數估計值為2.300 2±0.061 2，由此判斷其作混沌運動。

3)鑒于實測電壓時間序列數據受到多種因素(如采樣精度、噪音和誤差等)的影響，龐加萊截面和關聯維數估計值與理想情況有一定差距，要對其運動特性作出準確判斷，有必要綜合多種方法所表達的特征。

4)中性點不接地系統PT鐵磁諧振屬于零序性質，發生鐵磁諧振時對任一相電壓數據進行非線性動力學分析，均能對其類型作出準確判斷。

5)本文的研究為鐵磁諧振類型的在線自動識別和諧振過電壓的抑制策略研究奠定了一定的基礎，為實測鐵磁諧振數據的識別研究提供了新思路，豐富了鐵磁諧振的研究實例。

參考文獻

[1]Emin Z，Al Zahawi B A T，Auckland D W，et al.Ferroresonance in electromagnetic voltage transformers：a study based on nonlinear dynamics[J].IEE Proceedings-Generation，Transmission and Distribution，1997，144(4)：383-387.

[2]Mork B A，Stuehm D L.Application of nonlinear dynamics and chaos to ferroresonance in distribution systems[J].IEEE Transactions on Power Delivery，1994，9(2)：1009-1017.

[3]周麗霞，尹忠東，鄭立.配網PT鐵磁諧振機理與抑制的試驗研究[J].電工技術學報，2007，22(5)：153-158.

Zhou Lixia，Yin Zhongdong，Zheng Li.Research on principle of PT resonance in distribution power system and its suppression[J].Transactions of China Electrotechnical Society，2007，22(5)：153-158.

[4]張博，魯鐵成，杜曉磊.中性點接地系統鐵磁諧振非線性動力學分析[J].高電壓技術，2007，33(1)：31-35.

Zhang Bo，Lu Tiecheng，Du Xiaolei.Nonlinear dynamic analysis of ferroresonance in neutral-grounded power system[J].High Voltage Engineering，2007，33(1)：31-35.

[5]Val E M，Dudurych I，Redfern M A.Characterization of ferroresonant modes in HV substation with CB grading capacitors[J].Electric Power Systems Research，2007，77(11)：1506-1513.

[6]Araujo A E A，Soudack A C，Marti J R.Ferroresonance in power systems：chaotic behaviour[J].Generation，Transmission and Distribution，IEE Proceedings C，1993，140(3)：237-240.

[7]Emin Z，Al Zahawi B A T，Yu K T，et al.Quantification of the chaotic behavior of ferroresonant voltage transformer circuits[J].Circuits and Systems I：Fundamental IEEE Transactions on Theory and Applications，2001，48(6)：757-760.

[8]Mozaffari S，Henschel S，Soudack A C.Chaotic ferroresonance in power transformers[J].IEE Proceedings-Generation，Transmission and Distribution，1995，142(3)：247-250.

[9]Kieny C.Application of the bifurcation theory in studying and understanding the global behavior of a ferroresonant electric power circuit[J].IEEE Transactions on Power Delivery，1991，6(2)：866-872.

[10]司馬文霞，鄭哲人，楊慶，等.用參數不匹配混沌系統的脈沖同步方法抑制鐵磁諧振過電壓[J].電工技術學報，2012，27(6)：218-225.

Sima Wenxia，Zheng Zheren，Yang Qing，et al.Suppressing the ferroresonance overvoltage by impulsive synchronization of chaotic systems with parameter mismatched[J].Transactions of China Electrotechnical Society，2012，27(6)：218-225.

[11]Takens F.Detecting strange attractors in turbulence[M]//Lecture Notes in Mathematica.Berlin：Springer Press，1981：366-381.

[12]Packard N H，Crutchfield J P，Farmer J D，et al.Geometry from a time series[J].Physical Review Letters，1980，45：712-715.

[13]行鴻彥，龔平，徐偉.嵌入窗方法確定混沌系統重構參數的仿真研究[J].系統仿真學報，2013，25(6)：1219-1225.

Xing Hongyan，Gong Ping，Xu Wei.Simulation research of chaos system reconstruction parameters based on embedded window[J].Journal of System Simulation，2013，25(6)：1219-1225.

[14]Albano A M，Muench J，Schwartz C，et al.Singular-value decomposition and the Grassberger-Procaccia algorithm[J].Physical Review A，1988，38(6)：3017-3026.

[15]Fraser A M，Swinney H L.Independent coordinates for strange attractors from mutual information[J].Physical Review A，1986，33(2)：1134-1440.

[16]Grassberger P，Procaccia I.Measuring the strangeness of strange attractors[J].Physica D Nonlinear Phenomena，1983，9(1/2)：189-208.

[17]Broomhead D S，King G P.Extracting qualitative dynamics from experimental data[J].Physica D，1986，20(2/3)：217-236.

[18]Kennel M B，Brown R，Abarbanel H D I.Determining embedding dimension for phase-space reconstruction using a geometrical construction[J].Physical Review A，1992，45(6)：3403-3411.

[19]Cao L Y.Practical method for determining the minimum embedding dimension of a scalar time series[J].Physica D，1997，110(1/2)：43-50.

[20]Kugiumtzis D.State space reconstruction parameters in the analysis of chaotic time series—the role of the time window length[J].Physica D：Nonlinear Phenomena，1996，95(1)：13-28.

[21]Kim H S，Eykholt R，Salas J D.Nonlinear dynamics，delay times，and embedding windows[J].Physica D：Nonlinear Phenomena，1999，127(1/2)：48-60.

[22]Otani M，Jones A.Automated embedding and creep phenomenon in chaotic time series[EB/OL].http：//users.cs.cf.ac.uk/Antonia.J.Jones/UnpublishedPapers/Creep.pdf/，2000-10-14/2011-03-15.

[23]Franca L，Savi M A.Distinguishing periodic and chaotic time series obtained from an experimental nonlinear pendulum[J].Nonlinear Dynamics，2001，26(3)：253-271.

[24]劉秉正，彭建華.非線性動力學[M].北京：高等教育出版社，2005.

[25]Grassberge P，Procaccia I.Characterization of strange attractors[J].Physical Review Letters，1983，50(5)：346-349.

[26]Kantz H，Schreiber T.Nonlinear time series analysis[M].Cambridge：Cambridge University Press，1997.

[27]張雨.時間序列的混沌和符號分析及實踐[M].長沙：國防科技大學出版社，2007.

[28]Theiler J.Efficient algorithm for estimating the correlation dimension from a set of discrete points[J].Physical Review A，1987，36(9)：4456-4462.

[29]Grassberger P.An Optimized box-assisted algorithm for fractal dimensions[J].Physics Letters A，1990，148(1/2)：63-68.

[30]Schreiber T.Efficient neighbor searching in nonlinear time-series analysis[J].International Journal of Bifurcation and Chaos，1995，5(2)：349-358.

[31]杜林，李欣，司馬文霞，等.110 kV變電站過電壓在線監測系統及其波形分析[J].高電壓技術，2012，38(3)：535-543.

Du Lin，Li Xin，Sima Wenxia，et al.Overvoltage on-line monitoring system for 110 kV substation and its waveforms analysis[J].High Voltage Engineering，2012，38(3)：535-543.

黃艷玲女，1978年生，博士，講師，研究方向為電力系統過電壓。

E-mail：huangyanling@cqu.edu.cn(通信作者)

司馬文霞女，1965年生，教授，博士生導師，研究方向為電力系統的防雷與過電壓防護、特殊環境中外絕緣放電特性及機理。

E-mail：cqsmwx@cqu.edu.cn

作者簡介

中圖分類號：TM864

收稿日期2015-03-04改稿日期2015-06-21

國家自然科學基金(51177182)、中央高校基本科研業務費(0209005206004)和國家創新研究群體科學基金(51321063)資助。