基于設(shè)計驗算點擬合的響應(yīng)面法可靠性分析

陳鵬，朱翊洲，謝永誠

（上海核工程研究設(shè)計院有限公司，上海 200233）

0 引言

傳統(tǒng)的機(jī)械結(jié)構(gòu)設(shè)計多采用確定性分析方法，所有變量及參數(shù)不考慮隨機(jī)性和不確定性，為了保證在工程應(yīng)用中的可靠性，通常在設(shè)計中考慮一定的安全系數(shù)以增加設(shè)計保守性。近年來，隨著對設(shè)計精度和可靠性要求的不斷提高，很多機(jī)械結(jié)構(gòu)在設(shè)計中需要在考慮隨機(jī)不確定性的基礎(chǔ)上開展結(jié)構(gòu)概率可靠性研究[1-2]。

不同于確定性結(jié)構(gòu)分析，概率可靠性分析需要在概率分布的基礎(chǔ)上開展更為復(fù)雜的統(tǒng)計學(xué)分析，因此需要計算效率較高的結(jié)構(gòu)分析模型。目前工程中應(yīng)用最為廣泛的結(jié)構(gòu)可靠性方法是以一次2階矩法[3]為基礎(chǔ)的1階可靠性方法（FORM），該方法要求顯式的結(jié)構(gòu)失效功能函數(shù)，然而對于大多數(shù)工程問題，都難以直接給出以基本隨機(jī)變量為變量的功能函數(shù)，如壓力容器斷裂失效，螺紋局部疲勞失效等。當(dāng)結(jié)構(gòu)形式復(fù)雜時，基于ANSYS等有限元模型開展Monte-Carlo數(shù)值實驗[4]被認(rèn)為是精度最高的結(jié)構(gòu)可靠性分析方法，但當(dāng)結(jié)構(gòu)規(guī)模較大時，大量的數(shù)值實驗需要仿真時間難以估量的模擬實驗。

為了提高求解效率以適用于工程分析，對于復(fù)雜的機(jī)械機(jī)構(gòu)，目前使用最多的是用響應(yīng)面法對真實失效曲面進(jìn)行擬合，由于真實失效曲面的復(fù)雜性，在全域擬合得到精度較高的近似失效曲面是不現(xiàn)實的，傳統(tǒng)響應(yīng)面法在各隨機(jī)變量均值點附近開展擬合實驗，當(dāng)真實失效曲面非線性程度不高時，可以得到精度較高的擬合失效面和結(jié)構(gòu)可靠度，但對于大多數(shù)結(jié)構(gòu)可靠性問題，基于均值點擬合的響應(yīng)面只能保證均值點局部的曲面近似精度。為了提高擬合精度，需要在真實失效曲面上的失效點局部區(qū)域進(jìn)行擬合實驗，而真實失效點通常只能夠通過大量迭代不斷逼近，對于復(fù)雜問題，計算效率并不高。

本文基于復(fù)雜結(jié)構(gòu)的確定性有限元模型和響應(yīng)面法，擬合得到結(jié)構(gòu)失效的近似功能函數(shù)，根據(jù)1階可靠性方法中的設(shè)計驗算點法計算結(jié)構(gòu)失效的近似驗算點，并在該點局部區(qū)域重新開展近似實驗得到基于設(shè)計驗算點的結(jié)構(gòu)響應(yīng)面功能函數(shù)，在此基礎(chǔ)上根據(jù)矩法計算結(jié)構(gòu)可靠度，通過數(shù)值仿真和基于ANSYS的Monte-Carlo實驗結(jié)果，驗證了本文方法的計算精度。

1 基于響應(yīng)面法的功能函數(shù)擬合

對于工程中的結(jié)構(gòu)可靠性問題，描述結(jié)構(gòu)失效的功能函數(shù)通常可以表示為

式中：S為基于相關(guān)強(qiáng)度理論和設(shè)計規(guī)范的設(shè)計應(yīng)力強(qiáng)度值；σ為根據(jù)失效模式定義的目標(biāo)應(yīng)力值；x為影響應(yīng)力分量σ的所有隨機(jī)變量，當(dāng)S也是隨機(jī)變量時，功能函數(shù)對應(yīng)的基本隨機(jī)變量x應(yīng)包含S。

在工程中，對于大多數(shù)結(jié)構(gòu)都無法給出目標(biāo)應(yīng)力分量σ的顯式表達(dá)式，更多的只能采用有限元等復(fù)雜的求解方法。對于結(jié)構(gòu)可靠性問題，要得到結(jié)構(gòu)失效概率或可靠度，通常需要顯式的功能函數(shù)，因此基于多種多項式近似的響應(yīng)面法得到了廣泛應(yīng)用，其中應(yīng)用最多也最可靠的是Bucher不含交叉項的響應(yīng)面形式：

其中，H為試驗點系數(shù)矩陣。

在可靠性分析中，基本隨機(jī)變量在每個點的分布概率均不相同，為保證在真實失效曲面附近和變量分布概率更高的區(qū)域擬合精度更高，本文采用加權(quán)最小二乘法進(jìn)行擬合，權(quán)重矩陣和每個實驗點的權(quán)重可以分別表示為：

其中，wi為權(quán)矩陣W中的對角元素。

對于工程中大多數(shù)有限元問題，當(dāng)非線性程度較低時，在各基本隨機(jī)變量xˉi±3δxi范圍內(nèi)擬合的超曲面精度都是可以接受的，因此響應(yīng)面法也在ANSYS可靠性分析模塊等商業(yè)軟件中得到廣泛應(yīng)用。

對于大多數(shù)結(jié)構(gòu)可靠性問題，基本隨機(jī)變量均值點都不會在真實失效曲面上，因此在均值點附近擬合得到的響應(yīng)面都有較大的誤差，尤其是可靠度較高、功能函數(shù)非線性較高的可靠性問題。

經(jīng)典響應(yīng)面法認(rèn)為，第一次圍繞均值點選取實驗點擬合出響應(yīng)面后，再采用式（6）進(jìn)行迭代，可保證算出的實驗點中心近似落在真實失效面上。

式中：μ為均值點；xD為設(shè)計驗算點。

這種方法盡管可以得到相對較高的計算精度，但對于復(fù)雜結(jié)構(gòu)問題，迭代過程難以實現(xiàn)且效率較低。

2 基于設(shè)計驗算點的可靠性方法

由于基本變量平均值不在極限狀態(tài)曲面上，隨著均值點到極限狀態(tài)曲面距離的增大，展開后的線性極限狀態(tài)面可能較大程度地偏離原來的極限狀態(tài)曲面。對于這一問題，采用驗算點一次2階矩法可以解決，即在極限狀態(tài)曲面上選取一個最大可能失效點P*，并在該點處進(jìn)行Taylor展開。

首先，設(shè)計驗算點需要滿足狀態(tài)函數(shù)：

其中，μxi和σxi分別為基本變量xi的均值和標(biāo)準(zhǔn)差。

定義可靠性靈敏度系數(shù)αi：

1階可靠性方法中假定結(jié)構(gòu)功能函數(shù)的分布為正態(tài)分布，可利用功能函數(shù)的均值和標(biāo)準(zhǔn)差計算得到可靠度指標(biāo)，即2階矩可靠度指標(biāo)。然而對于工程中大多數(shù)結(jié)構(gòu)可靠性問題，基本隨機(jī)變量和功能函數(shù)可能并不服從正態(tài)分布，需要考慮更高階分布參數(shù)，在統(tǒng)計學(xué)中認(rèn)為前4階中心矩即可描述一個隨機(jī)變量的大部分分布信息。

在高階矩可靠性方法中，可以采用Edgeworth級數(shù)法[5]和點估計法[6-8]得到功能函數(shù)前4階矩μg、δg、α3g、α4g，一種可行的功能函數(shù)標(biāo)準(zhǔn)正態(tài)化變換方法是Ono和Idota提出的高階矩標(biāo)準(zhǔn)化方法（HOMST），即假定存在多項式變換：

其中，aj是待定系數(shù)，通過計算Sx（x）的有限階中心矩，分別得到等于標(biāo)準(zhǔn)正態(tài)分布u的對應(yīng)階中心矩。當(dāng)考慮功能函數(shù)4階矩時，將功能函數(shù)的標(biāo)準(zhǔn)化變量通過如下關(guān)系式表示為標(biāo)準(zhǔn)正態(tài)變量：

在本文結(jié)構(gòu)可靠性方法中，首先根據(jù)式（4）擬合得到均值點附近的響應(yīng)面功能函數(shù)（如式（2）），再根據(jù)基于設(shè)計驗算點的1階可靠性方法得到該功能函數(shù)下的驗算點，最后，以驗算點為中心再次進(jìn)行響應(yīng)面系數(shù)擬合，得到本文提出的基于驗算點擬合的響應(yīng)面功能函數(shù)。盡管該驗算點可能不在真實失效曲面上，但相對于均值點已經(jīng)有較高的精度提升。對于非正態(tài)分布假設(shè)，本文采用點估計法計算功能函數(shù)前4階中心矩。

3 仿真算例

對于核反應(yīng)堆中控制棒驅(qū)動機(jī)構(gòu)的承壓殼體結(jié)構(gòu)，在ANSYS中建立軸對稱模型，其中局部結(jié)構(gòu)如圖1所示。在運(yùn)行過程中，承壓殼體主要承受內(nèi)壓P和步躍沖擊載荷F的作用。本文對承壓殼體的可靠性分析分別考慮承壓殼體與頂蓋貫穿件異種金屬焊縫處的強(qiáng)度失效（評定截面ASN1），鉤爪殼體中部強(qiáng)度失效（評定截面ASN2），以及棒行程殼體和鉤爪殼體Canopy焊縫處的強(qiáng)度失效（評定截面ASN3）。

圖1 CRDM承壓殼體有限元模型和評定截面

考慮內(nèi)壓P、步躍沖擊載荷F及材料設(shè)計應(yīng)力強(qiáng)度Sm三個變量為隨機(jī)變量，則3個失效模式一次薄膜應(yīng)力可靠度的功能函數(shù)可以表示為

其中，x= ［SmP F］T為基本隨機(jī)變量，假定所有隨機(jī)變量均服從正態(tài)分布，分布參數(shù)如表1所示。

表1 各基本隨機(jī)變量分布參數(shù)

用本文方法得到的近似驗算點和驗算點處的真實功能函數(shù)值如表2所示。

表2 近似驗算點及真實功能函數(shù)值

從表2可以看出，本文在均值點響應(yīng)面功能函數(shù)基礎(chǔ)上得到的近似驗算點處，真實功能函數(shù)值都接近于0，可認(rèn)為都處于真實失效面上。

由于響應(yīng)面功能函數(shù)有顯式形式，采用Monte-Carlo仿真也能夠有極高的效率。分別基于有限元模型進(jìn)行5000次仿真，基于響應(yīng)面功能函數(shù)進(jìn)行107次仿真，得到3個失效截面可靠度，如表3所示。

表3 基于響應(yīng)面和有限元模型的MC仿真結(jié)果

表3中，RSM_E表示基于均值點擬合得到的響應(yīng)面結(jié)果，RSM_D表示本文基于近似驗算點擬合得到的響應(yīng)面結(jié)果，ANSYS表示在ANSYS中PDS模塊得到的收斂可靠度，Error表示本文結(jié)果與ANSYS結(jié)果的相對誤差。

可以看出，本文提出的基于近似驗算點擬合的響應(yīng)面與真實失效曲面更為接近，在計算可靠度時有更高的仿真精度。需要指出的是，由于有限元模型復(fù)雜，在ANSYS中用重要采樣法進(jìn)行5000次仿真需要約5 h的計算時間，而用近似響應(yīng)面功能函數(shù)進(jìn)行107次抽樣仿真只需要數(shù)分鐘時間，在精度相近的情況下計算效率提升明顯。

當(dāng)基本隨機(jī)變量為非正態(tài)分布時，假定Sm服從Weibull分布（115.92,12.15），內(nèi)壓服從對數(shù)正態(tài)分布（2.83,0.1），步躍載荷服從正態(tài)分布（60 000,6000），用式(16)計算可靠度指標(biāo)得到的4階矩可靠度求解得到的截面可靠度如表4所示。

表4 4階矩可靠度與MC法的對比

表4中，R_SM表示2階矩可靠度指標(biāo)計算得到的結(jié)構(gòu)可靠度，R_FM表示4階矩可靠度指標(biāo)計算得到的結(jié)構(gòu)可靠度，Error表示4階矩可靠度與ANSYS中進(jìn)行5000次Monte-Carlo仿真得到的可靠度間的相對誤差。

從表4中可以看出，當(dāng)基本隨機(jī)變量為非正態(tài)分布時，用4階矩可靠度指標(biāo)計算得到的結(jié)果可靠度有更高的精度。

4 結(jié)論

本文針對結(jié)構(gòu)復(fù)雜的可靠性功能函數(shù)，基于均值響應(yīng)面法和設(shè)計驗算點法，得到近似驗算點，再通過近似驗算點局部區(qū)域的響應(yīng)面擬合得到精度更高的響應(yīng)面功能函數(shù)。對于非正態(tài)分布的結(jié)構(gòu)可靠性問題，基于4階矩可靠度指標(biāo)法計算的結(jié)構(gòu)可靠度也有更高的精度，相對于傳統(tǒng)Monte-Carlo方法，近似的結(jié)構(gòu)失效功能函數(shù)也有更高的計算效率。對于工程中結(jié)構(gòu)復(fù)雜且可靠度要求較高的結(jié)構(gòu)可靠性問題，本文方法有很好的適用性和參考意義。