基于離散時間貝葉斯網絡的動態故障樹分析的改良方法

蘭　杰, 袁宏杰, 夏　靜

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

0　引　言

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

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

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

1　DTBN

DTBN是一個有向無環圖,主要包括兩部分:一部分是由有向邊和節點構成的網絡;另一部分是與每個節點相關的條件概率表(conditional probability table, CPT)。

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

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

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

2　動態故障樹轉化的誤差分析與補償

2.1　優先與門轉化分析

2.1.1傳統的優先與門轉化方法

優先與門有兩個輸入事件,當且僅當A失效,然后B失效,其輸出C才會失效。優先與門的輸入可以是基本事件,也可以是一個動態故障樹的子樹。優先門向DTBN的轉化圖如圖1所示。

圖1　優先與門轉化為DTBN

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

(1)

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

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

(2)

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

(3)

當A、B在時間段[T,∞)內失效,即A、B位于n+1狀態下的失效概率為

(4)

(5)

由優先與門性質可知,當且僅當A失效,然后B失效,其輸出C才會失效。傳統的貝葉斯轉化方法設定:如果節點A和節點B在同一個時間區間[(x-1)Δ,xΔ]內失效,即A、B同時處于狀態x時,節點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)的表。節點C的條件概率分布式為

(6)

當A、B中有任一個節點位于n+1狀態,即在時間區間[T,∞)內失效,C節點在時間區間[T,∞)內失效的條件概率為1。

對單一優先門,若只需求得節點C在任務時間[0,T]的失效概率和可靠度時,可簡化C的CPT的規模,設定節點A、B工作時間[0,∞)分為n+1段,節點C的時間區間分為[0,T]、[T,∞)兩段。此時節點C的CPT為2×(n+1)×(n+1)的表,規模大幅度變小,且計算結果與第3.1節方法的結果一樣。節點C在任務區間[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)

節點C在時間區間[T,∞)失效,則條件概率表分布式為

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

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

(9)

節點C在任務區間[0,T]內失效的條件概率表為如表1所示,同理也可得節點C在區間[T,∞)的CPT。

表1　優先與門節點C的CPT

2.1.2誤差分析

在傳統的BN轉化方法下,優先門輸出點C在任務區間[0,T]的失效概率Pc為

(10)

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

(11)

設ta為A的失效點,tb為B的失效點,根據復合梯形積分方法可得優先與門的輸出頂事件C在任務區間[0,T]的失效概率Pc1為

(12)

Pc1是優先與門的頂事件的失效概率的真值,則傳統的BN轉化方法的計算結果的誤差ε1為ε1=|Pc-Pc1|。

其中可用將任務區間[0,T]分為n段的方法,對式(12)進行轉化,則

(13)

由式(11)和式(13)得

(14)

利用夾逼定理可得

(15)

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

2.1.3改良方法

根據前面得出的傳統方法的計算誤差和復合梯形積分方法,利用得到的誤差直接修正節點C的CPT,改良優先與門的轉化方法。改良方法在動態門中不存在重復事件的情況下完全補償掉了計算誤差,進而減小了存在重復事件情形下的動態門誤差。

設定如果節點A和節點B在同一個時間區間[(x-1)Δ,xΔ]內失效時,節點C的失效概率不直接為0,即

PT,x=

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

(16)

特別對于服從指數分布的系統,即

(17)

則改良后節點C在任務區間[0,T]失效,條件概率表分布式為

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

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

(18)

改良節點C在任務區間[0,T]內失效的條件概率表如表2所示。

表2　改良后的優先與門節點C的CPT

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

2.2　備份門轉化分析

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

2.2.1傳統的備份門轉化方法

由CSP性質可知,當且僅當主件失效后,備份單元才進入工作狀態,備份單元在儲備狀態時失效率為零。CSP向DTBN轉化如圖2所示。

圖2　冷備份轉化為DTBN

當x

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

(19)

對于服從指數分布的系統,不妨設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)

節點B的CPT如表3所示。

表3　備份門節點B的CPT

節點C的條件概率分布式為

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

(24)

2.2.2誤差分析

在傳統轉化方法下,優先門輸出點C在任務區間[0,T]的失效概率Pc為

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

(25)

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

(26)

設ta為A的失效點,tb為B的失效點,根據復合梯形積分方法可得冷備份門的節點C在任務區間[0,T]的失效概率Pc1為

(27)

則Pc1是CSP輸出的失效概率真值,傳統BN轉化方法的計算結果誤差ε2為

(28)

2.2.3改良方法

根據傳統方法的計算誤差和復合梯形積分方法,利用得到的誤差直接修正節點B的CPT,改良備份門的轉化方法。在以前方法的基礎上,額外設定當x=y時,B的條件概率為

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

(29)

特別地，對于服從指數分布的系統為

(30)

則改良后節點B在任務區間[0,T]的條件概率表分布式為

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

(31)

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

3　實例分析

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

圖3　系統的動態故障樹

3.1　指數分布情形

不妨設定所有的部件在工作狀態下都服從指數分布。設其分布參數為:λ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。將動態故障樹轉化為BN,代入相關數據進行分析計算。

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

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

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

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

蒙特卡羅仿真方法適用于服從任何分布類型的動態故障樹分析,但其往往需要耗費大量的計算時間。運用蒙特卡羅方法對系統進行仿真的計算結果如表5所示。

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

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

當n取值100時,運用馬爾可夫方法、傳統方法和改良方法得到系統可靠度如圖4所示。

圖4　基于3種方法的系統可靠度變化曲線

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

3.2　非指數分布情形

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

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

對于部件服從非指數分布的情形,運用蒙特卡羅方法對系統進行仿真的計算結果如表7所示。

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

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

4　結　論

分析計算了DTBN方法的計算誤差,驗證了當離散時間分段n的取值越大時,動態門的計算結果誤差越小,但同時計算時間變長,節點CPT的規模變大。

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

通過某測量系統的實例分析,與馬爾可夫方法和蒙特卡羅仿真方法進行對比,驗證了改良方法的可行性。

改良方法繼承了DTBN方法的優點,補償了傳統方法的誤差,得到了更高的計算精度和計算效率。

參考文獻：

[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]. 西安交通大學學報, 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.