柳詩雨,呂震宙*,員婉瑩,肖思男
西北工業大學 航空學院,西安 710072
小失效概率情況下的全局可靠性靈敏度分析的高效方法
柳詩雨,呂震宙*,員婉瑩,肖思男
西北工業大學 航空學院,西安 710072
針對目前很多算法都無法準確、高效地計算小失效概率(10-4,甚至更小)情況下的全局可靠性靈敏度問題,本文提出了一種高效求解小失效概率情況下的全局可靠性靈敏度新算法。所提算法通過擴大標準差構造重要抽樣密度函數來進行空間分割(SP),再與無跡變換(UT)結合,利用函數在分割后的子空間內非線性程度的降低和無跡變換方法可以高效計算低非線性程度函數的前二階矩,來高效準確地計算小失效概率情況下的全局可靠性靈敏度。所提算法的優點有:重要抽樣密度函數的選擇可以使得空間分割時向重要區域偏移,并且在分割區域內功能函數的復雜性被降低,從而可以利用無跡變換方法高效計算失效概率,進而高效求得全局可靠性靈敏度。與已有的算法相比,算例說明了本文所提方法的優勢。
小失效概率;全局可靠性靈敏度;重要抽樣;空間分割;無跡變換
全局可靠性靈敏度分析是可靠性設計中非常重要的一項工作,它可以衡量輸入變量不確定性對失效概率的貢獻程度[1],因而近年來,國內外越來越多的學者將精力投入到全局可靠性靈敏度分析的研究中。Cui等[2-3]定義了全局可靠性靈敏度指標,分別用概率密度演化和態相關參數(SDP)方法來求解。在文獻[4]中,Wei等先采用Saltelli等提出的單層蒙特卡羅(MC)法[5-6]計算全局可靠性靈敏度,對于失效概率較大(10-2~10-1)的問題,單層MC法效率和精度均比較高;但對于小失效概率問題,為獲得較高的精度,就必須抽取大量的樣本以保證有足夠的樣本落入失效域中,因此效率較低。針對此問題,Wei等[4]又提出了采用單層重要抽樣(Importance Sampling,IS)法[7-9]和截斷重要抽樣(Truncated Importance Sampling,TIS)法[10]來提高抽樣的效率,其基本思想是通過改變抽樣概率密度函數,將抽樣中心轉移到設計點處,使得更多樣本點落入失效域從而加速收斂速度[11],取得了一定的效果,但是重要抽樣密度函數的確定依賴于設計點的求解,而求解設計點對于隱式功能函數來說較為困難。任博等[12]在IS法基礎上,提出一種改進估計量與IS法相結合的方法求解全局可靠性靈敏度指標。以上大多數方法都是基于數字模擬方法,對失效概率較小(10-4,甚至更小)的這類問題,仍需要隨機抽取大量樣本點才能使結果收斂,效率較低。
文獻[13]中提出了空間劃分的思想求解基于方差的全局可靠性靈敏度指標,該方法可以最大程度利用數字模擬樣本點從而減少計算量。受此文思想的啟發,本文提出空間分割(Space-Partition,SP)的數值積分法,該方法首先通過擴大的輸入變量的標準差來構造重要抽樣密度函數,以此重要抽樣密度函數來劃分空間,從而使得劃分空間向重要的失效域偏移,然后再在每個子空間中應用低階無跡變換(Unscented Transformation,UT)來求解全局可靠性靈敏度指標。與數字模擬方法相比,所提方法計算量不受失效概率大小影響,在保證計算結果同等精度條件下,對小失效概率(10-4,甚至更小)的模型,所提方法的計算效率大幅提升,本文最后采用數值算例及工程算例驗證了所提方法的合理性和高效性。
對于結構極限狀態函數Y=g(X),X=[X1X2… Xn]為n維輸入隨機變量,且Y>0時結構安全,Y<0時結構失效,此時輸入變量Xi的全局可靠性靈敏度可定義為[6]

式中:EXi(·)為Xi的期望算子;PfY為結構的無條件失效概率;PfY|X為當變量Xi取其實現值時
iY的條件失效概率;Si為單個輸入變量Xi的不確定性對失效概率PfY的貢獻。因此若Si越大,表示Xi的不確定性對PfY的貢獻越大。在工程實際中,可以通過減少Si比較大的變量Xi的不確定性減小失效概率PfY。
由式(1)可知,為計算全局可靠性靈敏度指標Si需計算無條件的失效概率PfY和條件失效概率平方的期望EXi(),EXi()的求解包含內層條件失效概率和外層期望的求解。和可以由式(2)和式(3)中所示的無條件指示函數IF(Y)和條件指示函數IF(Y |Xi)的期望求得。

式中:fX(x)為Xi的聯合概率密度函數;fX~i|Xi(x~i)為給定Xi條件下除Xi以外的剩余輸入變量X~i的聯合概率密度函數;E(·)為期望算子;且

以往采用矩方法(數值積分法的一種)來求解失效概率都是先求極限狀態函數Y=g(X)的各階矩,然后再由g(X)的各階矩來逼近g(X)的概率密度進而求得失效概率[14-16],以往矩方法不直接采用特征點和權函數來求E (IF)的原因是指示函數IF(Y)在整個變量空間中的性態比較復雜,很難采用特征點和權函數來捕捉IF(Y)的足夠信息以求得E (IF)。為了避免以往方法多次近似產生的多重誤差,本文采用空間分割的方法,由于局部區域IF(Y)的函數復雜性程度往往遠低于全域函數的復雜程度,這將使得在局部子空間中可以采用特征點和相應的權函數高精度求得E (IF)。相比于目前很多先求極限狀態函數的各階矩然后再逼近失效概率的方法,本文直接采用特征點和相應的權函數求解E (IF)的方法的誤差來源要少得多。
本文求解PfY的方法包括兩部分:其一是利用積分空間的可加性將求解PfY的整個變量空間積分等價變換成子空間積分的和,由于子空間中被積函數IF(Y)的復雜性遠低于整個空間的復雜性,從而可以利用特征點和權函數來求得子空間中的積分;其二是劃分子空間時構造了重要抽樣密度函數,按照重要抽樣密度進行子空間的劃分可以使得劃分的子空間向對PfY貢獻較大的區域移動,從而提高算法的效率精度。
2.2.1 重要抽樣密度函數的構造
對于工程中常見的小失效概率問題,常用到數字模擬法中的IS法,該方法的基本思想就是引入重要抽樣密度函數hX(x)代替原來的抽樣密度函數fX(x),從而使得對失效概率貢獻大的樣本點能以較大的概率出現,進而減小數字模擬法的估計方差。
本文借鑒IS法來劃分子空間,其基本思想是使得對失效概率貢獻大的區域能夠有更多的子空間,從而提高算法的精度和效率。構造重要抽樣密度函數的方法有很多,選擇其中較為容易實現的一種,該方法只需將原來抽樣密度函數fX(x)變量的標準差σf擴大r倍,就可構成重要抽樣密度函數hX(x),即hX(x ,σn)=fX(x ,rσf),然后按hX(x)進行等概率子空間的劃分,并由下述過程實現PfY的求解。
2.2.2 基于SP的PfY的計算
按照上述方法構造的hX(x)對每一維輸入變量進行區間劃分,設Xi(i=1,2,…,n)被劃分為Ni個子區間Aki=(aki-1,aki)(ki=1,2,…,Ni;i=1,2,…,n),則在此劃分下,失效概率PfY的計算公式為


因此,式(6)可改寫為


因此,式(8)又等價為

值W(j)X及相應的n維Sigma點x(j)(j=0,1,…,2n):


與以前采用特征點及權函數直接求解Y=g(X)的各階矩然后由Y的各階矩逼近Y的密度函數進而再近似求得失效概率的算法相比,本文方法的誤差來源少,只是在子空間做積分時產生誤差,但由于子空間中函數的非線性程度較低,因而可以得到較高的近似精度。


斯積分點及對應權重系數可查閱文獻[19]。

綜上,本文計算全局可靠性靈敏度指標Si的計算流程簡要總結如下:
步驟1 通過擴大輸入變量X的標準差來構造重要密度函數hX(x),隨后根據此密度函數hX(x)將整個積分空間等概率地劃分為(i=1,2,…,n)個子空間,在每個子空間內應用UT求解失效概率PfY。
步驟2 根據一維輸入變量Xi的重要抽樣密度函數hXi(xi)產生M 個高斯點(xi1,xi2,…,xiM)及其相應的權重(wi1,wi2,…,wiM)。
步驟3 固定變量Xi于xil(l=1,2,…,M)處,除Xi外變量X~i根據其重要抽樣密度hX~i(x~i)將 整 個 積 分 空 間 等 概 率 地 劃 分 為(j=1,2,…,i-1,i+1,…,n)個子空間,在每個子空間內應用UT求解條件失效概率
步驟4 將步驟3重復M次,得到變量Xi固定于高斯點(xi1,xi2…,xiM)處的條件失效概率的(l=1,2,…,M),對上 述(l=1,2,…,M)進行期望求解。
步驟5 將步驟1得到的PfY和步驟4得到的EXi()代入式(1)即可計算得到基本輸入變量Xi的全局可靠性靈敏度指標Si。
將每一維變量空間分別劃分Ni個子區間,Ni個子空間。在每個子空間內應用UT,需產生2n+1個Sigma點,由于選擇標準 UT,即W0=0,因此舍去第一個Sigma點,即:在每個子空間中僅需要調用模型2n次。因此,計算PfY的計算量為整個積分空間被劃分為,計算)時,外層采用M個高斯點,內層采用SP-UT法的計算量為,總計算量為,于是計算Si的總計算量為
綜上所述,對于低維小失效概率(10-4或更小)情況下的全局可靠性靈敏度的求解,為獲得較高的精度,MC法或IS法等數字模擬方法必須抽取大量的樣本以保證有足夠的樣本落入失效域中,因此效率較低。而本文所提方法的計算量受失效概率的大小影響很小,主要與函數的非線性程度和輸入變量的維度有關,因此計算效率較高,在小失效概率(10-4或更小)的全局可靠性靈敏度的計算上有明顯優勢。
采用一個數值算例和兩個工程算例來說明本文所提方法(Proposed)在計算小失效概率(10-5~10-3)問題上的合理性和高效性。在對比結果中,使用MC法在大樣本抽樣(9×106)情況下的結果作為近似精確解,以說明所提算法的精度。另外,引入文獻[4]提出的IS法進行比較,以說明所提算法的高效性。本文采用9點高斯積分。N用來表示每種方法的模型計算量。
圖1給出了在空間分割情況一致的條件下,本文所提直接利用特征值和權函數求E (IF)和通過求 g(x)的四階矩 (Fourth-Order Moment,FOM)來逼近E (IF)這兩種求解失效概率PfY的方法隨模型調用次數的變化規律對照圖。從圖1中可以看出,與直接計算指示函數的均值的方法相比,通過計算極限狀態函數g(x)的四階矩來逼近得到的失效概率與精確解相差很大,且收斂速度很慢,從這個意義上講,本文所提出的通過空間分割直接復雜函數IF(x)的均值的方法大大提高了失效概率的計算效率和精度。

圖1 本文所提方法和四階矩(FOM)方法求解失效概率PfY的收斂趨勢圖Fig.1 Schematic comparison of convergence rates of failure probability PfYcalculated by proposed and fourth-order moment(FOM)methods
圖2和圖3分別給出了不同標準差擴大倍數r的情況下,失效概率PfY和全局可靠性靈敏度指標S1隨著空間分割數目r的增加而變化的曲線。從圖2可以明顯看出,當采用原概率密度函數(r=1)時,失效概率始終為0,即無失效情況出現,而通過重要抽樣密度函數(r=2~5)來劃分空間,每一維空間劃分數目從5個開始即可出現失效。從兩張圖中都可以看出,當r取2~5時,隨著空間分割數的增多,所提方法都能較快較好地收斂到其精確解。從這個意義講,通過擴大標準差構造重要抽樣函數使得按重要抽樣密度函數抽取的樣本點有更大的概率落入失效域,從而減小估計值的方差,提高了計算精度和速度。
雖然從理論上講,擴大系數取大于1的值就可以提高數字模擬法的效率,但是研究表明,如果r值取得過大,就會使得過多的樣本點又落入了對失效概率貢獻不大的失效域,這也不利于計算效率的提高[20-21]。然而最優的放大系數的搜索涉及額外的優化計算,并且這種尋優過程所耗費的計算量有可能較大。因此如何合理地確定標準差放大系數,還有必要做進一步的理論和數值研究。本文算例中均采用r=2,建議標準差擴大倍數r取2~5。

圖2 標準差放大系數r對失效概率PfY估計值的影響Fig.2 Influence of amplification ratios of standarddeviation r on failure probability PfY

圖3 標準差放大系數r對全局可靠性靈敏度指標S1的影響Fig.3 Influence of amplification ratios of standard deviation r on global reliability sensitivity indices S1
表1給出了采用三種方法求解算例1的全局可靠性靈敏度指標的結果。從表1的結果可以直觀地看出,本文所提方法計算得到的全局可靠性靈敏度指標Si的值與文獻[4]的IS法及精確解MC法的計算結果相吻合,從而說明了本文方法的正確性。而在計算量上,本文所提方法優越很多,只需要調用功能函數1 170次,而MC法計算量非常龐大,為4.5×107。IS法較MC法效率提高很多,但仍然需要2×104次的模型計算量。為了進一步說明本文所提方法的高效性,圖4給出了本文所提方法和IS法計算S1指標的收斂趨勢圖,從圖4中可以看出,本文所提方法計算S1指標的收斂速度比IS法快得多。

表1 算例1全局可靠性靈敏度指標計算結果Table 1 Estimates for global reliability sensitivity indices of Example 1

圖4 本文所提方法和IS法計算算例1的全局可靠性靈敏度指標S1收斂趨勢圖Fig.4 Schematic comparison of convergence rates of global reliability sensitivity indices S1calculated by proposed and IS methods in Example 1
算例2 如圖5所示一矩形截面懸臂梁受到均布載荷,以其自由端撓度不超過L/325為約束建立極限狀態函數為g(ω,b)=L/325-ωbL4/8EI式中:ω、b、L、E、I分別為單位載荷、截面尺寸、梁的長度、彈性模量和截面慣性矩,其中L和E為已知常量,L=6m,E=26GPa,I=b4/12,將極限狀態函數簡化為g(ω,b)=0.018 461 54-74.769 12ω/b3其中:ω~N(1 000,1002),b~N(250,252)。
本工程算例分析了輸入變量對矩形截面懸臂梁模型失效概率影響的全局可靠性靈敏度指標Si。從表2的結果可以直觀地看出,在計算工程算例中,本文所提方法計算精度與MC法和IS法一致,而在計算量上,則遠優于MC法和IS法,這也證明了本文所提方法同樣適用于工程算例。
為了進一步說明本文所提方法的高效性,圖6給出了IS法和本文所提方法計算S2指標的收斂趨勢圖,從圖中可以看出,本文所提方法計算S2指標的收斂速度比IS法快得多。

圖5 矩形截面懸臂梁Fig.5 Rectangular-section cantilever beam

表2 算例2全局可靠性靈敏度指標計算結果Table 2 Estimates for global reliability sensitivity indices of Example 2

圖6 本文所提方法和IS法計算算例2的全局可靠性靈敏度指標S2收斂趨勢圖Fig.6 Schematic comparison of convergence rates of global reliability sensitivity indices S2calculated by proposed and IS methods in Example 2
算例3 如圖7所示的十桿結構,其中水平桿和豎直桿的長度均為L;每根桿的截面積為Ai=0.001(i=1,2…,10);彈性模量為E;P1、P2和P3為作用在圖上所示位置的外載荷,P2=P3=100kN。設L、E和P1服從正態分布,其均值為L=1m,E=100GPa,P1=800kN,變異系數均為0.05。當2節點縱向位移超過0.04m時認為結構失效,因此極限狀態函數為
由MC法求得該結構失效概率為4.6×10-5,屬于小失效概率情況。下面采用本文所提方法、IS法和MC法計算各變量的全局可靠性靈敏度指標,計算結果如圖8所示。從圖8的全局可靠性靈敏度指標對比結果可以看出,對于工程中的隱式非線性功能函數情況,文獻[6]所提的IS法的計算結果出現了明顯的偏差,原因可能是重要抽樣密度函數的確定依賴于設計點的求解,而求解設計點對于隱式功能函數來說較為困難。而本文所提方法的結果與MC法基本一致,證明本文所提方法在計算工程中的隱式非線性功能函數的全局可靠性靈敏度指標上仍然保證很好的計算精度。

圖7 平面十桿桁架結構Fig.7 Planer 10-bar structure

圖8 算例3全局可靠性靈敏度指標計算結果Fig.8 Estimates for global reliability sensitivity indices of Example 3
1)將空間分割(SP)和無跡變換(UT)方法結合,提出了高效計算全局可靠性靈敏度指標的方法。該方法計算量受失效概率的大小影響很小,而與輸入變量的維度和函數的非線性程度相關,因此在計算小失效概率(10-5~10-3)情況下的全局可靠性靈敏度指標上,特別是在輸入變量維度較低時很有優勢。
2)數值算例和工程算例都表明,與IS法相比,該方法不僅在求解全局可靠性靈敏度指標時具有更高的收斂速度,而且還適用于工程中的隱式非線性功能函數情況。
[1] 宋述芳,呂震宙.系統可靠性靈敏度分析方法及其應用研究[J].機械強度,2007,29(1):53-57.SONG S F,LU Z Z.Reliability sensitivity analysis method for structural system and its application[J].Journal of Mechanical Strength,2007,29(1):53-57(in Chinese).
[2] CUI L J,LU Z Z,WANG W.Importance measures of basic random variable and their probability density evolution solutions[J].Journal of Nanjing University of Aeronautics and Astronautics,2011,43(2):165-171.
[3] CUI L J,LU Z Z,ZHAO X.Moment-independent importance measure of basic random variable and its probability density evolution solution[J].Science China-Technological Sciences,2010,53(4):1138-1145.
[4] WEI P F,LU Z Z,HAO W R,et al.Efficient sampling methods for global reliability sensitivity analysis[J].Computer Physics Communications,2012,183(8):1728-1743.
[5] HOMMA T,SALTELLI A.Importance measures in global sensitivity analysis of nonlinear models[J].Reliability Engineering &System Safety,1996,52(1):1-17.
[6] SALTELLI A,RATTO M,ANDRES T,et al.Global sensitivity analysis: The primer [M ]. Chichester:Wiley,2008.
[7] AU S K,BECK J L.A new adaptive importance sampling scheme for reliability calculations[J].Structural Safety,1999,21(2):135-158.
[8] HARBITZ A.An efficient sampling method for probability of failure calculation[J].Structural Safety,1986,3(2):109-115.
[9] MELCHER R E.Importance sampling in structural systems[J].Structural Safety,1989,6(1):3-10.
[10] GROOTEMAN F.Adaptive radial-based importance sampling method for structural reliability[J].Structural Safety,2008,30(6):533-542.
[11] HASOFER A M,LIND N C.Exact and invariant secondmoment code format[J].Journal of the Engineering Mechanics Division,1974,100(1):111-121.
[12] 任博,呂震宙,呂召燕,等.失效概率全局靈敏度分析的改進重要抽樣法[J].機械強度,2014,36(2):193-200.REN B,LU Z Z,LU Z Y,et al.An improved importance sampling method for failure probability based global sensitivity analysis[J].Journal of Mechanical Strength,2014,36(2):193-200(in Chinese).
[13] ZHAI Q,YANG J,ZHAO Y.Space-partition method for the variance-based sensitivity analysis:Optimal partition scheme and comparative study[J].Reliability Engineering &System Safety,2014,131:66-82.
[14] ZHAO Y,LU Z.Applicable range of the fourth-moment method for structural reliability[J].Journal of Asian Architecture and Building Engineering,2007,6(1):151-158.
[15] ZHAO Y G,ONO T.Moment methods for structural reliability[J].Structural Safety,2001,23(1):47-75.
[16] ZHAO Y G,ONO T.On the problems of the fourth moment method[J].Structural Safety,2004,26(3):343-347.
[17] RICHE J H.Reliability estimation using unscented transformation[C]/3rd International Workshop on Dependable Control of Discrete Systems.Piscataway,NJ:IEEE Press,2011.
[18] JULIER S J,UHLMANN J K.Unscented filtering and nonlinear estimation[J].Proceedings of the IEEE,2004,92(3):401-422.
[19] ZHANG X,PANDEY M D.Structural reliability analysis based on the concepts of entropy,fractional moment and dimensional reduction method[J].Structural Safety,2013,43(0):28-40.
[20] 吳斌,歐進萍,張紀剛,等.結構動力可靠度的重要抽樣法[J].計算力學學報,2001,18(4):478-482.WU B,OU J P,ZHANG J G,et al.Importance sampling techniques in dynamical structural reliability[J].Chinese Journal of Computational Mechanics,2001,18(4):478-482(in Chinese).
[21] 陸立新,吳斌,歐進萍.結構動力可靠度的重要抽樣法的探討[J].世界地震工程,2007,23(3):74-78.LU L X,WU B,OU J P.Study on an importance sampling method of structural dynamic reliability[J].World Earthquake Engineering,2007,23(3):74-78 (in Chinese).
Efficient method for global reliability sensitivity analysis with small failure probability
LIU Shiyu,LYU Zhenzhou*,YUN Wanying,XIAO Sinan
School of Aeronautics,Northwestern Polytechnical University,Xi’an 710072,China
At present,there are many methods for the estimation of global reliability sensitivity.However,these methods cannot efficiently and accurately estimate the global reliability probability in case of small failure probability(10-4or smaller).In this work,a highly efficient method to compute the global reliability sensitivity is proposed for the small failure probability.The proposed method combines the space-partition(SP)with unscented transformation(UT)which can obtain the first two moments of lowly nonlinear response function.The importance sampling density function,which is constructed by increasing the standard deviation,is employed to partition the input space into a series of subspaces,and thus the subspaces partitioned by the constructed importance sampling density function can move to the important area for the failure probability.Because the complexity of response function is reduced in the partitioned subspace,in which UT can estimate effectively the failure probability,the proposed method can estimate the global reliability sensitivity indices efficiently.All the above contribute to the efficiency and accuracy of the proposed method to compute the global reliability sensitivity.In this paper,the proposed method is compared with the existing methods and examples,and it is shown that the proposed method outperforms the others.
small failure probability;global reliability sensitivity;importance sampling;space-partition;unscented transformation
2015-09-16;Revised:2015-11-22;Accepted:2016-01-25;Published online:2016-03-09 14:47
URL:www.cnki.net/kcms/detail/11.1929.V.20160309.1447.006.html
s:National Natural Science Foundation of China(51475370);the Fundamental Research Funds for the Central Universities(3102015BJ(II)CG009)
V19;TB114.3
A
1000-6893(2016)09-2766-09
10.7527/S1000-6893.2016.0029
2015-09-16;退修日期:2015-11-22;錄用日期:2016-01-25;網絡出版時間:2016-03-09 14:47
www.cnki.net/kcms/detail/11.1929.V.20160309.1447.006.html
國家自然科學基金(51475370);中央高校基本科研業務費專項資金(3102015BJ(II)CG009)
*通訊作者.Tel.:029-88460480 E-mail:zhenzhoulu@nwpu.edu.cn
柳詩雨,呂震宙,員婉瑩,等.小失效概率情況下的全局可靠性靈敏度分析的高效方法[J].航空學報,2016,37(9):2766-2774.LIU S Y,LYU Z Z,YUN W Y,et al.Efficient method for global reliability sensitivity analysis with small failure probability[J].Acta Aeronautica et Astronautica Sinica,2016,37(9):27662-774.
柳詩雨 女,碩士研究生。主要研究方向:飛行器設計及可靠性工程。
Tel.:029-88460480
E-mail:shiyu_liu22@163.com
呂震宙 女,博士,教授,博士生導師。主要研究方向:飛行器設計及可靠性工程。
Tel.:029-88460480
E-mail:zhenzhoulu@nwpu.edu.cn
員婉瑩 女,博士研究生。主要研究方向:飛行器設計及可靠性工程。
Tel.:029-88460480
E-mail:wanying_yun@163.com
肖思男 男,博士研究生。主要研究方向:飛行器設計及可靠性工程。
Tel.:029-88460480
E-mail:bruce1209@163.com
*Corresponding author.Tel.:029-88460480 E-mail:zhenzhoulu@nwpu.edu.cn