可靠性全局靈敏度求解的Meta-IS-AK算法

周蘇婷， 呂震宙， 凌春燕，王燕萍

西北工業大學 航空學院，西安 710072

全局靈敏度分析可以從輸入變量的整個分布范圍來衡量輸入變量的不確定性對工程設計中所感興趣的輸出性能統計特征的貢獻程度[1-2]，通常在可靠性分析中，更關注結構是否失效以及輸入變量對于結構失效概率的影響程度。為了衡量輸入變量對于結構失效概率的影響并對輸入變量的重要性進行排序，Cui等[3]提出了可靠性全局靈敏度指標。可靠性全局靈敏度指標定義為某一輸入變量固定時，結構的無條件失效概率與條件失效概率的差異在輸入變量的整個分布范圍內的平均值，該指標反映了輸入變量固定時對失效概率的平均影響。與可靠性局部靈敏度指標[4-5]不同，可靠性全局靈敏度指標可以從輸入變量的整個不確定性范圍來衡量其對于失效概率的平均影響。

本文提出了利用貝葉斯公式轉換以及元重要抽樣算法[6-7]結合自適應Kriging代理的可靠性全局靈敏度，并組織了求解可靠性全局靈敏度指標的高效算法。基于貝葉斯公式可以將可靠性全局靈敏度指標的原始定義等價表示為輸入變量的無條件概率密度函數(Probability Density Fanction, PDF)與條件概率密度函數之間的差異與失效概率的乘積形式，從而避免條件失效概率的求解[8-9]?；诖说葍r形式，只需要一組樣本，便可以得到各輸入變量的可靠性全局靈敏度指標。利用貝葉斯算法進行可靠性全局靈敏度分析的研究已有很多，Wang等[8]利用貝葉斯公式，將可靠性全局靈敏度指標的定義式進行了等價轉換，并采用蒙特卡羅法來計算可靠性全局靈敏度指標。該方法適用范圍雖廣，但在計算工程實際中的小失效概率問題時計算量過于龐大。Wang等[9]則結合重要抽樣方法來計算可靠性全局靈敏度指標的貝葉斯等價形式。該方法只需一組樣本便可以計算得到所有輸入變量的可靠性全局靈敏度指標，且計算量獨立于輸入變量的維數，明顯地提高了計算的效率。但是重要抽樣方法[10-11]需要構造重要抽樣密度函數，一般來說需要將重要抽樣密度函數的中心放在設計點處，這就使得該方法依賴于其他方法來尋找設計點，并且重要抽樣方法對于多設計點的情況和多失效域的情況并不適用。Yun等[12]基于可靠性全局靈敏度指標的貝葉斯公式轉換后的等價形式，利用子集模擬方法[13-16]結合重要抽樣的方法進行可靠性全局靈敏度的計算。該方法在利用子集模擬結合重要抽樣的方法計算得到無條件失效概率的同時，通過重復利用子集模擬重要抽樣方法抽取的樣本，并利用Metropolis-Hastings(M-H)[17]準則來實現失效樣本的轉換，進而估計各輸入變量的失效條件概率密度函數。該方法同樣只需要一組樣本便可計算得到可靠性全局靈敏度指標，且對于小失效概率問題是有效可行的。但該方法抽樣過程中所結合的重要抽樣方法仍然依賴于設計點的選取，因而該方法依舊很難適用于多設計點和多失效域的情況。

對于實際工程中普遍存在大量復雜的功能函數，自適應Kriging(Adaptive Kriging, AK)代理模型法是處理這類問題的高效算法[18]。可靠性分析中涉及到的自適應Kriging算法有自適應Kriging代理模型結合Monte Carlo模擬(Monte Carlo Simulation, MCS)法(AK-MCS)[19]、自適應Kriging代理模型結合重要抽樣(Important Sampling, IS)法(AK-IS)[20]以及元重要抽樣算法(Meta-IS)[6-7]等。AK-MCS算法將自適應Kriging代理模型與Monte Carlo模擬法相結合，其首先由輸入變量的聯合概率密度函數產生MCS備選樣本池，然后利用學習函數在備選樣本池中逐步挑選對失效面擬合貢獻較大的點來更新Kriging模型，最終能夠確保Kriging模型能在一定的置信水平下識別樣本池內樣本的功能函數值的正負號，從而較準確地求得失效概率。但對于工程上常見的小失效概率問題，AK-MCS算法的樣本池容量龐大，進而使得代理過程十分耗時。將重要抽樣法與自適應Kriging過程相結合形成的AK-IS算法可以大大降低備選樣本池的規模，從而提高自適應學習的效率，減少自適應學習過程所消耗的時間。但是對于多設計點及多失效域問題，基于一次二階矩法構造重要抽樣密度函數的AK-IS算法將失效，因而對于多設計點及多失效域的問題需要另尋解決途徑。元模型重要抽樣算法(Meta-IS)通過元模型構造重要抽樣密度函數，進而來抽取重要抽樣樣本點。該方法在降低AK-MCS樣本池規模的同時避免了設計點的求解，適用于多設計點和多失效域問題的可靠性分析。然而，Meta-IS算法在抽取重要抽樣樣本后仍需要求解重要抽樣樣本點處真實的功能函數值來估計失效概率，所以該算法仍有較大的計算量。如果在元重要抽樣的基礎上嵌入自適應Kriging模型，形成Meta-IS-AK算法，則有可能極大程度地提高失效概率求解的效率，本文將采用這種思路求解失效概率，并在此基礎上組織可靠性全局靈敏度的算法。

利用可靠性全局靈敏指標的貝葉斯轉換形式，只要能求得無條件失效概率和失效域條件下輸入變量的概率密度函數即可直接求得可靠性全局靈敏度，而求解無條件失效概率的數字模擬法可以同時完成失效概率的計算以及失效域條件下輸入變量概率密度函數的求解，為此本文將Meta-IS與自適應Kriging模型結合起來，以便完成可靠性全局靈敏度的高效求解。本文所組織的Meta-IS-AK算法主要分3步來執行。第1步是由元重要抽樣的迭代策略得到逐步逼近最優重要抽樣函數的樣本點，在迭代收斂后便可以得到第1步的Kriging模型以及相應的重要抽樣樣本點。第2步是基于重要抽樣樣本構建自適應Kriging模型，該步的Kriging模型是在第1步得到的Kriging模型的基礎上逐步更新得到的，其目的是為了構建能夠對第1步中重要抽樣樣本點失效與否做出準確預測的Kriging模型，從而高效地求解出無條件失效概率。然而Meta-IS-AK算法抽取的失效域內的樣本點的密度函數為最優重要抽樣密度函數，因而無法直接用于估計輸入變量失效域內的條件概率密度函數，為此第3步則對失效域內服從于重要抽樣密度函數的樣本點進行轉換，該轉換將利用Metropolis-Hastings準則由重要抽樣密度函數的失效樣本點來得到原概率密度函數在失效域內的樣本點，進而求得各個輸入變量在失效域中的條件概率密度函數。

本文主要由以下幾部分組成。第1節簡要地回顧了已有的利用貝葉斯公式轉換的可靠性全局靈敏度的形式，第2節詳細介紹了本文所提的計算可靠性全局靈敏度的高效算法。第3節分別給出了數值與工程算例，驗證了所提算法的高效性，第4節則得出結論。

1 基于貝葉斯公式的可靠性全局靈敏度指標的轉換形式

(1)

利用貝葉斯公式，可將條件失效概率Pf|Xi進行如下轉換[8-9,12,23-24]：

(2)

將式(2)代入到式(1)中，即可得到基于貝葉斯公式的可靠性全局靈敏度指標表達式：

(3)

2 求解可靠性全局靈敏度指標的高效算法

2.1 失效概率Pf的估計

(4)

步驟1以gK1(x)作為gK2(x)的初始模型，即令gK2(x)=gK1(x)。

步驟3計算g(xu)，將(xu,g(xu))添加到訓練集T中，形成更新的訓練集T={T∪(xu,g(xu))}。

步驟4由T更新gK2(x)，并返回步驟2。

(5)

(6)

(7)

2.2 失效域條件下輸入變量的概率密度函數估計

步驟2計算下列比值r,并由r確定下一個馬爾可夫鏈樣本

i=1,2,…,NF-1

(8)

i=1,2,…,NF-1

(9)

式中：u為[0,1]區間上服從均勻分布的隨機數。

2.3 Meta-IS-AK算法求解可靠性全局靈敏度指標流程圖

由2.1和2.2節的內容可給出Meta-IS-AK算法求解可靠性全局靈敏度的流程如圖1所示。

圖1 Meta-IS-AK算法求解的流程圖Fig.1 Flowchart of estimating Meta-IS-AK algorithm

3 算例分析

為了驗證所提算法在求解可靠性全局靈敏度的效率和精度，本節給出了3個算例，分別采用準Monte Carlo算法(QMC)、AK-MCS算法、Meta-IS算法及Meta-IS-AK算法求解可靠性全局靈敏度指標。其中QMC算法結果將作為參照解。

3.1 不連續多失效域算例

考慮一個多失效域的算例，其功能函數為

(10)

(11)

圖2 算例3.1失效邊界Fig.2 Failure boundary of Example 3.1

表1 算例3.1可靠性全局靈敏度指標的計算結果

Table 1Results of reliability global sensitivity index for Example 3.1

算法X1/10-3X2/10-3NcallQMC14.5[0.30]14.5[0.46]1×105AK-MCS14.2[0.23]14.1[0.35]440Meta-IS14.2[0.41]14.0[0.65]42+5 000Meta-IS-AK14.2[0.43]13.9[0.65]42+254

注：數據右上標為可靠性全局靈敏度指標估計值的標準差,Ncall為功能函數的調用次數。

3.2 鋼架結構算例

對于一鋼架結構，有4個失效模式：

(12)

因為是串聯系統，所以系統的功能函數g為g=min{g1,g2,g3,g4}。在各個失效模式的功能函數中，Mi(i=1,2,3)與S是獨立的，均值和標準差分別為:μMi=2(i=1,2,3),μS=1,σMi=2(i=1,2,3),σS=0.25。各輸入變量所對應的可靠性全局靈敏度指標的計算結果如表2所示。

表2 算例3.2可靠性全局靈敏度指標的計算結果

由表2中結果可以看出，Meta-IS-AK算法得到的計算結果與QMC算法計算所得結果基本一致，且模型的調用次數低于QMC算法與AK-MCS算法以及Meta-IS算法。這說明Meta-IS-AK算法能夠很好地適用于多個失效模式的情況。對于此算例，輸入變量的重要性排序為S>M3>M1>M2。此外，輸入變量S對失效概率的影響較其他輸入變量的影響大很多。因此在進行可靠性設計時通過調整輸入變量S的不確定性可以更有效地滿足可靠性要求。此外，M2對失效概率影響的靈敏度指標較其他的輸入變量小很多，因此可以忽略M2的不確定性對于失效概率的影響，可以通過將M2固定在均值處來簡化失效概率的求解模型。

3.3 簡化的翼盒結構算例(多設計點算例)

圖3(a)是簡化翼盒模型的示意圖。這種翼盒結構由64個桿和42個板組成。64個桿根據它們的方向分為3組，x、y和z方向上桿的長度分別是2L、L和 3L。所有桿的截面積均為A，所有板的厚度均為TH，E和P分別是所有板和桿的彈性模量和外部載荷。泊松比為0.3。假設輸入變量是獨立的正態變量，其分布參數如表3所示。

圖3 翼盒結構Fig.3 Wing box structure

表3 算例3.3輸入變量的分布參數

Table 3Distribution parameters of input variables of Example 3.3

變量分布均值變異系數A/m2正態1×10-40.1L/m正態0.20.1E/Pa正態7.1×10100.12P/N正態1 5000.1TH/m正態2.5×10-30.15

由于翼盒結構承擔了來自機翼的大部分垂直方向的載荷，其自由端的位移可能相對較大。大的位移可能會引起桿的變形并導致桿或板的破壞，這在飛行過程中是不允許的。因此計算翼盒結構的位移是很有必要的，且位移的計算可以通過有限元分析完成。圖3(b)顯示了在ANSYS 14.0中構建的翼盒結構的有限元模型，其中輸入變量固定在它們的均值處，此外，翼盒的變形如圖3(c)所示。

本文中，該翼盒的隱式功能函數是根據翼盒的最大位移Dmax不超過臨界值構造的(這里臨界值等于0.01)，其表達式為

g(X)=0.01-|Dmax|

(13)

式中：Dmax=D(L,A,E,P,TH)，是關于輸入變量的隱式函數，由有限元分析確定。算例3.3的可靠性全局靈敏度計算結果如圖4所示，計算結果的標準差如表4所示。

圖4 算例3.3可靠性全局靈敏度指標的計算結果Fig.4 Results of reliability global sensitivity index for Example 3.3

表4 算例3.3可靠性全局靈敏度指標結果

Table 4Reliability global sensitivity index forExample 3.3

變量SDi/10-5QMCAK-MCSMeta-ISMeta-IS-AKA/m20.270.200.170.18L/m0.260.110.140.11E/Pa5.363.875.045.05P/N0.350.420.330.32TH/m8.178.856.376.30Ncall5×10512033+5 00033+59

由表4中的結果可以看出，Meta-IS-AK算法得到的計算結果精度較高且模型的調用次數低于QMC算法與AK-MCS算法以及Meta-IS算法，這說明了Meta-IS-AK算法能夠很好地適用于隱式功能函數的情況。對于算例3.3，輸入變量的的不確定性對于失效概率影響的重要性排序為TH>E>P>L>A。從圖4的結果可以看出，輸入變量TH對失效概率的影響較其他輸入變量的影響大很多，因此在進行可靠性設計時調整輸入變量TH的不確定性可以更有效地滿足可靠性要求。從圖4的結果還可以看出，L和A的不確定性對于失效概率的影響較其他的輸入變量小很多，因此忽略L和A的不確定性對于失效概率的影響可以簡化失效概率的求解模型。

4 結 論

1) 可靠性全局靈敏度指標可以有效地衡量輸入變量對于結構失效概率的影響。為了高效地計算可靠性全局靈敏度指標，本文提出了可靠性全局靈敏度指標計算的Meta-IS-AK算法。

2) Meta-IS-AK算法在計算可靠性全局靈敏度指標時充分利用了可靠性全局靈敏度指標貝葉斯算法的維度獨立性，也使得該指標的求解轉化為輸出樣本的分類問題，這就使可靠性全局靈敏度指標可以通過嵌入式代理模型結合數字模擬法來求解。相比于傳統的重要抽樣函數，元重要抽樣法構造的準最優抽樣密度函數能夠很好地適用于多設計點和多失效域的情況。在此基礎上構造重構的功能函數的代理模型能夠準確地對準重要抽樣密度函數抽取的樣本進行分類，進一步地提高了失效概率的計算效率。M-H準則無需額外調用功能函數便可以將準重要抽樣密度函數的失效樣本轉化成原概率密度函數的失效樣本，大大提高了求解輸入變量的失效條件概率密度函數的效率。

3) 所提方法能夠很好地適用于多設計點和多失效域的情況，也適用于隱函數問題的可靠性靈敏度分析。所提方法的計算效率遠遠高于QMC算法，同時較AK-MCS算法以及Meta-IS算法也有提高。所給算例充分驗證了以上結論。