基于離散時間貝葉斯網(wǎng)絡(luò)的動態(tài)故障樹分析的改良方法

蘭　杰, 袁宏杰, 夏　靜

(北京航空航天大學(xué)可靠性與系統(tǒng)工程學(xué)院, 北京 100191)

0　引　言

故障樹分析作為極其重要的系統(tǒng)可靠性分析方法,已廣泛應(yīng)用于電子通信、醫(yī)療、能源、軍事和航空航天等領(lǐng)域。傳統(tǒng)靜態(tài)故障樹分析方法只重視部件間的與、或、表決等靜態(tài)特性,而忽略了部件間的優(yōu)先關(guān)系、順序關(guān)系、功能相關(guān)性等動態(tài)失效特性,因此在分析具有動態(tài)失效特征的系統(tǒng)可靠性時存在較大誤差[1]。為了彌補(bǔ)靜態(tài)故障樹分析方法在動態(tài)特性上建模功能的不足,文獻(xiàn)[2]提出了動態(tài)故障樹分析方法,以Markov理論為基礎(chǔ),建立了動態(tài)故障樹模型。但Markov方法只適用于指數(shù)分布[3-4],而且存在狀態(tài)空間爆炸問題,難以適用于大型動態(tài)故障樹的分析。為避開Markov模型來解決動態(tài)故障樹,文獻(xiàn)[5]提出了基于復(fù)合梯形積分的計算方法,但該方法仍然存在組合爆炸的問題。

貝葉斯網(wǎng)絡(luò)技術(shù)(Bayesian networks, BN)因具有系統(tǒng)建模、推理和診斷的優(yōu)點(diǎn)[6-8],廣泛應(yīng)用于可靠性分析領(lǐng)域。2001年,文獻(xiàn)[9]提出了靜態(tài)故障樹轉(zhuǎn)化為BN的方法,并給出BN方法的實(shí)例。2005年,文獻(xiàn)[10-11]提出了基于BN解決動態(tài)故障樹的方法,在建模和分析能力方面有著顯著優(yōu)勢,有效地避免了狀態(tài)爆炸問題。文獻(xiàn)[12-13]提出了動態(tài)BN分析方法,并給出了優(yōu)先與門、溫備份門等動態(tài)門轉(zhuǎn)化為動態(tài)BN的方法。2007年,為了克服靜態(tài)或均勻離散帶來的精度差、效率低等問題,文獻(xiàn)[14]提出了動態(tài)離散化的方法。2010年,文獻(xiàn)[15]提出了聯(lián)合動態(tài)離散化算法和基于事件的BN建模的方法,來提高服從任何分布的復(fù)雜系統(tǒng)的計算精度,但計算和迭代過程復(fù)雜。

傳統(tǒng)的離散時間BN(discrete time Bayesian networks, DTBN)方法的計算時間、精度和效率受離散的時間分段n影響極大。針對此問題,本文基于復(fù)合梯形積分方法[5],分析討論了傳統(tǒng)方法計算的誤差來源,并提出了新的動態(tài)門轉(zhuǎn)化方法。與傳統(tǒng)方法和Markov方法相比,新的轉(zhuǎn)化方法補(bǔ)償了傳統(tǒng)方法的計算誤差,提高了結(jié)果的計算精度和計算效率,適用于服從各種常見分布的復(fù)雜系統(tǒng)。最后以某測量系統(tǒng)為例,建立了系統(tǒng)的動態(tài)故障樹和BN,驗證了改良方法的可行性和高效性。

1　DTBN

DTBN是一個有向無環(huán)圖,主要包括兩部分:一部分是由有向邊和節(jié)點(diǎn)構(gòu)成的網(wǎng)絡(luò);另一部分是與每個節(jié)點(diǎn)相關(guān)的條件概率表(conditional probability table, CPT)。

將動態(tài)故障樹中每個基本事件、靜態(tài)或動態(tài)邏輯門和頂事件,用單獨(dú)的BN節(jié)點(diǎn)來代替;用有向邊來代表兩個節(jié)點(diǎn)之間存在邏輯關(guān)系。如果節(jié)點(diǎn)A到節(jié)點(diǎn)B之間有一條A指向B的有向邊,則節(jié)點(diǎn)A是B的父節(jié)點(diǎn),B是A的子節(jié)點(diǎn)。對于給定節(jié)點(diǎn)Vi,擁有父節(jié)點(diǎn)pa(Vi),則:P(Vi|V1,V2,…,Vi-1)=P(Vi|pa(Vi))。

每一個節(jié)點(diǎn)都有其相關(guān)的條件概率表CPT。由BN的條件獨(dú)立性可知,每一個節(jié)點(diǎn)的CPT可以由P(Vi|pa(Vi))得到,它描述了節(jié)點(diǎn)之間的邏輯關(guān)系。

把整個任務(wù)工作時間[0,t]分成n個均勻的區(qū)間,每個區(qū)間的長度為Δ=t/n,整個任務(wù)時間被劃分n+1個部分:[0,Δ],[Δ,2Δ],[2Δ,3Δ], …，[(n-1)Δ,nΔ],[nΔ,+∞]。設(shè)節(jié)點(diǎn)E=[(x-1)Δ,xΔ],則節(jié)點(diǎn)E在[(x-1)Δ,xΔ]區(qū)間內(nèi)失效,或稱節(jié)點(diǎn)E處于狀態(tài)x。若頂事件ET失效,則必然發(fā)生在這(n+1)個時間區(qū)間中。頂事件ET在任務(wù)時間[0,T]內(nèi)失效的概率為

2　動態(tài)故障樹轉(zhuǎn)化的誤差分析與補(bǔ)償

2.1　優(yōu)先與門轉(zhuǎn)化分析

2.1.1傳統(tǒng)的優(yōu)先與門轉(zhuǎn)化方法

優(yōu)先與門有兩個輸入事件,當(dāng)且僅當(dāng)A失效,然后B失效,其輸出C才會失效。優(yōu)先與門的輸入可以是基本事件,也可以是一個動態(tài)故障樹的子樹。優(yōu)先門向DTBN的轉(zhuǎn)化圖如圖1所示。

圖1　優(yōu)先與門轉(zhuǎn)化為DTBN

將[0,t]分為n段,每段區(qū)間為Δ=t/n,則當(dāng)A的概率密度為fA(t)時,A在狀態(tài)x的失效概率為

(1)

式中,A可以服從任何概率分布。

不妨設(shè)A的失效概率服從指數(shù)分布,失效率為λa,則概率密度fA(t)=λae-λat。B的失效概率也服從指數(shù)分布,失效率為λb,概率密度fB(t)=λbe-λbt?？傻肁、B在狀態(tài)x下的失效概率為

(2)

Pb,x=P(B=[(x-1)Δ,xΔ])=(eλbΔ-1)e-λbxΔ

(3)

當(dāng)A、B在時間段[T,∞)內(nèi)失效,即A、B位于n+1狀態(tài)下的失效概率為

(4)

(5)

由優(yōu)先與門性質(zhì)可知,當(dāng)且僅當(dāng)A失效,然后B失效,其輸出C才會失效。傳統(tǒng)的貝葉斯轉(zhuǎn)化方法設(shè)定:如果節(jié)點(diǎn)A和節(jié)點(diǎn)B在同一個時間區(qū)間[(x-1)Δ,xΔ]內(nèi)失效,即A、B同時處于狀態(tài)x時,節(jié)點(diǎn)C的失效概率直接為0，或者為1[16]。

若將A、B、C的工作時間[0,T]分為n段,則A的CPT為1×(n+1)行列的表,B的CPT為1×(n+1)行列的表,C的CPT為(n+1)×(n+1)×(n+1)的表。節(jié)點(diǎn)C的條件概率分布式為

(6)

當(dāng)A、B中有任一個節(jié)點(diǎn)位于n+1狀態(tài),即在時間區(qū)間[T,∞)內(nèi)失效,C節(jié)點(diǎn)在時間區(qū)間[T,∞)內(nèi)失效的條件概率為1。

對單一優(yōu)先門,若只需求得節(jié)點(diǎn)C在任務(wù)時間[0,T]的失效概率和可靠度時,可簡化C的CPT的規(guī)模,設(shè)定節(jié)點(diǎn)A、B工作時間[0,∞)分為n+1段,節(jié)點(diǎn)C的時間區(qū)間分為[0,T]、[T,∞)兩段。此時節(jié)點(diǎn)C的CPT為2×(n+1)×(n+1)的表,規(guī)模大幅度變小,且計算結(jié)果與第3.1節(jié)方法的結(jié)果一樣。節(jié)點(diǎn)C在任務(wù)區(qū)間[0,T]失效,則條件概率表分布式為

PT,x,y=P(C=[0,T]|

A=[(x-1)Δ,xΔ],B=[(y-1)Δ,yΔ])=

(7)

PT,x,n+1=P(C=[0,T]|A=

[(x-1)Δ,xΔ],B=[T,∞))=

PT,n+1,y=P(C=[0,T]|A=[T,∞),

B=[(y-1)Δ,yΔ])=

PT,n+1,n+1=P(C=[0,T]|A=[T,∞),

B=[T,∞))=0

(8)

節(jié)點(diǎn)C在時間區(qū)間[T,∞)失效,則條件概率表分布式為

P(C=[T,∞)|A=[(x-1)Δ,xΔ],

B=[(y-1)Δ,yΔ])=1-PT,x,y

(9)

節(jié)點(diǎn)C在任務(wù)區(qū)間[0,T]內(nèi)失效的條件概率表為如表1所示,同理也可得節(jié)點(diǎn)C在區(qū)間[T,∞)的CPT。

表1　優(yōu)先與門節(jié)點(diǎn)C的CPT

2.1.2誤差分析

在傳統(tǒng)的BN轉(zhuǎn)化方法下,優(yōu)先門輸出點(diǎn)C在任務(wù)區(qū)間[0,T]的失效概率Pc為

(10)

由表1可得,式(10)可轉(zhuǎn)化為

(11)

設(shè)ta為A的失效點(diǎn),tb為B的失效點(diǎn),根據(jù)復(fù)合梯形積分方法可得優(yōu)先與門的輸出頂事件C在任務(wù)區(qū)間[0,T]的失效概率Pc1為

(12)

Pc1是優(yōu)先與門的頂事件的失效概率的真值,則傳統(tǒng)的BN轉(zhuǎn)化方法的計算結(jié)果的誤差ε1為ε1=|Pc-Pc1|。

其中可用將任務(wù)區(qū)間[0,T]分為n段的方法,對式(12)進(jìn)行轉(zhuǎn)化,則

(13)

由式(11)和式(13)得

(14)

利用夾逼定理可得

(15)

當(dāng)n取值小時,傳統(tǒng)貝葉斯轉(zhuǎn)化方法的誤差會很大;當(dāng)n→∞時,誤差變小并趨近于0。驗證了當(dāng)n取值越大時傳統(tǒng)貝葉斯轉(zhuǎn)化方法的精度越高,誤差越小,同時可知當(dāng)n越大時,計算時間變長,各個部件的CPT的規(guī)模變大。

2.1.3改良方法

根據(jù)前面得出的傳統(tǒng)方法的計算誤差和復(fù)合梯形積分方法,利用得到的誤差直接修正節(jié)點(diǎn)C的CPT,改良優(yōu)先與門的轉(zhuǎn)化方法。改良方法在動態(tài)門中不存在重復(fù)事件的情況下完全補(bǔ)償?shù)袅擞嬎阏`差,進(jìn)而減小了存在重復(fù)事件情形下的動態(tài)門誤差。

設(shè)定如果節(jié)點(diǎn)A和節(jié)點(diǎn)B在同一個時間區(qū)間[(x-1)Δ,xΔ]內(nèi)失效時,節(jié)點(diǎn)C的失效概率不直接為0,即

PT,x=

P(C=[0,T]|A=[(x-1)Δ,xΔ],B=[(x-1)Δ,xΔ])=

(16)

特別對于服從指數(shù)分布的系統(tǒng),即

(17)

則改良后節(jié)點(diǎn)C在任務(wù)區(qū)間[0,T]失效,條件概率表分布式為

PT,x,y=P(C=[0,T]|

A=[(x-1)Δ,xΔ],B=[(y-1)Δ,yΔ])=

(18)

改良節(jié)點(diǎn)C在任務(wù)區(qū)間[0,T]內(nèi)失效的條件概率表如表2所示。

表2　改良后的優(yōu)先與門節(jié)點(diǎn)C的CPT

對于動態(tài)門中不存在重復(fù)事件的情況,改良后的轉(zhuǎn)化方法,在n值很小時能保證結(jié)果的正確性,即ε1=0。在動態(tài)門中存在重復(fù)事件時,會存在比傳統(tǒng)方法更加微小的誤差。

2.2　備份門轉(zhuǎn)化分析

根據(jù)備份件在儲備狀態(tài)時能量消耗的不同,在備份門可分為冷備份門(cold spare gate,CSP)、溫備份門(warm spare gate,WSP)和熱備份門(hot spare gate,HSP)。當(dāng)主件處于工作狀態(tài)時,備份單元處于儲備狀態(tài);當(dāng)且僅當(dāng)主件失效后,備份單元才進(jìn)入工作狀態(tài)。對于冷備份門,備份單元在進(jìn)入工作狀態(tài)前視為無能量耗損,即失效率為零。

2.2.1傳統(tǒng)的備份門轉(zhuǎn)化方法

由CSP性質(zhì)可知,當(dāng)且僅當(dāng)主件失效后,備份單元才進(jìn)入工作狀態(tài),備份單元在儲備狀態(tài)時失效率為零。CSP向DTBN轉(zhuǎn)化如圖2所示。

圖2　冷備份轉(zhuǎn)化為DTBN

當(dāng)x

Px,y=P(B=[(y-1)Δ,yΔ]|A=[(x-1)Δ,xΔ])=

(19)

對于服從指數(shù)分布的系統(tǒng),不妨設(shè)A的失效率為λa,B的失效率為λb,概率密度分別為:fA(t)=λae-λat,fB(t)=λbe-λbt,則

Px,n+1=P(B=[T,∞)|A=[(x-1)Δ,xΔ])=

(20)

B的條件概率分布式為

P(B=[(y-1)Δ,yΔ]|A=[(x-1)Δ,xΔ])=

(21)

P(B=[(y-1)Δ,yΔ]|A=[T,∞))=0

(22)

P(B=[T,∞)|A=[T,∞))=1

(23)

節(jié)點(diǎn)B的CPT如表3所示。

表3　備份門節(jié)點(diǎn)B的CPT

節(jié)點(diǎn)C的條件概率分布式為

P(C=[(z-1)Δ,zΔ]|B=[(y-1)Δ,yΔ])=

(24)

2.2.2誤差分析

在傳統(tǒng)轉(zhuǎn)化方法下,優(yōu)先門輸出點(diǎn)C在任務(wù)區(qū)間[0,T]的失效概率Pc為

Pc=P(C=[0,T])=

(25)

由表3可知,式(25)可轉(zhuǎn)化為

(26)

設(shè)ta為A的失效點(diǎn),tb為B的失效點(diǎn),根據(jù)復(fù)合梯形積分方法可得冷備份門的節(jié)點(diǎn)C在任務(wù)區(qū)間[0,T]的失效概率Pc1為

(27)

則Pc1是CSP輸出的失效概率真值,傳統(tǒng)BN轉(zhuǎn)化方法的計算結(jié)果誤差ε2為

(28)

2.2.3改良方法

根據(jù)傳統(tǒng)方法的計算誤差和復(fù)合梯形積分方法,利用得到的誤差直接修正節(jié)點(diǎn)B的CPT,改良備份門的轉(zhuǎn)化方法。在以前方法的基礎(chǔ)上,額外設(shè)定當(dāng)x=y時,B的條件概率為

Pxx=P(B=[(x-1)Δ,xΔ]|A=[(x-1)Δ,xΔ])=

(29)

特別地，對于服從指數(shù)分布的系統(tǒng)為

(30)

則改良后節(jié)點(diǎn)B在任務(wù)區(qū)間[0,T]的條件概率表分布式為

P(B=[(y-1)Δ,yΔ]|A=[(x-1)Δ,xΔ])=

(31)

對于備份門中不存在重復(fù)事件的情況下,改良后的轉(zhuǎn)化方法,在n值很小時能保證結(jié)果的正確性,即ε1=0。改良后的轉(zhuǎn)化方法比傳統(tǒng)方法擁有更高的計算精度。

3　實(shí)例分析

以某測量系統(tǒng)為例,驗證改良方法的可行性和高效性。測量系統(tǒng)的動態(tài)故障樹模型如圖3所示,其中M為密封膠,L、J是兩個關(guān)鍵電子元器件,D為紅外傳感器,E為接口,F為備用接口。

圖3　系統(tǒng)的動態(tài)故障樹

3.1　指數(shù)分布情形

不妨設(shè)定所有的部件在工作狀態(tài)下都服從指數(shù)分布。設(shè)其分布參數(shù)為:λA=2×10-3h-1、λB=5×10-3h-1、λC=8×10-3h-1、λD=3×10-3h-1、λE=10-3h-1和λF=4×10-3h-1。將動態(tài)故障樹轉(zhuǎn)化為BN,代入相關(guān)數(shù)據(jù)進(jìn)行分析計算。

使用Matlab軟件分別計算傳統(tǒng)DTBN方法和新的改良方法[17]。根據(jù)各個部件失效率的數(shù)量級,不妨設(shè)定工作時間T=1 000 h,兩種方法的結(jié)果對比如表4所示。

表4　兩種方法的失效概率與耗時

由表4可知，當(dāng)n取值相同時,兩種方法各個部件的CPT規(guī)模相同,改良方法的精度比傳統(tǒng)方法的精度更高,誤差更小。

對于所有部件服從指數(shù)分布的系統(tǒng),使用馬爾可夫方法計算得到系統(tǒng)的失效概率F(1 000)=0.557 203?？芍牧挤椒梢杂行Р⒑芸斓氐玫较到y(tǒng)的失效概率。在n很小時,改良方法就能得到系統(tǒng)失效概率的數(shù)量級。

蒙特卡羅仿真方法適用于服從任何分布類型的動態(tài)故障樹分析,但其往往需要耗費(fèi)大量的計算時間。運(yùn)用蒙特卡羅方法對系統(tǒng)進(jìn)行仿真的計算結(jié)果如表5所示。

表5　基于蒙特卡羅方法的系統(tǒng)失效概率與耗時

由表5可知，改良方法與蒙特卡羅方法相比,要得到同樣的精度,前者耗時更短。

當(dāng)n取值100時,運(yùn)用馬爾可夫方法、傳統(tǒng)方法和改良方法得到系統(tǒng)可靠度如圖4所示。

圖4　基于3種方法的系統(tǒng)可靠度變化曲線

由圖4可知，改良方法和馬爾可夫方法的曲線誤差區(qū)別極小,基本相互重疊,說明改良方法具有精度高、效率高等優(yōu)點(diǎn)。

3.2　非指數(shù)分布情形

當(dāng)部件L、E的失效分布服從兩參數(shù)威布爾分布時,其分布參數(shù)分別為:mB=0.95、ηB=3×102h;mE=1.05、ηE=2×102h。部件M、J、D和F服從指數(shù)分布,λA=2×10-3h-1、λC=8×10-3h-1、λD=3×10-3h-1和λF=4×10-3h-1。使用Matlab軟件,將動態(tài)故障樹轉(zhuǎn)化為BN,代入相關(guān)數(shù)據(jù)進(jìn)行分析計算。設(shè)定工作時間T=1 000 h,兩種方法的結(jié)果對比如表6所示。

表6　兩種方法的失效概率與耗時

對于部件服從非指數(shù)分布的情形,運(yùn)用蒙特卡羅方法對系統(tǒng)進(jìn)行仿真的計算結(jié)果如表7所示。

表7　基于蒙特卡羅方法的系統(tǒng)失效概率與耗時

由表7可知，改良方法與蒙特卡羅方法相比,要得到失效概率的數(shù)量級,前者耗時更短。

4　結(jié)　論

分析計算了DTBN方法的計算誤差,驗證了當(dāng)離散時間分段n的取值越大時,動態(tài)門的計算結(jié)果誤差越小,但同時計算時間變長,節(jié)點(diǎn)CPT的規(guī)模變大。

針對傳統(tǒng)方法的誤差,提出改良方法補(bǔ)償?shù)袅藗鹘y(tǒng)方法的計算誤差。改良方法擁有更高的計算精度和計算效率,解決了DTBN方法的計算時間、精度和效率受離散時間分段n極大影響的問題。在n很小時,改良方法就能得到系統(tǒng)失效概率的數(shù)量級。

通過某測量系統(tǒng)的實(shí)例分析,與馬爾可夫方法和蒙特卡羅仿真方法進(jìn)行對比,驗證了改良方法的可行性。

改良方法繼承了DTBN方法的優(yōu)點(diǎn),補(bǔ)償了傳統(tǒng)方法的誤差,得到了更高的計算精度和計算效率。

參考文獻(xiàn)：

[1] DUGAN J B, BAVUSO S J, BOYD M A. Dynamic fault-tree for fault-tolerant computer systems[J]. IEEE Trans.on Reliability, 1992, 41(3): 363-376.

[2] DUGAN J B, SULLIVAN K J, COPPIT D. Developing a low-cost high-quality software tool for dynamic fault-tree analysis[J]. IEEE Trans.on Reliability, 2000, 49(1): 49-59.

[3] SMOTHERMAN M K, ZEMIUDEH K. A non-homogenedous Makrov model for phased-mission reliability analysis[J]. IEEE Trans.on Reliability, 1989, 38(5): 585-590.

[4] DUGAN J B, Galileo: a tool for dynamic fault tree analysis[C]∥Proc.of the 11th International Conference on Computer Performance Evaluation: Modelling Techniques and Tools, 2000: 328-331.

[5] AMARI S, DILL G, HOWALD E. A new approach to solve dynamic fault trees[C]∥Proc.of the Annual IEEE Reliability and Maintainability Symposium, 2003: 374-379.

[6] AAMODT C A A, BANGS O, JENSEN F V, et al. Bayesian networks with applications in reliability analysis[J]. IEEE Trans.on Parallel & Distributed Systems, 2002, 22(3): 501-513.

[7] LANGSETH H, PORTINALE L, LANGSETH H,et al. Applications of Bayesian networks in reliability analysis[J].Bayesian Network Technologies Applications & Graphical Models,2007.

[8] SIMON C, WEBER P, EVSUKOFF A. Bayesian networks inference algorithm to implement dempster shafer theory in reliability analysis[J]. American Review of Respiratory Disease, 2008, 93(7): 950-963.

[9] BOBBIO A, PORTINALE L, MINICHINO M, et al. Improving the analysis of dependable systems by mapping fault trees into Bayesian networks[J].Reliability Engineering & System Safety,2001,71(3): 249-260.

[10] BOUDALI H, DUGAN J B. A new Bayesian network approach to solve dynamic fault trees[C]∥Proc.of the Reliability and Maintainability Symposium, 2005: 451-456.

[11] BOUDALI H, DUGAN J B. A discrete-time Bayesian network reliability modeling and analysis framework[J]. Reliability Engineering and System Safety, 2005, 87(3): 337-349.

[12] MONTANI S, PORTINALE L, BOBBIO A. Dynamic Bayesian networks for modeling advanced fault tree features in dependability analysis[C]∥Proc.of the ESREL, 2005:1414-1422.

[13] MONTANI S, PORTINALE L, BOBBIO A, et al. A tool for automatically translating dynamic fault trees into dynamic Bayesian networks[C]∥Proc.of the Reliability & Maintainability Symposium Annual, 2006: 434-441.

[14] NEIL M, TAILOR M, MARQUEZ D. Inference in hybrid Bayesian networks using dynamic discretization[J]. Statistics and Computing, 2007, 17(17): 219-233.

[15] MARQUEZ D, NEIL M, FENTON N. Improved reliability modeling using Bayesian networks and dynamic discretization[J]. Reliability Engineering and System Safety,2010,95(4):412-425.

[16] 周忠寶, 周經(jīng)倫, 孫權(quán), 等. 基于離散時間貝葉斯網(wǎng)絡(luò)的動態(tài)故障樹分析方法[J]. 西安交通大學(xué)學(xué)報, 2007, 41(6): 732-736.

ZHOU Z B, ZHOU J L, SUN Q, et al. Dynamic fault tree analysis method based on discrete-time Bayesian networks[J]. Journal of Xi’an Jiaotong University, 2007, 41(6): 732-736.

[17] MURPHY K P. The bayes net toolbox for Matlab[J]. Computing Science and Statistics,2001,33(2):1024-1034.