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

(輸配電裝備及系統(tǒng)安全與新技術(shù)國(guó)家重點(diǎn)實(shí)驗(yàn)室(重慶大學(xué))　重慶　400044)

?

實(shí)測(cè)鐵磁諧振時(shí)間序列的非線性動(dòng)力學(xué)分析

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

(輸配電裝備及系統(tǒng)安全與新技術(shù)國(guó)家重點(diǎn)實(shí)驗(yàn)室(重慶大學(xué))重慶400044)

摘要在僅有實(shí)測(cè)單一電壓時(shí)間序列的情況下，為準(zhǔn)確識(shí)別出電力系統(tǒng)鐵磁諧振類型，運(yùn)用坐標(biāo)延遲方法重構(gòu)與原系統(tǒng)拓?fù)湟饬x上等價(jià)的相空間。采用相平面、龐加萊截面和關(guān)聯(lián)維數(shù)3種基于重構(gòu)相空間的非線性動(dòng)力學(xué)分析工具表達(dá)電壓時(shí)間序列的運(yùn)動(dòng)特征，由此識(shí)別其所屬的鐵磁諧振類型。對(duì)3例典型的中性點(diǎn)不接地系統(tǒng)電磁式電壓互感器鐵磁諧振實(shí)測(cè)時(shí)間序列進(jìn)行了動(dòng)力學(xué)特征分析。其中兩例時(shí)間序列的相平面軌跡、龐加萊截面均呈現(xiàn)周期運(yùn)動(dòng)特性，關(guān)聯(lián)維數(shù)估計(jì)值分別為1.015 0±0.003 9、1.006 1±0.000 6，由此判別它們作周期運(yùn)動(dòng)；另一例時(shí)間序列的相平面軌跡、龐加萊截面均未呈現(xiàn)周期或準(zhǔn)周期運(yùn)動(dòng)特性，關(guān)聯(lián)維數(shù)估計(jì)值為2.300 2±0.061 2，由此判別其作混沌運(yùn)動(dòng)。實(shí)測(cè)電壓時(shí)間序列受到多種因素的影響，綜合以上3種分析方法所表達(dá)的特征，對(duì)所發(fā)生的鐵磁諧振類型作出的判斷更準(zhǔn)確和更有說(shuō)服力。

關(guān)鍵詞：鐵磁諧振時(shí)間序列非線動(dòng)力學(xué)相空間重構(gòu)關(guān)聯(lián)維數(shù)電磁式電壓互感器

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引言

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

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

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

因一個(gè)變量隨時(shí)間的變化隱含了整個(gè)系統(tǒng)的運(yùn)動(dòng)規(guī)律，根據(jù)F.Takens提出的嵌入定理[11]，從一維時(shí)間序列可重構(gòu)一個(gè)與原動(dòng)力系統(tǒng)在拓?fù)湟饬x上等價(jià)的相空間。因此，可對(duì)單一變量的時(shí)間序列運(yùn)用相空間重構(gòu)方法構(gòu)造相空間，通過(guò)分析重構(gòu)相空間軌線的分布或結(jié)構(gòu)得到原系統(tǒng)的動(dòng)力學(xué)特征。

本文在僅有某一中性點(diǎn)不接地系統(tǒng)PT鐵磁諧振實(shí)測(cè)電壓時(shí)間序列的情況下，運(yùn)用相空間重構(gòu)方法獲得系統(tǒng)相空間。在此基礎(chǔ)上，運(yùn)用相平面、龐加萊截面和關(guān)聯(lián)維數(shù)3種基于相空間重構(gòu)的非線性動(dòng)力學(xué)分析工具表達(dá)系統(tǒng)的動(dòng)力學(xué)特征，以識(shí)別鐵磁諧振類型。

1非線性動(dòng)力學(xué)分析方法

1.1相空間重構(gòu)

系統(tǒng)相空間是由狀態(tài)變量xi(i=1,2,…,n)張成的空間Rn， 相空間的一個(gè)點(diǎn)對(duì)應(yīng)系統(tǒng)的一個(gè)狀態(tài)。N.H.Parkard等[12]于1980年提出了基于單一變量時(shí)間序列的兩種相空間重構(gòu)方法：導(dǎo)數(shù)重構(gòu)法和坐標(biāo)延遲重構(gòu)法。因?qū)崪y(cè)數(shù)據(jù)混有噪聲，對(duì)其求取的微分誤差較大，因此實(shí)際應(yīng)用中主要使用坐標(biāo)延遲重構(gòu)法進(jìn)行一維時(shí)間序列的相空間重構(gòu)。

1.1.1坐標(biāo)延遲重構(gòu)法

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

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

(1)

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

由此，重構(gòu)相空間的軌線分布或結(jié)構(gòu)便反映了系統(tǒng)的運(yùn)動(dòng)特征。將重構(gòu)相空間的某兩個(gè)向量的取值序列分別作為橫、縱坐標(biāo)繪制在平面中就構(gòu)成重構(gòu)相平面。系統(tǒng)作周期運(yùn)動(dòng)時(shí)，其相平面軌跡是一條不斷重復(fù)的閉合曲線[2]；而作混沌運(yùn)動(dòng)時(shí)，相平面軌跡永不重復(fù)[1]。

1.1.2坐標(biāo)延遲重構(gòu)參數(shù)選取

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

本文在確定相空間重構(gòu)參數(shù)時(shí)采用第一種觀點(diǎn)，即τ和m的值采用不同的方法確定。采用互信息法獲取最佳τ值，因互信息刻畫的是兩個(gè)變量間的普遍(或總體)相關(guān)性，適用于非線性系統(tǒng)。本文采用的計(jì)算過(guò)程無(wú)需事先確定最優(yōu)m值，但通過(guò)計(jì)算一組m值下的關(guān)聯(lián)維數(shù)可間接獲取最優(yōu)m值(即G-P法確定最優(yōu)m值，實(shí)現(xiàn)方法見(jiàn)1.3.1節(jié))。

1.2龐加萊截面

龐加萊截面是相平面軌跡的一種簡(jiǎn)化表達(dá)，是對(duì)相平面軌跡以一定的頻率采樣獲取的相點(diǎn)圖。

若對(duì)相平面軌跡采樣的頻率為系統(tǒng)外加激勵(lì)函數(shù)的頻率(如電源頻率)，則對(duì)作周期運(yùn)動(dòng)的系統(tǒng)，在龐加萊截面上表現(xiàn)為一個(gè)不動(dòng)點(diǎn)或是多個(gè)孤立點(diǎn)[1]；對(duì)作準(zhǔn)周期運(yùn)動(dòng)的系統(tǒng)，龐加萊截面是一條分布著無(wú)限個(gè)不同點(diǎn)的封閉曲線[6，7]；而當(dāng)系統(tǒng)作混沌運(yùn)動(dòng)時(shí)，龐加萊截面呈現(xiàn)為局限于平面某一特定區(qū)域的一些隨機(jī)分布的相點(diǎn)集合[1]。

1.3關(guān)聯(lián)維數(shù)

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

1.3.1關(guān)聯(lián)維數(shù)G-P算法

給定一尺度r和嵌入維數(shù)m，定義相空間距離小于r的點(diǎn)對(duì)數(shù)與所有點(diǎn)對(duì)數(shù)的比例為關(guān)聯(lián)積分Cm(r)，即

(2)

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

(3)

對(duì)于無(wú)限長(zhǎng)序列(n→∞)和足夠小的r值，若Cm(r)與r呈冪律關(guān)系，即Cm(r)∝rD2， 則定義關(guān)聯(lián)維數(shù)D2為

(4)

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

min{‖Vi-Vj‖}

(5)

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

因此，計(jì)算關(guān)聯(lián)維數(shù)采取的常規(guī)步驟為：

1)在式(5)的范圍內(nèi)取一尺度序列ri(i=1,2,…,l)，計(jì)算與ri對(duì)應(yīng)的關(guān)聯(lián)積分序列Cm(ri)。

2)求出序列Cm(ri)和ri的對(duì)數(shù)序列l(wèi)nCm(r)i和lnri，以序列l(wèi)nri為橫坐標(biāo)、序列l(wèi)nCm(r)i為縱坐標(biāo)繪制雙對(duì)數(shù)曲線。

3)將雙對(duì)數(shù)曲線線性區(qū)(或無(wú)標(biāo)度區(qū))的斜率作為D2的估計(jì)值。

因雙對(duì)數(shù)曲線或多或少呈現(xiàn)S形，很難通過(guò)目測(cè)或經(jīng)驗(yàn)確定其線性區(qū)域。文獻(xiàn)[26]指出一個(gè)可信的無(wú)標(biāo)度區(qū)必須是清晰可辨的，表現(xiàn)在lnCm(r)對(duì)lnr曲線的局部斜率dsm對(duì)lnr曲線存在隨m增大而趨于飽和的平坦區(qū)域。

為使估算所得D2盡量準(zhǔn)確，本文計(jì)算關(guān)聯(lián)維數(shù)采取的步驟與常規(guī)的不同之處為：

1)計(jì)算一組m(m=2~10)的對(duì)數(shù)序列l(wèi)nCm(r)i。

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

3)若dsm對(duì)lnr曲線存在隨m值增大而趨于飽和的平坦區(qū)域，則將各m值下dsm對(duì)lnr曲線的平坦區(qū)域作為無(wú)標(biāo)度區(qū)，并將此區(qū)域?qū)?yīng)的數(shù)據(jù)點(diǎn)對(duì)(lnCm(r)，lnr)進(jìn)行直線擬合得到序列D2_m，由此獲取D2估計(jì)值。同時(shí)，D2_m趨于飽和時(shí)對(duì)應(yīng)的m值可作為該時(shí)間序列重構(gòu)的最優(yōu)m值，即G-P法確定最優(yōu)m值。

局部斜率dsm的計(jì)算式為

(6)

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

對(duì)于隨機(jī)過(guò)程的時(shí)間序列，其D2隨m的增大而增大，不會(huì)趨于飽和；而確定性系統(tǒng)的D2隨m的增大而逐漸趨于飽和。不同的D2飽和值對(duì)應(yīng)于不同的系統(tǒng)運(yùn)動(dòng)狀態(tài)[27]：D2=1，系統(tǒng)作周期運(yùn)動(dòng)；D2=2，系統(tǒng)作準(zhǔn)周期或擬周期運(yùn)動(dòng)；D2>2或D2不為整數(shù)，系統(tǒng)作混沌運(yùn)動(dòng)。

1.3.2關(guān)聯(lián)維數(shù)G-P優(yōu)化算法

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

2典型鐵磁諧振時(shí)間序列的非線性動(dòng)力學(xué)分析

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

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

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

1)時(shí)間序列1

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

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

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

圖3　時(shí)間序列1的重構(gòu)相平面圖Fig.3　Reconstructed phase plane plot of time series 1

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

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

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

圖6　時(shí)間序列1的dsm對(duì)lnr曲線Fig.6　Curves dsm versus lnr of time series 1

m=2～10時(shí)，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時(shí)，D2_m趨于穩(wěn)定，由此估算所得D2=1.015 0±0.003 9。該序列重構(gòu)所需的最小嵌入維數(shù)為3。

時(shí)間序列1的D2相對(duì)于整數(shù)1的誤差低于2%，非常小，因此其D2理論值應(yīng)為整數(shù)1，作周期運(yùn)動(dòng)。這與從相平面軌跡和龐加萊截面得出的結(jié)論一致。

對(duì)與時(shí)間序列1同一時(shí)間同一位置獲取的B相、C相數(shù)據(jù)也進(jìn)行了動(dòng)力學(xué)特征分析，其運(yùn)動(dòng)模式與時(shí)間序列1的相同，D2分別為1.007 6±0.008 7，1.028 8±0.011 5。

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

2)時(shí)間序列2

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

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

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

圖8　時(shí)間序列2的重構(gòu)相平面圖Fig.8　Reconstructed phase plane plot of time series 2

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

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

對(duì)時(shí)間序列1和2的分析表明，測(cè)量誤差的存在導(dǎo)致龐加萊截面上應(yīng)重疊的相點(diǎn)未完全重疊，給分析帶來(lái)困難。因此，對(duì)周期比較長(zhǎng)的時(shí)間序列，運(yùn)用龐加萊截面時(shí)需非常謹(jǐn)慎，時(shí)間序列數(shù)據(jù)要足夠長(zhǎng)，至少要包含2個(gè)運(yùn)動(dòng)周期。

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

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

圖11　時(shí)間序列2的dsm對(duì)lnr曲線Fig.11　Curves dsm versus lnr of time series 2

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

時(shí)間序列2的D2相對(duì)于整數(shù)1的誤差低于0.7%，非常小，因此其D2理論值應(yīng)為整數(shù)1，作周期運(yùn)動(dòng)。這與從相平面軌跡和龐加萊截面得出的結(jié)論一致。

對(duì)與時(shí)間序列2同一時(shí)間同一位置獲取的B相、C相數(shù)據(jù)也進(jìn)行了動(dòng)力學(xué)特征分析，其運(yùn)動(dòng)模式與時(shí)間序列2的相同，D2分別為1.001 3±0.001 0，1.018 4±0.005 3。

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

3)時(shí)間序列3

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

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

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

圖13　時(shí)間序列3的重構(gòu)相平面圖Fig.13　Reconstructed phase plane plot of time series 3

與圖13對(duì)應(yīng)的龐加萊截面如圖14所示。圖14的相點(diǎn)呈不均勻的離散分布，未形成一條近似閉曲線，這表明時(shí)間序列3不是作準(zhǔn)周期運(yùn)動(dòng)。因圖14中的相點(diǎn)數(shù)較少，未呈現(xiàn)自相似的分層結(jié)構(gòu)，因此不能判斷時(shí)間序列3是否作混沌運(yùn)動(dòng)。

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

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

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

圖16　時(shí)間序列3的dsm對(duì)lnr曲線Fig.16　Curves dsm versus lnr of time series 3

m=2～10時(shí)，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時(shí)，D2_m趨于穩(wěn)定，由此估算可得D2=2.300 2±0.061 2。該序列重構(gòu)所需的最小嵌入維數(shù)為5。

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

對(duì)與時(shí)間序列3同一時(shí)間同一位置獲取的B相、C相數(shù)據(jù)也進(jìn)行了動(dòng)力學(xué)特征分析，其運(yùn)動(dòng)模式與時(shí)間序列3的相同，D2分別為2.223 6±0.043 6，2.408 3±0.082 4。

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

3結(jié)論

1)為準(zhǔn)確識(shí)別實(shí)測(cè)鐵磁諧振一維電壓時(shí)間序列所屬的諧振類型，本文運(yùn)用相空間重構(gòu)方法構(gòu)造系統(tǒng)相空間。在此基礎(chǔ)上，采用相平面、龐加萊截面和關(guān)聯(lián)維數(shù)3種非線性動(dòng)力學(xué)分析方法，以不同的形式表達(dá)其動(dòng)力學(xué)特征，由此對(duì)其所屬的鐵磁諧振類型作出了準(zhǔn)確判斷。

2)對(duì)3例典型的中性點(diǎn)不接地系統(tǒng)PT鐵磁諧振實(shí)測(cè)時(shí)間序列進(jìn)行了動(dòng)力學(xué)特征分析。其中兩例時(shí)間序列的相平面軌跡、龐加萊截面均呈現(xiàn)周期運(yùn)動(dòng)特性，關(guān)聯(lián)維數(shù)估計(jì)值分別為1.015 0±0.003 9、1.006 1±0.000 6，由此判斷它們作周期運(yùn)動(dòng)。另一例時(shí)間序列的相平面軌跡、龐加萊截面均未呈現(xiàn)周期或準(zhǔn)周期運(yùn)動(dòng)特性，關(guān)聯(lián)維數(shù)估計(jì)值為2.300 2±0.061 2，由此判斷其作混沌運(yùn)動(dòng)。

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

4)中性點(diǎn)不接地系統(tǒng)PT鐵磁諧振屬于零序性質(zhì)，發(fā)生鐵磁諧振時(shí)對(duì)任一相電壓數(shù)據(jù)進(jìn)行非線性動(dòng)力學(xué)分析，均能對(duì)其類型作出準(zhǔn)確判斷。

5)本文的研究為鐵磁諧振類型的在線自動(dòng)識(shí)別和諧振過(guò)電壓的抑制策略研究奠定了一定的基礎(chǔ)，為實(shí)測(cè)鐵磁諧振數(shù)據(jù)的識(shí)別研究提供了新思路，豐富了鐵磁諧振的研究實(shí)例。

參考文獻(xiàn)

[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]周麗霞，尹忠東，鄭立.配網(wǎng)PT鐵磁諧振機(jī)理與抑制的試驗(yàn)研究[J].電工技術(shù)學(xué)報(bào)，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]張博，魯鐵成，杜曉磊.中性點(diǎn)接地系統(tǒng)鐵磁諧振非線性動(dòng)力學(xué)分析[J].高電壓技術(shù)，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]司馬文霞，鄭哲人，楊慶，等.用參數(shù)不匹配混沌系統(tǒng)的脈沖同步方法抑制鐵磁諧振過(guò)電壓[J].電工技術(shù)學(xué)報(bào)，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，F(xiàn)armer J D，et al.Geometry from a time series[J].Physical Review Letters，1980，45：712-715.

[13]行鴻彥，龔平，徐偉.嵌入窗方法確定混沌系統(tǒng)重構(gòu)參數(shù)的仿真研究[J].系統(tǒng)仿真學(xué)報(bào)，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]劉秉正，彭建華.非線性動(dòng)力學(xué)[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]張雨.時(shí)間序列的混沌和符號(hào)分析及實(shí)踐[M].長(zhǎng)沙：國(guó)防科技大學(xué)出版社，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變電站過(guò)電壓在線監(jiān)測(cè)系統(tǒng)及其波形分析[J].高電壓技術(shù)，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年生，博士，講師，研究方向?yàn)殡娏ο到y(tǒng)過(guò)電壓。

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

司馬文霞女，1965年生，教授，博士生導(dǎo)師，研究方向?yàn)殡娏ο到y(tǒng)的防雷與過(guò)電壓防護(hù)、特殊環(huán)境中外絕緣放電特性及機(jī)理。

E-mail：cqsmwx@cqu.edu.cn

作者簡(jiǎn)介

中圖分類號(hào)：TM864

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

國(guó)家自然科學(xué)基金(51177182)、中央高?；究蒲袠I(yè)務(wù)費(fèi)(0209005206004)和國(guó)家創(chuàng)新研究群體科學(xué)基金(51321063)資助。