基于AK-SDMCS 的含缺陷核管道可靠性分析

薛國峰，汪大洋，張永山，蘇邁佳

1. 佛山科學技術學院 交通與土木建筑學院，廣東 佛山 528225

2. 廣州大學 土木工程學院，廣東 廣州 510006

3. 特倫托大學 土木、環境與機械工程學院，意大利 特倫托 38100

管道系統是核電站冷卻劑壓力邊界的重要組成部分，是核電站安全和正常運行的重要保障[1-3]。由于制造工藝及焊接等因素影響，管道系統不可避免地存在初始缺陷[4-5]。在運行荷載及溫度、疲勞、腐蝕、老化等因素作用下，管道中可能萌生新的缺陷，已存在的初始缺陷則可能進一步擴展，直至發生泄漏、管道斷裂等事故[6]。因此，非常有必要開展核電站含缺陷管道的安全評定[7]。

現行規范主要采用確定性方法進行管道安全評定。如英國R6 評定方法將含缺陷管道分為脆性斷裂、塑性失穩和彈塑性斷裂共3 種失效模式，根據管道及缺陷的幾何尺寸、材料、荷載等因素判斷失效模式，進而進行安全評定[8-9]。美國機械工程師學會（ American Society of Mechanical Engineers, ASME）管道與壓力容器規范中的評定方法以R6 為基礎，分別給出基于失效模式判定的含缺陷鐵素體和奧氏體管道安全評定方法以及基于失效評定圖的含缺陷管道安全評定方法[10]。但應注意到，管道系統的材料、幾何形狀、構造、工作環境及承受的荷載等均具有很強的隨機性，缺陷形狀、尺寸、位置和檢測手段等也具有一定的隨機性。在這些復雜因素綜合作用下，核電站管道系統的工作狀態及失效均表現出較強的隨機性。因此，有必要采用概率統計方法合理考慮這些隨機性，準確計算含缺陷管道的失效概率，從而準確評價管道的安全性和可靠性，為管道系統設計、優化、檢測及維護提供理論基礎[11-12]。

現有的含缺陷核管道可靠性分析方法主要包括一次可靠性方法(first order reliability method,FORM)、二次可靠性方法(second order reliability method, SORM)、 蒙特卡羅模擬(Monte Carlo simulation, MCS)及重要抽樣法(importance sampling,IS)等[9,13-20]。FORM 和SORM 計算效率高，但對非線性功能函數計算誤差較大[21]。數值模擬法(MCS、IS 等)計算精度高，但由于需要生成大量模擬樣本并計算結構響應而計算量過大[22-23]。近年來，主動學習Kriging 代理模型可靠性分析方法由于同時結合了代理模型法效率高和數值模擬法精度高的優點而受到較多關注[24-27]。其基本思想是基于數值模擬法的樣本池，采用主動學習方法以少量樣本來訓練代理模型，進而采用此代理模型代替真實結構功能函數來計算失效概率，從而實現了以少量結構功能函數調用（即少量結構分析次數）獲得高精度失效概率的目的，提高了可靠性分析的計算效率和精度。經典的主動學習方法主要為基于主動學習 Kriging 代理模型的Monte Carlo 模擬方法(an active learning reliability method combining Kriging and Monte Carlo Simulation, AK-MCS)[25]，其主要基于MCS 樣本池展開。但由于MCS 樣本池容量通常很大，在其上進行主動學習時需要多次調用Kriging 模型計算響應。尤其是對于小失效概率問題，主動學習本身所需的計算量就非常大，阻礙了主動學習方法的應用。為了降低主動學習所需的計算量，人們提出基于方差縮減技術的主動學習可靠性分析方法如基于主動學習 Kriging 代理模型的重要抽樣法 (a combined Importance sampling and active Kriging reliability method, AK-IS)和基于主動學習Kriging 代理模型的子集模擬方法(an active learning method combining Kriging and subset simulation, AK-SS)等[28-30]。但這些方法仍存在一些共性問題，主要包括：1)無法保持AK-MCS 的通用性；2)需要前處理獲取失效域信息；3)抽樣區域大于重要區域，使得在對失效概率無較大貢獻的區域內浪費了部分功能函數調用。因此，本文提出了基于參數空間球分解蒙特卡羅模擬的主動學習代理模型可靠性分析方法（an active learning reliability method combining adaptive Kriging and spherical decomposition-MCS, AK-SDMCS）[31]，其通過將參數空間沿著徑向由內向外劃分為一系列超球環，從而將整個空間內進行的主動學習分解為一系列超球環內的主動學習，通過“化整為零、分而治之”，提高了主動學習的計算效率。AKSDMCS 既保持了AK-MCS 的通用性和魯棒性，同時克服了在求解小失效概率問題時主動學習本身計算量過大的困難。

因此，本文圍繞含缺陷核管道可靠性分析開展研究。首先，基于ASME 鍋爐與壓力容器規范建立含環向內表面缺陷核管道的可靠性分析模型，包括塑性失穩、彈塑性撕裂失穩和脆性斷裂3 種主要失效模式及對應的極限狀態方程，歸納凝練了含缺陷核管道可靠性分析應考慮的隨機因素：運行壓力、缺陷尺寸、材料斷裂韌性和流變應力等；其次，基于參數空間徑向球分解思想，將整個參數空間內生成的樣本池上進行的主動學習（AK-MCS）分解為劃分所得的一系列超球環內樣本池上的主動學習，降低主動學習的計算量，克服AK-MCS 對小失效概率計算量過大的缺點，建立AK-SDMCS 方法；最后，采用AK-SDMCS 計算含缺陷核管道的失效概率，建立高效適用的含缺陷核管道可靠性分析方法。本文還給出3 個不同數量級（10-5~10-3）的算例，通過與MCS、AKMCS 對比，驗證了所提AK-SDMCS 進行含缺陷核管道可靠性分析的計算精度和效率。

1 含缺陷核管道可靠性分析模型

1.1 含缺陷管道截面參數

本文主要研究含環向內表面缺陷核管道的可靠性分析，截面如圖1 所示。

圖1 含缺陷核管道截面示意

圖1 中t為管道壁厚，mm；R1、R2、Rm分別為管道內半徑、外半徑和平均半徑，mm；缺陷通常簡化為弧形，a為 缺陷深度，mm； θ為缺陷半長度角(弧度制)；c=θ·R1為缺陷半長，mm；l=2c為缺陷長度，mm； β為中性角(弧度制)。

1.2 含缺陷核管道失效模式

含缺陷管道的失效模式通常分為塑性失穩、彈塑性撕裂失穩和脆性斷裂3 類，主要通過綜合考慮了荷載、缺陷尺寸、材料性質等因素的篩選參數SC進行判別[32]，如圖2 所示。

圖2 判斷含缺陷核管道的失效模式

1) 當SC< 0.2 時，判定為塑性失穩失效模式，需要進行極限塑性荷載分析(limit load (fully plastic) analysis)。

2) 當0.2≤SC<1.8 時，判定為彈塑性撕裂失穩失效模式，需要進行彈塑性斷裂力學分析(elasticplastic fracture mechanics analysis, EPFM)。

3) 當SC≥1.8 時，判定為脆性斷裂失效模式，應采用線彈性斷裂力學方法分析(linear elastic fracture mechanics analysis, LEFM)，一般以應力強度因子為判據。

篩選參數SC按下式計算：

應力強度因子KI為

式中Fm和Fb為計算參數，可參考ASME 規范計算，二者均為管道內徑和缺陷尺寸的函數；Q為缺陷形狀無量綱參數。

式中： σy為材料屈服強度； σf為材料流變應力，一般定義為屈服強度和拉伸強度的平均值。中性角β計算為

其中

塑性失穩時的極限狀態方程（對應SC< 0.2）為

式中SC1為塑性失穩條件下的許用彎曲應力。

式中：SFb和SFm為結構系數，一般根據管道所處的不同工況取不同值；為失效彎曲應力。

塑性失穩時的極限狀態方程（對應0.2 ≤SC<1.8）為

式中SC2為塑性失穩條件下的許用彎曲應力。

式中Z是根據管道尺寸計算得到的乘子，可參考ASME 規范計算。

脆性斷裂失效時的極限狀態方程（對應SC≥1.8）為

式中Kc=(JIcE′/1 000)0.5為材料的臨界斷裂韌性。應注意的是，這里的應力強度因子計算公式不同于計算SC時采用的應力強度因子，主要區別在于取消了缺陷形狀無量綱參數Q，并增加了結構系數SFb和SFm。應力強度因子具體計算公式為

其中

式中Fm和Fb的計算公式同前文，均可參考ASME 規范計算。

1.3 隨機參數

含缺陷核管道可靠性分析考慮的隨機因素主要為運行壓力、缺陷尺寸、材料斷裂韌性和流變應力等[13]。

1.3.1 運行壓力的隨機性

管道所受的流體內壓p通常假定為正態分布隨機變量。根據內壓計算薄膜應力 σm和彎曲應力σb，通過對數據統計則可獲得相應的概率分布。應注意到 σm和 σb均依賴于內壓p，顯然二者并非相互獨立隨機變量。實際處理中，可通過有限元分析等手段建立內壓與薄膜應力和彎曲應力的確定性函數關系，如σm=f(p)和σb=h(p)，進而由內壓的概率分布推導出薄膜應力和彎曲應力的概率分布。在隨機模擬場合，只需生成內壓的隨機樣本，再把這些隨機樣本代入確定性函數f(p)和h(p)，則可獲得對應的薄膜應力和彎曲應力樣本[15]。

1.3.2 缺陷尺寸的隨機性

為了便于計算，環向內表面缺陷通常簡化為弧形。由檢測精度和表征誤差引起的缺陷深度和長度的隨機特征一般均可表示為正態分布，相關分布參數可通過對大量實測數據統計獲得[33-34]。

1.3.3 流變應力和斷裂韌性的隨機性

由于制造工藝、檢測手段等多種復雜因素影響，材料性質通常并非常數。即使是同一批次生產，材料性質仍表現出較大的隨機性。對于含缺陷核管道可靠性分析，材料的流變應力和斷裂韌性為主要影響參數，對應的概率分布可通過對實測數據統計分析獲得，通常假定為對數正態分布[33]。

1.4 含缺陷管道的極限狀態函數

綜上，含環向內表面缺陷核管道的極限狀態函數可表示為三分段函數，如下式所示：

式中x為基本隨機變量，這里主要包括薄膜應力σm、彎曲應力 σb、缺陷深度a、缺陷長度c(或缺陷半長度角θ)、流變應力 σf和斷裂韌性JIC等。基于上述極限狀態函數，含缺陷管道的失效概率可表示為

式中fX(x)為基本隨機變量的聯合概率密度。

2 含缺陷核管道可靠性分析方法

2.1 Kriging 模型

Kriging 模型由Krige 和Matheron 等首先提出并應用于地質統計學[35]，后應用于計算機模擬試驗[36]、全局優化[37]等領域。在結構可靠性分析領域，Kriging 模型主要用作一種代理模型[38-39]。此外，很多主動學習可靠性分析方法均基于Kriging 模型建立[24,40-41]。

Kriging 模型為基于統計理論的半參數插值模型，由參數化的線性回歸模型和非參數化的隨機過程2 部分組成[37]。線性回歸模型表示原函數的整體趨勢，隨機過程則用于修正回歸模型的偏差。Kriging 模型的應用主要包括2 步：首先基于給定的實驗設計確定模型參數；其次通過最佳線性無偏預測計算結構在新輸入樣本下的響應。給定實驗設計：

式中：fT(x)β和f(x)=[f1(x),f2(x),···,fm(x)]T分別為回歸模型和基函數；β=[β1,β2,···,βm]T為對應的回歸系數。本文采用常數回歸模型，即取f(x)=1，β=β，從而將Kriging 模型簡化為

式中z(x)為零均值平穩高斯過程，其協方差為

式中： σ2為隨機過程的方差；R(x,w)為任意兩點x和w的自相關函數，通常采用高斯自相關函數。

式中：xi和wi分別為點x和w的第i個分量， θi為第i維的相關參數。

令R為p×p維對稱相關矩陣， 其中Rij=R(x(i),x(j)),i,j=1,2,···,p，F為p×1維單位向量，則回歸系數 β和隨機過程的方差 σ2可按以下式進行估計：

相關參數 θ可采用極大似然法按下式估計：

對于新輸入樣本x(0)，Kriging 模型可給出結構響應的最佳線性無偏預測和相應的Kriging方差，分別為

式中

容易驗證，對實驗設計中每一個樣本x(i)，均有即Kriging 模型為精確插值模型。對于實驗設計之外的輸入樣本，Kriging 方差一般不為0，可表示模型預測值的局部不確定性。因此，Kriging 方差常被用于改進實驗設計和提高代理模型精度。本文采用MATLAB工具箱DACE 進行Kriging 模型的構造與預測[42]。

2.2 AK-SDMCS 簡介

AK-SDMCS 的基本思想是首先將參數空間沿著徑向由內向外劃分為一系列超球環（即參數空間球分解）；其次，在這些超球環內分別采用主動學習方法構造Kriging 模型；最后，采用最終獲得的Kriging 模型計算結構響應和失效概率。

2.2.1 SDMCS

首先，將參數空間D=沿徑向劃分為一系列互不重疊且完備的子區域Di,i=1,2,···,m，使得：

式中：Di={x|Ri-1≤∥x∥

定義子區域Di上的截斷概率密度為

將式(5)代入式(4)可得失效概率估計：

對應的估計方差為

參數空間的球分解也采用自適應方式進行。首先將空間沿著徑向劃分為2 個超球環，然后將外側超球環進一步劃分為2 個超球環，如此重復迭代直到最外側超球環中所需的樣本量小于事先規定的樣本量，即可自動確定半徑序列R0

由于標準正態空間內樣本到原點的距離服從自由度為N的卡方分布(|x|2)， 因此有=(1-10-j)。當最外側超球環(子區域)內所需的樣本量滿足≤?N0時，表明個樣本即可充分探索該區域（注：=NMCSθi，NMCS可由事先設定的變異系數，結合MCS 計算公式獲得。實際應用中，失效概率通常未知，可采用當前失效概率估計值替代）。此時，可停止空間分解，相應的超球環數量即為m。

2.2.2 AK-SDMCS

在SDMCS 基礎上，在每個超球環內分別采用主動學習方法構造Kriging 代理模型，并采用最終獲得的Kriging 模型計算失效概率，即可建立AK-SDMCS。

AK-SDMCS 的初始實驗設計采用拉丁超立方抽樣(latin hypercube sampling, LHS)生成。模擬區域選為[u±VMσ]，其中u和 σ分別為基本隨機變量的均值和標準差，VM為標準差縮放系數。初始實驗設計的樣本容量為max{15,N+2}。由于LHS 具有更好的空間填充特性，因此，訓練所得的初始Kriging 模型將較好地刻畫真實功能函數的全局特征。

其次，考慮到經典的AK-MCS 采用的停止條件過于保守[43]，本文提出由AK-MCS 的停止條件和文獻[44]中的停止條件組成的復合停止條件。

其中第1 個停止條件表示失效域的穩定性。

式中為第j輪主動學習迭代中的失效樣本數量。第2 個停止條件表示相鄰2 次主動學習中候選樣本響應符號發生變化的比例。

式中：NS為候選樣本總量；為第j輪迭代訓練所得的代理模型； ×表示異或運算；ε1,ε2為2 個閾值，均取0.001[44]。

2.3 基于AK-SDMCS 的含缺陷核管道可靠性分析方法

基于前文建立的含缺陷核管道可靠性分析模型，根據基本隨機變量概率分布生成樣本，采用AK-SDMCS 計算失效概率，即可建立基于AKSDMCS 的含缺陷核管道可靠性分析方法，基本流程如圖3 所示。

圖3 含缺陷核管道可靠性分析基本流程

方法基本步驟如下：

1) 初始化：主要包括確定管道參數和基本隨機變量的概率分布以及AK-SDMCS 基本參數：Kriging 模型的超參數、初始實驗設計抽樣區域和容量N0、失效概率估計的目標變異系數 δ0、子區域內的初始樣本池容量?N0、主動學習停止條件的2 個閾值ε1和ε2，參數空間球分解的半徑序列為

2) 空間球分解和子區域Dm-1內抽樣：將外層超球環劃分為2 個子區域。令m=m+1，在次外層超球環Dm-1內生成?N0個樣本。

4) 空間球分解的停止條件：計算最外側超球環所需的樣本量=NMCSθi。若>?N0，則返回步驟2)；反之，若≤?N0，則繼續執行步驟5)。

5) 生成最外側超球環Dm內的個樣本。

3 算例分析

本節通過3 個算例驗證所提方法的有效性、計算精度和效率。3 個算例具有不同數量級的失效概率，用于考察可靠性分析方法對失效概率量級的魯棒性。AK-SDMCS 的計算參數取值如下：目標變異系數取為δ0=5%；Kriging 模型超參數初值取100，搜索范圍取[10-3, 103]；子區域內初始樣本量和后續樣本增量均為?N0=104；初始LHS 抽樣區域為[μ±6σ]；初始實驗設計容量取N0=12；主動學習停止條件的2 個閾值取ε1=ε2=0.001。作為對比，分別采用MCS 和AK-MCS 對這3 個算例進行可靠性分析。失效概率估計的目標變異系數均取 5%。AK-MCS 的計算參數為：初始實驗設計容量取50，采用隨機抽樣方式在區域[μ±6σ]內生成；Kriging 模型超參數初值取10，搜索范圍取[0.1, 20]。

對所有算例均采用MCS 計算結果作為參考解。AK-MCS 和AK-SDMCS 均獨立運行10 次，以獲取平均結果和經驗變異系數。通過比較所需的功能函數調用次數Ncall、失效概率變異系數和相對誤差來驗證所提方法的精度和效率。

3.1 算例1 (10-3 量級)

管道模型參數如表1 所示。假定管道處于正常工況，結構系數分別取SFb=2.3和SFm=2.7。材料的流變應力 σf和斷裂韌性JIC、工作壓力p、管道內壁環向表面缺陷的深度a和半長度角θ均假定為隨機變量，對應的概率分布參數如表2 所示。為便于驗證，這里參考文獻[14]將含缺陷截面處的薄膜應力和彎曲應力分別取為σm=9.625p+3.2，σb=18.38p+10.83。

表1 算例1 管道模型參數

表2 算例1 基本隨機變量分布參數

分別采用MCS、AK-MCS 和AK-SDMCS 計算該管道的失效概率，結果如表3 所示。可以發現，AK-MCS 和AK-SDMCS 均具有較低的變異系數（最大約3.41%）和相對誤差（最大約3.57%），所得的失效概率估計具有較高的精度。而在結構功能函數調用次數方面，二者均比MCS 減少了4 個數量級，驗證了主動學習可靠性分析方法較高的計算效率。和AK-MCS 相比，AK-SDMCS 僅需29.3%（不到1/3）的結構功能函數調用次數，即可獲得相近精度的失效概率估計，計算效率更高。

表3 算例1 可靠性分析結果

圖4 為本文方法在一次典型計算中的失效概率收斂曲線，其中紅色虛線為MCS 所得的失效概率參考值。AK-MCS 在迭代約130 次后接近參考值，共調用結構功能函數333 次；而AK-SDMCS則在迭代約40 次后就接近失效概率參考值，所需的結構功能函數調用次數僅為95，收斂速度更快，計算效率更高。

圖4 算例1 失效概率收斂曲線

此外，10 次獨立運行中，AK-MCS 和AKSDMCS 的樣本池容量均值分別為2.43×105和70 440。可以看出，AK-SDMCS 顯著減少了樣本池容量（為AK-MCS 的29.0%）。由于主動學習過程需要基于樣本池中的樣本重復進行Kriging 模型訓練和預測，故AK-SDMCS 也降低了代理模型主動學習的計算量。通過減少樣本池容量、參數空間劃分以及降低結構功能函數調用次數，AKSDMCS 進一步提高了可靠性分析的計算效率。

3.2 算例2 (10-4 量級)[14]

管道模型參數如表4 所示。結構系數同算例1，分別取SFb=2.3和SFm=2.7。基本隨機變量的概率分布參數如表5 所示。假定含缺陷截面處的薄膜應力和彎曲應力分別為σm=10.66p+5.3，σb=18.25p+7.64。

本算例中，AK-SDMCS 方法10 次獨立運行生成的樣本池容量均值為1.97×105。而AK-MCS 則需要生成容量為3.2×106的樣本池（按5%變異系數估算）。這時，Kriging 代理模型本身的訓練和預測就需要大量的計算資源，個人電腦（CPU：Intel Core i7-10510U；內存：32 GB；硬盤：250 GB固態）已無法進行計算。其根本原因在于，主動學習的每一輪迭代都需要基于當前實驗設計點更新Kriging 模型，并采用更新后的Kriging 計算樣本池內所有樣本的響應值和方差。這導致主動學習本身就需要大量計算資源，使得AK-MCS 無法進行，故表6 僅列出MCS 和AK-SDMCS 的計算結果。可以看出，AK-SDMCS 所得的失效概率估計變異系數與MCS 接近，相對誤差為3.96%，具有較高的計算精度。在計算量方面，AK-SDMCS所需的功能函數調用次數比MCS 減少了5 個數量級，顯著提高了可靠性分析的計算效率。圖5為AK-SDMCS 在一次典型計算中的失效概率收斂曲線。

表6 算例2 可靠性分析結果

圖5 算例2 失效概率收斂曲線

3.3 算例3 (10-5 量級)[15-16]

本算例管道模型的外徑377 mm、壁厚10 mm，其余確定性參數均與算例1 相同。基本隨機變量的概率分布參數如表7 所示。假定含缺陷截面處的薄膜應力和彎曲應力分別為σm=9.625p+5.4 , σb=25.24p+4.7。

表7 算例3 基本隨機變量分布參數

表8 列出了可靠性分析結果。10 次獨立運行中， AK-SDMCS 的樣本池容量均值為3.62×105。類似于算例2，AK-MCS 由于樣本池容量太大而無法進行，再次驗證了AK-MCS 不適用于小失效概率問題。本文所提AK-SDMCS 給出的失效概率估計變異系數與MCS 接近，相對誤差僅為0.62%，所需的結構功能函數調用次數僅為158.3，比MCS 減少了5 個數量級，計算效率有了明顯提高。圖6 為AK-SDMCS 在一次典型計算中的失效概率收斂曲線。

表8 算例3 可靠性分析結果

圖6 算例3 失效概率收斂曲線

比較3 個算例可以發現， AK-SDMCS 對10-5~10-3數量級范圍的失效概率均具有較高的計算效率，所需的結構功能函數調用次數為101.0~158.3，并沒有隨失效概率量級的增加而大幅增加。表明AK-SDMCS 方法對失效概率量級不敏感，具有較好的魯棒性。

4 結論

本文基于ASME 鍋爐與壓力容器規范建立了含缺陷核管道的可靠性分析模型，包括失效模式、主要隨機參數以及極限狀態方程。進而，基于Kriging 代理模型及參數空間球分解思想，建立了適用于含缺陷管道可靠性分析的主動學習代理模型方法（AK-SDMCS），并通過3 個不同數量級的算例驗證了所提方法的適用性、精度和效率。結果表明：

1) 所提方法只需少量功能函數調用次數即可給出高精度的失效概率估計（對10-5~10-3數量級范圍的失效概率，所需功能函數調用次數在102量級），具有較高的精度和效率。

2) 與主動學習方法AK-MCS 相比，所提方法既保留了AK-MCS 通用性和魯棒性的優點，同時通過空間球分解思想顯著降低了主動學習的計算量（特別是對小失效概率問題），克服了AK-MCS求解小失效概率時計算量過大的缺點。

3) 所提方法對失效概率數量級具有一定的魯棒性，所需的結構功能函數調用次數并未隨著失效概率數量級的增加而顯著提高。

4) 克服了現有方法面臨的計算量大、精度低等缺點，可為含缺陷核管道可靠性分析提供高效適用的途徑。