基于蒙特卡洛的動車組差動保護(hù)可靠性研究

孟建軍，趙文濤

(1.蘭州交通大學(xué)機(jī)電技術(shù)研究所，甘肅 蘭州，730070；2.甘肅省物流及運(yùn)輸裝備信息化工程技術(shù)研究中心，甘肅 蘭州，730070；3.甘肅省物流與運(yùn)輸裝備行業(yè)技術(shù)中心，甘肅 蘭州，730070)

1 引言

保護(hù)對于電力系統(tǒng)而言是非常重要的一環(huán)，而動車組作為高可靠性要求的產(chǎn)品保護(hù)更是必不可少。差動保護(hù)系統(tǒng)是動車組牽引電力系統(tǒng)的主保護(hù)，其安全可靠運(yùn)行對于整個(gè)動車組的運(yùn)營至關(guān)重要。然而我國鐵路部門對其采用的檢修方案仍是“定期維修”，造成了人力、物力、財(cái)力的巨大浪費(fèi)[1]，急需相關(guān)可靠性理論作為其維修方案的參考。

針對動車組差動保護(hù)系統(tǒng)組件數(shù)量多、離散程度高、結(jié)構(gòu)復(fù)雜[2]等特點(diǎn)，選用故障樹分析法來對其進(jìn)行可靠性分析，為表示動車組差動保護(hù)系統(tǒng)中的時(shí)序邏輯，引入動態(tài)邏輯門，故采用動態(tài)故障樹來對動車組差動保護(hù)系統(tǒng)進(jìn)行可靠性的量化與評估。

1991年，Dugan等最早提出了動態(tài)故障樹分析法的概念[3，4]，而國內(nèi)的熊小平等重點(diǎn)考慮了冗余的動態(tài)特性，建立動態(tài)故障樹對變電站保護(hù)系統(tǒng)進(jìn)行了可靠性分析[5]。對于動態(tài)故障樹的求解，1992年，Souza等首次將其轉(zhuǎn)化為馬爾可夫模型使用解析法求解[6，7]，近年來周廣林等對其進(jìn)行了優(yōu)化與擴(kuò)展，融合了三角模糊數(shù)并在礦井作業(yè)中得到了應(yīng)用[8]，但仍存在狀態(tài)爆炸、對于組件的故障率只能設(shè)為常數(shù)而實(shí)際中故障率隨時(shí)間變化等問題。蒙特卡洛法更多的是一種思想，利用大量的樣本數(shù)據(jù)進(jìn)行仿真，對其結(jié)果實(shí)行概率分析，進(jìn)而解決問題的思想[9]，而系統(tǒng)規(guī)模和狀態(tài)維數(shù)的爆炸性增長對蒙特卡洛法的運(yùn)算復(fù)雜度影響較小[10，11]。故本文對動態(tài)故障樹模型采用蒙特卡洛仿真方法求解，能夠很好的解決上述問題。

2 動車組差動保護(hù)系統(tǒng)簡介

基于基爾霍夫定理實(shí)現(xiàn)的動車組差動保護(hù)裝置用于動車組牽引系統(tǒng)的過壓、過流保護(hù)，牽引主電路正常工作時(shí)，牽引變壓器原邊輸入端和輸出端的電流互感器檢測電流差為零，即差動電流為零或在整定范圍以內(nèi)；當(dāng)差動電流值超過整定值時(shí)，觸發(fā)差動保護(hù)將牽引變壓器與故障電流斷開，實(shí)現(xiàn)牽引電力系統(tǒng)的保護(hù)[12]。其接線原理如圖1所示。

圖1 差動保護(hù)系統(tǒng)原理示意圖

3 動車組差動保護(hù)系統(tǒng)動態(tài)故障樹的建立

針對動車組差動保護(hù)系統(tǒng)，收集研究其故障機(jī)理[13]，如差動保護(hù)系統(tǒng)使用的電流互感器鐵芯類型為開口式，長期工作會導(dǎo)致其發(fā)生銹蝕、質(zhì)量下降、交合處氣隙增大、鐵芯磁導(dǎo)率下降等使得線路誤差增大超過限定值，進(jìn)而導(dǎo)致差動保護(hù)誤動[14]。確定以動車組差動保護(hù)系統(tǒng)故障為頂事件，以電流互感器鐵芯故障、電流互感器繞組故障、分接開關(guān)故障為中間事件建立三個(gè)故障子樹，然而在電流互感器中，無論是一次繞組還是二次繞組，若發(fā)生過電壓的情況，則一定會導(dǎo)致繞組局部過熱、局部場強(qiáng)增強(qiáng)、局部場強(qiáng)集中等故障的發(fā)生，這是經(jīng)典故障樹無法表達(dá)出來的時(shí)序性，故在上述故障樹中加入動態(tài)邏輯門繪制如圖2-6所示的動態(tài)故障樹。

圖2 差動保護(hù)系統(tǒng)動態(tài)故障樹

圖3 電流互感器鐵芯故障A子樹

圖4 電流互感器繞組故障B子樹

圖5 B子樹補(bǔ)充a子樹

圖6 分接開關(guān)故障C子樹

4 動車組差動保護(hù)系統(tǒng)可靠性分析

目前在可靠性分析理論中，認(rèn)為所收集的元件故障運(yùn)維數(shù)據(jù)是完全數(shù)據(jù)，即均符合收集數(shù)據(jù)要求的各項(xiàng)標(biāo)準(zhǔn)。但這與實(shí)際情況有很大的出入。收集到的故障數(shù)據(jù)并不都可以體現(xiàn)元件故障的特性。一些故障數(shù)據(jù)可能僅代表其壽命不小于該數(shù)值，無法得知確切壽命，這樣的數(shù)據(jù)稱為截尾數(shù)據(jù)。若并不是所有元件都從同一時(shí)間開始數(shù)據(jù)收集，而且在收集過程中，一些元件由于某種原因尚未故障中途撤離，如部分元件信息中途丟失、在現(xiàn)場收集時(shí)某些元件尚未失效等，這樣的截尾數(shù)據(jù)則成為隨機(jī)截尾數(shù)據(jù)。對于動車組差動保護(hù)這樣的壽命周期較長的安全系統(tǒng)，由于元件的定期維修或更換等原因，收集到的故障樣本數(shù)據(jù)更符合上述隨機(jī)截尾數(shù)據(jù)的特征，且為使本文所作的研究更具有普遍意義，認(rèn)為所使用的數(shù)據(jù)均為隨機(jī)截尾數(shù)據(jù)。

鑒于篇幅限制，本文將以分接開關(guān)故障C子樹中的電氣故障動態(tài)子樹為例，對動車組差動保護(hù)系統(tǒng)的可靠性研究進(jìn)行詳細(xì)的闡述。

4.1 元件故障分布參數(shù)估計(jì)

為更好地利用先驗(yàn)知識和現(xiàn)場組件故障數(shù)據(jù)統(tǒng)計(jì)的期望，本文采用了馬爾可夫蒙特卡洛方法，結(jié)合貝葉斯估計(jì)法和最大似然估計(jì)優(yōu)點(diǎn)的混合算法。不同的元件根據(jù)各自的先驗(yàn)知識與現(xiàn)場故障數(shù)據(jù)統(tǒng)計(jì)樣本，選取不同的分布模型進(jìn)行故障分布參數(shù)估計(jì)。對于電氣故障子樹而言，其各個(gè)底事件均為電氣常發(fā)故障，故障分布更符合二參數(shù)威布爾分布，故認(rèn)為電氣故障子樹中的底事件故障分布均為二參數(shù)威布爾分布。

二參數(shù)威布爾分布的故障概率分布函數(shù)為

(1)

其概率分布密度函數(shù)為

(2)

本文認(rèn)為所采集到的故障數(shù)據(jù)均為隨機(jī)結(jié)尾故障數(shù)據(jù)，可整理為

(t1，δ1)，(t2，δ2)，(t3，δ3)…(tn，δn)

(3)

δi=0為截尾數(shù)據(jù)，而δi=1為正常故障數(shù)據(jù)，本文以底事件62外部短路為例，其數(shù)據(jù)樣本如表1所示．

表1 底事件62外部短路的截尾故障樣本

其中*數(shù)據(jù)代表截尾數(shù)據(jù)。

鑒于直接采用貝葉斯方法的后驗(yàn)分布進(jìn)行統(tǒng)計(jì)推斷非常困難，本文建立一個(gè)穩(wěn)定分布與后驗(yàn)分布一致的馬爾科夫鏈，當(dāng)其收斂時(shí)便可以將其模擬值作為從后驗(yàn)分布中抽取的樣本。

電氣故障子樹中的底事件故障率均會隨著時(shí)間緩慢增加，如日久消耗等原因。極限情況為其故障率與時(shí)間呈線性關(guān)系，此時(shí)參數(shù)α=2；若其故障率與時(shí)間無關(guān)，此時(shí)參數(shù)α=1。故有理由相信對于參數(shù)α，先驗(yàn)概率在1.5兩側(cè)以同樣的速度遞減。此時(shí)α的先驗(yàn)邊緣概率估計(jì)函數(shù)為

(4)

故參數(shù)α的先驗(yàn)概率分布取以1.5為中點(diǎn)的三角分布。其抽樣公式如下

α=1+0.5(R1+R2)

(5)

其中R1，R2為在電腦上產(chǎn)生的[0，1]上均勻分布的隨機(jī)數(shù)。

對電氣故障動態(tài)子樹的底事件62外部短路發(fā)生時(shí)間樣本進(jìn)行統(tǒng)計(jì)分析，可得其故障時(shí)間期望為1000小時(shí)到1700小時(shí)之間，故該底事件的故障時(shí)間期望函數(shù)為

(6)

故參數(shù)β的估計(jì)值抽樣公式如下

(7)

其中R3為在電腦上產(chǎn)生的[0，1]上均勻分布的隨機(jī)數(shù)。

故電氣故障子樹的底事件62外部短路二參數(shù)威布爾分布的參數(shù)聯(lián)合先驗(yàn)概率為

(8)

前文提到當(dāng)構(gòu)建的馬爾科夫鏈?zhǔn)諗繒r(shí)進(jìn)行取樣便可模擬解決貝葉斯方法的實(shí)現(xiàn)困難之處，接受概率如下

(9)

其中π(α，β)為貝葉斯后驗(yàn)概率，正比于S(α，β)與g(α，β)的乘積，系數(shù)為標(biāo)準(zhǔn)化常量。

其中S(α，β)為二參數(shù)威布爾分布的最大似然函數(shù)，公式為

(10)

首先根據(jù)先驗(yàn)知識選取α為三角分布初始值為1.5，β為伽馬分布初始值為1350，迭代次數(shù)選為10000次，根據(jù)式(5)和(7)抽取抽樣值α′β′,以接受概率式(9)接受抽樣值，即若抽樣值與初始值的最大似然函數(shù)之比大于1則接受抽樣值，否則以比值的概率接受抽樣值，否則拒絕抽樣值。直到迭代結(jié)束。圖7為參數(shù)α的接受迭代抽樣值，圖8為參數(shù)β的接受迭代抽樣值。

圖7 參數(shù)α的迭代結(jié)果直方圖

圖8 參數(shù)β的迭代結(jié)果直方圖

在10000次迭代中，有5736次參數(shù)抽樣值可以被接受，故舍棄被拒絕的抽樣值，對迭代過程中所有參數(shù)的接受抽樣值取算數(shù)平均值即為α和β的估計(jì)值。結(jié)果為α=1.39，β=1431.1。

其它底事件的失效概率分布參數(shù)估計(jì)過程與62外部短路的估計(jì)方法相同，結(jié)果如表2所示。

表2 電氣故障動態(tài)子樹各底事件的概率分布函數(shù)參數(shù)估計(jì)值

4.2 動車組差動保護(hù)系統(tǒng)動態(tài)故障樹的蒙特卡洛仿真分析

本文所建立的動車組差動保護(hù)系統(tǒng)動態(tài)故障樹具有底事件數(shù)量大、故障發(fā)生的時(shí)間長、概率分布函數(shù)復(fù)雜等特點(diǎn)，故采用仿真固定時(shí)間步長法，即按照一個(gè)固定的仿真時(shí)間間隔作為時(shí)間增量來進(jìn)行仿真，以底事件的概率密度函數(shù)抽取底事件發(fā)生時(shí)間，與仿真時(shí)間做對比判斷元件的狀態(tài)，若仿真時(shí)間大于底事件發(fā)生時(shí)間則認(rèn)為該事件發(fā)生，在每個(gè)時(shí)間間隔推進(jìn)點(diǎn)上對各個(gè)事件的狀態(tài)進(jìn)行評估，更新整個(gè)動態(tài)故障樹的仿真狀態(tài)。

通過上面的計(jì)算已經(jīng)得到了各底事件發(fā)生時(shí)間分布的參數(shù)估計(jì)值，采用反函數(shù)隨機(jī)數(shù)抽取法來對發(fā)生故障的時(shí)間進(jìn)行抽樣，本文所設(shè)定的二參數(shù)威布爾分布的抽取公式為

(11)

其中r為在電腦上產(chǎn)生的[0，1]上均勻分布的隨機(jī)數(shù)。

將電氣故障動態(tài)故障子樹根據(jù)其時(shí)序邏輯轉(zhuǎn)化為布爾運(yùn)算函數(shù)

T=x52+x53+x54+x61+x62+…+x70+x71

(12)

輸入為底事件是否發(fā)生，輸出為頂事件是否發(fā)生。

輸入抽樣的各個(gè)底事件發(fā)生時(shí)間，按照固定步長法進(jìn)行仿真，每經(jīng)過一個(gè)時(shí)間間隔則與底事件發(fā)生時(shí)間進(jìn)行對比判斷底事件狀態(tài)，帶入布爾運(yùn)算函數(shù)中計(jì)算頂事件狀態(tài)。當(dāng)頂事件發(fā)生時(shí)記錄仿真時(shí)間，重復(fù)以上步驟直到設(shè)定的仿真次數(shù)。

在對電氣故障動態(tài)子樹頂事件進(jìn)行仿真時(shí)，設(shè)定仿真次數(shù)為10000，時(shí)間步長為2小時(shí)。

統(tǒng)計(jì)在整個(gè)仿真過程中頂事件發(fā)生的次數(shù)及仿真時(shí)間，對仿真結(jié)果進(jìn)行二參數(shù)威布爾分布曲線擬合，頂事件發(fā)生時(shí)間直方圖與函數(shù)擬合曲線如圖9所示。

圖9 電氣故障子樹頂事件發(fā)生概率分布函數(shù)直方圖和擬合曲線

故得到電氣故障子樹頂事件發(fā)生概率分布函數(shù)為

(13)

從圖中可以看出，頂事件發(fā)生次數(shù)在運(yùn)行1000小時(shí)左右最高，超過 3000小時(shí)次數(shù)逐漸趨近于0。

而整個(gè)動車組差動保護(hù)系統(tǒng)動態(tài)故障樹頂事件發(fā)生概率分布函數(shù)為：

(14)

4.3 動車組差動保護(hù)系統(tǒng)元件的重要度分析

從布爾運(yùn)算函數(shù)(12)來看，似乎在電氣故障動態(tài)故障子樹中，14個(gè)底事件對于頂事件發(fā)生的貢獻(xiàn)度是一樣的，而實(shí)際中不同底事件對于頂事件的重要度不相同也是不可能相同的，特定元件在系統(tǒng)的架構(gòu)、位置的不同可能造成其影響系統(tǒng)安全可靠運(yùn)行程度的不同，對動車組差動保護(hù)系統(tǒng)動態(tài)故障樹進(jìn)行重要度分析對系統(tǒng)設(shè)計(jì)、診斷及最優(yōu)化分析時(shí)有相當(dāng)大的作用，可根據(jù)各元件重要度決定系統(tǒng)檢查、維護(hù)及故障檢測執(zhí)行的先后順序，或是在系統(tǒng)改進(jìn)時(shí)改進(jìn)重要度較大的元件。

對于電氣故障動態(tài)故障子樹的重要度分析同樣使用蒙特卡洛仿真方法，即數(shù)學(xué)模型算法。系統(tǒng)內(nèi)包含諸多子樹，底事件數(shù)量較大，系統(tǒng)的故障時(shí)間函數(shù)中參數(shù)包含每個(gè)元件的故障時(shí)間，重要度便可使用系統(tǒng)發(fā)生故障時(shí)間對組件發(fā)生故障時(shí)間的偏微分來表示

(15)

在保持其它底事件的發(fā)生時(shí)間抽樣條件不變的情況下，僅改變某單個(gè)底事件的發(fā)生時(shí)間期望，然后使用蒙特卡洛仿真方法求解頂事件發(fā)生的時(shí)間期望，通過對比頂事件發(fā)生時(shí)間期望相對于初始狀態(tài)的改變來計(jì)算該底事件的重要度。

重要度的求解較為簡單，過程便不多贅述，動車組差動保護(hù)系統(tǒng)底事件重要度仿真結(jié)果如表3所示：

表3 動車組差動保護(hù)系統(tǒng)底事件重要度

由上述重要度仿真分析結(jié)果可知，對于動車組差動保護(hù)系統(tǒng)而言，底事件發(fā)生特性變化對頂事件發(fā)生概率分布重要度較大的是底事件1銹蝕、39諧振以及62外部短路。總體而言，繞組對于動車組差動保護(hù)系統(tǒng)重要度最大，分接開關(guān)的重要度最小。設(shè)法降低上述重要度較高的底事件發(fā)生概率可以有效提高系統(tǒng)的可靠性，在對動車組差動保護(hù)系統(tǒng)進(jìn)行系統(tǒng)檢查、維護(hù)及故障檢測時(shí)可參考本文重要度的仿真分析結(jié)果，優(yōu)先對繞組進(jìn)行檢測，重點(diǎn)關(guān)注是否發(fā)生了銹蝕、諧振以及外部短路等故障。

5 結(jié)束語

本文結(jié)合動車組差動保護(hù)系統(tǒng)的故障機(jī)理，完成了動態(tài)故障樹的構(gòu)建，采用了馬爾可夫蒙特卡洛方法結(jié)合貝葉斯估計(jì)法和最大似然估計(jì)優(yōu)點(diǎn)的混合算法來對動態(tài)故障樹各個(gè)底事件的概率分布函數(shù)的參數(shù)值進(jìn)行估計(jì)，構(gòu)造頂事件的布爾函數(shù)，利用固定時(shí)間步長仿真法進(jìn)行蒙特卡洛仿真分析求得頂事件發(fā)生的概率分布函數(shù)，證明了動車組差動保護(hù)系統(tǒng)可靠性較高，通過改變單一變量的方法系統(tǒng)各底事件對于頂事件的重要度進(jìn)行分析，得出了繞組的重要度較高，在進(jìn)行系統(tǒng)檢查、維護(hù)及故障檢測時(shí)可參考本文重要度的仿真分析結(jié)果，優(yōu)先對繞組進(jìn)行檢測，重點(diǎn)關(guān)注是否發(fā)生了銹蝕、諧振以及外部短路等故障的結(jié)論，可作為制訂動車組差動保護(hù)系統(tǒng)檢修策略的參考。