基于自適應(yīng)MCMC和子集模擬的非能動(dòng)系統(tǒng)熱工水力可靠性評(píng)估

蔣立志，蔡 琦，張永發(fā)

(海軍工程大學(xué) 核科學(xué)技術(shù)學(xué)院，湖北 武漢 430033)

與能動(dòng)系統(tǒng)相比，非能動(dòng)系統(tǒng)驅(qū)動(dòng)力較小，設(shè)計(jì)參數(shù)、初始條件等因素波動(dòng)帶來(lái)的不確定性以及對(duì)相關(guān)物理現(xiàn)象認(rèn)識(shí)的不充分，均可能造成熱工水力過(guò)程失效、無(wú)法實(shí)現(xiàn)設(shè)計(jì)功能。對(duì)非能動(dòng)系統(tǒng)進(jìn)行熱工水力可靠性(TH-R)評(píng)估能從可靠性角度增強(qiáng)對(duì)非能動(dòng)系統(tǒng)熱工水力過(guò)程不確定性的認(rèn)識(shí)，相關(guān)分析結(jié)果可作為系統(tǒng)性能評(píng)估與優(yōu)化的參考。最初非能動(dòng)系統(tǒng)TH-R評(píng)估主要采用以蒙特卡羅(MC)為基礎(chǔ)的模擬方法[1-2]，需大量調(diào)用系統(tǒng)熱工水力模型進(jìn)行計(jì)算，以獲取足夠的輸入、輸出數(shù)據(jù)，但以RELAP5為代表的熱工水力程序運(yùn)行耗時(shí)較長(zhǎng)，帶來(lái)極高的運(yùn)算代價(jià)和較低的評(píng)估效率。重要抽樣(IS)方法的引入[3-4]一定程度緩解了運(yùn)算代價(jià)的問(wèn)題，但在高維問(wèn)題中重要抽樣密度函數(shù)的構(gòu)造十分困難且存在維數(shù)災(zāi)難問(wèn)題[5]。子集模擬(SS)方法[6]適用于高維小失效率模型，它的引入[7-8]在評(píng)估效率和穩(wěn)健性之間帶來(lái)了一種均衡，其基本思想是將較小失效概率分解為較大失效概率的乘積并依次估計(jì)這些條件失效概率。SS中條件失效概率的估計(jì)需對(duì)各中間失效事件的條件概率密度函數(shù)進(jìn)行精確采樣，通常采用修正Metropolis-Hastings(MMH)[9]算法實(shí)現(xiàn)。MMH算法中建議分布的選擇對(duì)算法的效率和失效率估計(jì)精度影響很大，最優(yōu)建議分布的選擇較為困難，且隨失效率水平的降低算法性能會(huì)出現(xiàn)退化。文獻(xiàn)[10-11]討論了Metropolis-Hastings(MH)算法建議分布的最優(yōu)問(wèn)題和自適應(yīng)馬爾可夫鏈蒙特卡羅(MCMC)的合理框架，但相關(guān)優(yōu)化結(jié)論針對(duì)特定目標(biāo)分布假設(shè)，并不完全適用于子集模擬的中間條件分布。文獻(xiàn)[12]基于標(biāo)準(zhǔn)正態(tài)空間內(nèi)的條件采樣(CS)算法提出一種適用于SS的自適應(yīng)MCMC方法，稱為自適應(yīng)條件采樣(aCS)方法。該方法在當(dāng)前鏈狀態(tài)和候選樣本間引入相關(guān)性，避免了建議分布選擇的復(fù)雜工作[13]，能在條件采樣過(guò)程中直接和靈活地調(diào)整算法參數(shù)，使候選樣本接受率在目標(biāo)值附近保持穩(wěn)定，以獲取更高的算法效率和失效率估計(jì)精度。

本文引入aCS-SS方法，提出一種基于自適應(yīng)MCMC和SS的非能動(dòng)系統(tǒng)TH-R評(píng)估方法，以某型核動(dòng)力裝置二次側(cè)非能動(dòng)余熱排出試驗(yàn)系統(tǒng)為例進(jìn)行驗(yàn)證。

1 基本理論

1.1 SS的基本思想

假設(shè)系統(tǒng)輸入變量的聯(lián)合概率密度函數(shù)為fX(x)，則系統(tǒng)的失效率PF為：

(1)

式中：x=(x1,x2,…,xn)∈Rn為輸入?yún)?shù)向量，Rn為n維輸入?yún)?shù)空間；F為失效域；IF(x)為示性函數(shù)，x∈F時(shí)IF(x)=1，否則IF(x)=0。采用MC方法時(shí)系統(tǒng)失效率的估計(jì)值為：

(2)

(3)

式(3)中條件失效概率的估計(jì)值為：

(4)

式中，fX(x|Fj-1)=fX(x)IFj-1(x)/P(Fj-1)為x落入Fj-1內(nèi)的條件概率密度函數(shù)。由于F0是一個(gè)必然事件，P1的計(jì)算可采用標(biāo)準(zhǔn)的MC模擬進(jìn)行。當(dāng)j≥2時(shí)Pj的估計(jì)需根據(jù)條件分布fX(x|Fj-1)進(jìn)行采樣，高效的條件采樣是SS的核心。

1.2 基于MMH算法的SS

3) 重復(fù)上述步驟，直至樣本數(shù)量滿足要求。

1) 設(shè)定各層的樣本數(shù)N及條件失效概率p0，應(yīng)該盡量保證p0N為正整數(shù)。

1.3 CS的基本原理

上述基于CS算法的MCMC策略其穩(wěn)態(tài)分布漸進(jìn)收斂于條件分布fX(x|Fj-1)，文獻(xiàn)[12]給出了詳細(xì)的證明過(guò)程。CS算法的效率和性能僅依賴于相關(guān)系數(shù)ρi的選擇，與MMH算法相比避免了建議分布形式、參數(shù)的選擇，避免了預(yù)候選樣本的接受/拒絕過(guò)程，算法形式更加簡(jiǎn)潔、靈活，運(yùn)算代價(jià)更低。但針對(duì)不同的系統(tǒng)類型和失效率水平，如何合理選擇相關(guān)系數(shù)ρi依然是一個(gè)困難的問(wèn)題，較高和較低的相關(guān)系數(shù)均會(huì)導(dǎo)致條件樣本的相關(guān)性過(guò)高，降低算法效率和估計(jì)精度。

1.4 基于aCS的子集模擬

文獻(xiàn)[10-11]的相關(guān)研究表明，目標(biāo)分布為一維正態(tài)分布時(shí)MH算法存在一近似的最優(yōu)接受率a*≈0.44，文獻(xiàn)[15]的算例表明SS中MMH算法的最優(yōu)接受率在0.3～0.5之間。CS算法中候選樣本的接受/拒絕僅依賴于系統(tǒng)極限狀態(tài)函數(shù)的響應(yīng)，因此從條件分布fX(x|Fj-1)采樣可等效為一維條件采樣問(wèn)題，文獻(xiàn)[12]提出一種aCS方法，能使候選樣本接受率較穩(wěn)定地保持在目標(biāo)值附近。利用aCS算法從條件分布fX(x|Fj-1)采樣的詳細(xì)流程如下。

2) 自適應(yīng)迭代過(guò)程：令迭代次數(shù)r=1,2,…,Ns/Na，重復(fù)執(zhí)行以下步驟。

(1) 根據(jù)標(biāo)準(zhǔn)差初值σ0i和當(dāng)前尺度參數(shù)λr計(jì)算相關(guān)系數(shù)：

σi=min(λrσ0i,1)

(5)

(6)

(7)

(8)

2 算例分析及結(jié)果討論

本文以某型核動(dòng)力裝置二次側(cè)非能動(dòng)余熱排出試驗(yàn)系統(tǒng)[16-18]作為算例驗(yàn)證aCS-SS算法的有效性，系統(tǒng)原理如圖1所示。采用文獻(xiàn)[16]建立的熱工水力模型并選擇啟動(dòng)工況5為運(yùn)行工況進(jìn)行分析，該模型使用RELAP5/Mod3.2作為建模工具。由于缺少充足的試驗(yàn)數(shù)據(jù)和運(yùn)行經(jīng)驗(yàn)，通過(guò)專家分析及工程判斷方法識(shí)別影響系統(tǒng)熱工水力行為的32個(gè)不確定性輸入?yún)?shù)，參數(shù)分布類型列于表1。由于試驗(yàn)系統(tǒng)未模擬反應(yīng)堆一回路，因此選擇關(guān)注的系統(tǒng)輸出為系統(tǒng)運(yùn)行期間(20 000 s)冷凝器管側(cè)壓力峰值pMAX：若運(yùn)行期間pMAX不大于管道所能承受的壓力上限值pm，則認(rèn)為蒸汽發(fā)生器模擬體中的熱量能被成功排出至冷卻水箱，即g(x)=pm-pMAX<0時(shí)系統(tǒng)功能失效，否則系統(tǒng)功能正常。

圖1 非能動(dòng)余熱排出系統(tǒng)原理簡(jiǎn)圖Fig.1 Simplified figure of passive residual heat removal system

表1 不確定性參數(shù)分布類型Table 1 Distribution type of uncertainty parameter

注：N表示正態(tài)分布，U表示均勻分布

為比較不同失效率水平下的算法性能，設(shè)置4組壓力上限值pm1

1) MMH-SS：采用MMH算法產(chǎn)生條件樣本，建議分布選取為正態(tài)分布，其標(biāo)準(zhǔn)差根據(jù)種子樣本計(jì)算得到。

2) CS-SS：采用CS算法產(chǎn)生條件樣本，為分析不同相關(guān)系數(shù)對(duì)失效率估計(jì)精度的影響，各分量相關(guān)系數(shù)取3種配置：ρi=0.7、ρi=0.8、ρi=0.9，分別記為CS-SS-a、CS-SS-b和CS-SS-c。

3) aCS-SS-a：采用aCS算法產(chǎn)生條件樣本，標(biāo)準(zhǔn)差初值取σ0i=1，尺度參數(shù)初值取λ1=0.8。

4) aCS-SS-b：采用aCS算法產(chǎn)生條件樣本，標(biāo)準(zhǔn)差初值根據(jù)種子樣本計(jì)算得到，尺度參數(shù)初值取λ1=0.8。

由于a*≈0.44是在目標(biāo)分布為一維正態(tài)分布的假設(shè)下獲得的結(jié)論，并不完全適用于本文所評(píng)估的系統(tǒng)。經(jīng)過(guò)測(cè)試，系統(tǒng)失效率為1.94×10-3時(shí)MMH-SS算法的候選樣本接受率約為0.5，且失效率估計(jì)精度較高，因此將aCS-SS的目標(biāo)接受率設(shè)定為a*=0.5。設(shè)定N=1 000、p0=0.1，在不同失效率水平下，利用4種策略各進(jìn)行100次失效率估計(jì)，各策略下采用相同的100個(gè)隨機(jī)數(shù)種子。最終得到整個(gè)SS過(guò)程中候選樣本的平均接受率隨失效率水平的變化，如圖2所示。失效率估計(jì)平均誤差隨失效率水平的變化如圖3所示，失效率估計(jì)平均變異系數(shù)(SS變異系數(shù)的近似估計(jì)方法參考文獻(xiàn)[7,13])隨失效率水平的變化如圖4所示。

圖2 平均接受率隨失效率水平的變化Fig.2 Average acceptance rate plotted against varying failure probability level

圖3 失效率估計(jì)平均誤差隨失效率水平的變化Fig.3 Average error of failure probability estimate plotted against varying failure probability level

圖4 失效率估計(jì)平均變異系數(shù)隨失效率水平的變化Fig.4 Average coefficient of variation of failure probability estimate plotted against varying failure probability level

由圖2可知：隨系統(tǒng)失效率水平的不斷降低，4種估計(jì)方法的平均接受率均有所下降，但aCS-SS-a和aCS-SS-b能使接受率在目標(biāo)值附近保持穩(wěn)定水平；CS-SS的平均接受率與相關(guān)系數(shù)的取值有直接關(guān)系，相關(guān)系數(shù)越大平均接受率越高，反之亦然。結(jié)合圖3、4的結(jié)果可知：隨失效率水平的降低，4種方法的平均估計(jì)精度均降低、平均變異系數(shù)均增大，這是由于較低失效率水平下失效區(qū)域非常狹小，有效的條件采樣變得更加困難；失效率處于較高水平時(shí)4種方法精度差別不大，但失效率處于較低水平時(shí)aCS-SS-b的平均精度最高、平均變異系數(shù)最?。籆S-SS算法中相關(guān)系數(shù)ρi取不同值時(shí)失效率平均估計(jì)精度和平均變異系數(shù)變化劇烈，對(duì)于本文算例ρi=0.8時(shí)失效率平均估計(jì)精度較高、變異系數(shù)較大，但在實(shí)際工程應(yīng)用中對(duì)于不同的研究對(duì)象很難準(zhǔn)確地找到一個(gè)比較合理的相關(guān)系數(shù)取值；與aCS-SS-a相比，aCS-SS-b的平均精度更高、變異系數(shù)更小，這是由于aCS-SS-b中標(biāo)準(zhǔn)差初值σ0i由種子樣本估計(jì)得到，這一過(guò)程本質(zhì)上是對(duì)失效區(qū)域的一種自適應(yīng)探索，當(dāng)各輸入?yún)?shù)對(duì)系統(tǒng)性能的影響不相同時(shí)aCS-SS-b的性能必然高于aCS-SS-a，但當(dāng)各輸入?yún)?shù)對(duì)系統(tǒng)性能的影響基本相同時(shí)二者的性能無(wú)明顯差別。

SS的失效率估計(jì)精度取決于各子集邊界的準(zhǔn)確劃分和條件失效域的高效采樣。采樣過(guò)程中候選樣本接受率過(guò)低或過(guò)高，均會(huì)導(dǎo)致條件失效樣本間相關(guān)性較高、重復(fù)樣本較多，最終導(dǎo)致各子集邊界劃分出現(xiàn)偏差以及最后一層子集的失效率估計(jì)誤差增大。在失效率為2.90×10-6的水平下，分別利用MMH-SS、CS-SS-b、aCS-SS-a和aCS-SS-b進(jìn)行1次失效率估計(jì)，第2～6層子集對(duì)應(yīng)的候選樣本接受率變化列于表2。由表2的變化規(guī)律可看出：隨各子集的逐漸劃分，候選樣本接受率逐漸降低，這是由于后續(xù)子集對(duì)應(yīng)的區(qū)域越來(lái)越狹小、采樣也越來(lái)越困難；MMH-SS和CS-SS算法中各子集接受率的降低較為劇烈，這是因?yàn)楦髯蛹臈l件采樣使用的均是相同的策略或相關(guān)系數(shù)；aCS-SS中接受率的下降較為緩慢，即使在最后一層子集也能保持接近0.5的水平，這是由于aCS-SS的自適應(yīng)調(diào)節(jié)過(guò)程可令各子集對(duì)應(yīng)的相關(guān)系數(shù)逐漸增大，特別是aCS-SS-b中對(duì)失效區(qū)域的自適應(yīng)探索過(guò)程還可使各輸入?yún)?shù)對(duì)應(yīng)的相關(guān)系數(shù)隨參數(shù)重要度的不同自動(dòng)調(diào)節(jié)。

表2 各子集候選樣本的接受率Table 2 Acceptance rate of candidate sample for each subset

3 小結(jié)

本文引入一種基于aCS的自適應(yīng)MCMC算法，提出一種基于自適應(yīng)MCMC和SS的非能動(dòng)系統(tǒng)熱工水力可靠性評(píng)估方法，以某型核動(dòng)力裝置二次側(cè)非能動(dòng)余熱排出試驗(yàn)系統(tǒng)為算例進(jìn)行驗(yàn)證，對(duì)4種估計(jì)方法的性能進(jìn)行了比較。計(jì)算結(jié)果表明：aCS-SS避免了MMH-SS中建議分布形式和參數(shù)的選擇，算法形式更加簡(jiǎn)潔、靈活，算法中的參數(shù)可自適應(yīng)調(diào)整以保持穩(wěn)健性；aCS-SS能更好地使候選樣本接受率在目標(biāo)值附近保持穩(wěn)定水平，在較低失效率水平(10-5、10-6)下aCS-SS的估計(jì)精度和穩(wěn)健性均高于目前常用的MMH-SS算法。隨非能動(dòng)技術(shù)的不斷發(fā)展、運(yùn)行經(jīng)驗(yàn)的不斷積累以及各類非能動(dòng)安全系統(tǒng)的不斷改進(jìn)，反應(yīng)堆非能動(dòng)系統(tǒng)的物理過(guò)程失效會(huì)進(jìn)一步降低，本文方法對(duì)于較低失效率水平下非能動(dòng)系統(tǒng)TH-R評(píng)估問(wèn)題是非常有前景的。