邱 濤，王增臣，游令非

(1.中國(guó)電子科技集團(tuán)公司 第二十八研究所，江蘇 南京 210000；2.北京航空航天大學(xué) 可靠性與系統(tǒng)工程學(xué)院，北京 100191；3.北京航空航天大學(xué) 可靠性與環(huán)境工程技術(shù)國(guó)防科技重點(diǎn)實(shí)驗(yàn)室，北京 100191)

0 引言

在機(jī)械結(jié)構(gòu)(航空、航天、武器)服役過(guò)程中，工作環(huán)境愈加惡劣，不確定影響因素日趨增加，建立顯式的結(jié)構(gòu)可靠性模型比較困難。對(duì)于復(fù)雜機(jī)械結(jié)構(gòu)，認(rèn)知局限和客觀隨機(jī)性的存在使其在研制過(guò)程中充滿不確定性[1]。當(dāng)可靠性模型的輸入變量信息較少時(shí)，很難得到輸入變量的精確概率分布，可以采用分布未知、邊界可知的區(qū)間變量進(jìn)行描述[2-3]。

對(duì)于概率不確定模型，有一次二階矩法(First Order Second Moment, FORM)[4-5]、HLRF(Hasofer-Lind and Rackwitz-Fiessler)[6-7]和aHLRF[8]。對(duì)于區(qū)間不確定模型，唐忠等[9]基于偏彈性理論，構(gòu)造了一種含區(qū)間變量的二階多項(xiàng)式功能函數(shù)；Wang等[10]考慮不確定信息對(duì)系統(tǒng)結(jié)構(gòu)安全的影響，提出一種基于非概率集合理論的區(qū)間可靠性分析方法；孫靜怡等[11]為抑制區(qū)間算法導(dǎo)致的區(qū)間擴(kuò)張問(wèn)題，提出一種區(qū)間仿射響應(yīng)面法。

對(duì)于概率—區(qū)間混合的可靠性模型，Elishakoff等[12]基于凸模型理論，提出一種解決混合可靠性問(wèn)題的分析方法；Luo等[13]提出一種迭代算法，用于尋找不確定性空間中的最壞情況點(diǎn)和最可能故障點(diǎn)，以計(jì)算可靠性指標(biāo)；Qiu等[14]將區(qū)間方法引入傳統(tǒng)可靠性理論，獲得系統(tǒng)失效概率區(qū)間；Du[15-16]針對(duì)概率—區(qū)間混合的可靠性問(wèn)題，提出一種概率分析與區(qū)間優(yōu)化嵌套處理的雙層優(yōu)化模型，但在工程實(shí)際中功能函數(shù)大多為隱式，采用該方法會(huì)大大降低計(jì)算效率；姜潮等[17-18]基于Bucher建立線性響應(yīng)面，提高了混合可靠性分析的計(jì)算效率；Alibrandi等[19]基于一次二階矩法，通過(guò)較好地近似設(shè)計(jì)點(diǎn)來(lái)提高計(jì)算效率，但對(duì)于非線性程度較高的功能函數(shù)，計(jì)算精度會(huì)降低；Xie等[20]在考慮區(qū)間變量相關(guān)的情況下，用單步可靠性序列迭代算法求解雙層優(yōu)化模型；Gao等[21]提出一種統(tǒng)一區(qū)間隨機(jī)可靠性抽樣方法，用于評(píng)估涉及多個(gè)不精確隨機(jī)場(chǎng)的工程結(jié)構(gòu)的安全性。根據(jù)以上研究表明，當(dāng)功能函數(shù)為隱式表達(dá)式時(shí)，存在計(jì)算效率降低的問(wèn)題；當(dāng)功能函數(shù)非線性程度較高時(shí)，會(huì)導(dǎo)致計(jì)算結(jié)果不收斂或使計(jì)算精度下降。

針對(duì)以上問(wèn)題，基于響應(yīng)面法和一次二階矩法，本文提出一種高效的混合可靠性分析算法。當(dāng)功能函數(shù)為隱式表達(dá)式時(shí)，為減少原始極限狀態(tài)函數(shù)的調(diào)用次數(shù)，提出一種基于二次插值抽樣的響應(yīng)面法；當(dāng)功能函數(shù)非線性程度較高時(shí)，針對(duì)計(jì)算結(jié)果不收斂或計(jì)算精度下降的問(wèn)題，提出一種基于二參數(shù)尋優(yōu)設(shè)計(jì)點(diǎn)的概率分析算法。

1 概率—區(qū)間混合可靠性模型

結(jié)構(gòu)的可靠性受材料屬性、幾何尺寸、環(huán)境載荷等多種不確定因素影響，其功能函數(shù)可表示為

g(Z)=g(X,Y)。

(1)

式中Z=(X,Y)，X=[X1,X2,…,Xn]T表示隨機(jī)變量，Y=[Y1,Y2,…,Ym]T表示區(qū)間變量。當(dāng)g(Z)>0時(shí)，結(jié)構(gòu)安全；當(dāng)g(Z)<0時(shí)，結(jié)構(gòu)發(fā)生失效；當(dāng)g(Z)=0時(shí)，結(jié)構(gòu)處于極限狀態(tài)。

根據(jù)圖1可知，失效概率Pf為一區(qū)間范圍，可表示為[15]：

(2)

(3)

2 基于二次插值抽樣的響應(yīng)面法

由于功能函數(shù)含有區(qū)間變量，不能直接進(jìn)行可靠度求解，首先采用高維模型表示(High Dimensional Model Representation, HDMR)模型將混合可靠性模型解耦為單層可靠性模型。在擬合單層可靠性模型時(shí)，為減少原始極限狀態(tài)函數(shù)的調(diào)用次數(shù)，提出二次插值抽樣法(Quadratic Interpolation Sampling, QIS)進(jìn)行抽樣。

2.1 單層可靠性模型

HDMR是一種近似函數(shù)模型，低階的HDMR模型同樣具有較高的擬合精度[22]。

基于高維模型表示方法，g(Z)可表示為

(4)

式中：l=n+m；g0為0階分量函數(shù)；gr(Zr)為1階分量函數(shù)；grs(Zr,Zs)為2階分量函數(shù)；grst(Zr,Zs,Zt)為3階分量函數(shù)；g12…l(Z1,Z2,…,Zl)為所有殘余的耦合輸入變量對(duì)響應(yīng)值的影響。

因?yàn)镠DMR展式中的1階分量函數(shù)gr(Zr)是泰勒展式中僅含變量Zr的集合，所以1階分量函數(shù)gr(Zr)可以為非線性。通過(guò)引入HDMR展式的1階分量函數(shù)，將功能函數(shù)解耦為隨機(jī)和區(qū)間變量相互分離的單元函數(shù)，近似表達(dá)為

(5)

式中：c為擬合近似表達(dá)式的參考點(diǎn)；g(Xi)為僅包含隨機(jī)變量Xi的單元函數(shù)；g(Yj)為僅包含區(qū)間變量Yj的單元函數(shù)。

2.2 二次插值抽樣法

在抽取樣本點(diǎn)時(shí)，提出了QIS，第1次插值用于確定參考點(diǎn)，使參考點(diǎn)更加靠近功能函數(shù)；第2次插值用于確定擬合樣本點(diǎn)，降低抽樣系數(shù)對(duì)樣本點(diǎn)位置的影響。

(6)

s=1,2,…,l；

(7)

s=l+1,l+2,…,2l。

(8)

雖然CS能在參考點(diǎn)附近很好地?cái)M合響應(yīng)面，但是抽取的樣本點(diǎn)與抽樣系數(shù)f的關(guān)系較大，導(dǎo)致計(jì)算結(jié)果對(duì)f較為敏感。如果f過(guò)大，則樣本點(diǎn)距功能函數(shù)較遠(yuǎn)，導(dǎo)致擬合精度下降；如果f過(guò)小，則樣本點(diǎn)較為密集，導(dǎo)致擬合矩陣奇異。

(9)

(10)

r=l+1,l+2,…,2l。

(11)

3 可靠性分析

因混合可靠性模型被解耦為隨機(jī)和區(qū)間變量分離的單層可靠性計(jì)算模型，可對(duì)含隨機(jī)變量的單元函數(shù)進(jìn)行概率分析，對(duì)含區(qū)間變量的單元函數(shù)進(jìn)行區(qū)間分析。

3.1 基于二參數(shù)尋優(yōu)設(shè)計(jì)點(diǎn)的概率分析

在概率分析時(shí)，基于FORM，通過(guò)增加兩個(gè)參數(shù)η和λ來(lái)調(diào)整每一次迭代的搜索步長(zhǎng)和搜索方向，提出基于二參數(shù)尋優(yōu)設(shè)計(jì)點(diǎn)的一次二階矩法tHLRF。

(Uk+1-Uk)。

(12)

式中:Uk為第k次迭代的設(shè)計(jì)驗(yàn)算點(diǎn)(Most Probabile Point, MPP)；Uk+1為第k+1次迭代的MPP；G(Uk)為功能函數(shù)G(U)在Uk處的梯度；為功能函數(shù)G(U)在Uk處的一階泰勒展開(kāi)。

(13)

=ηG(Uk)。

(14)

式中η∈[0,1)。

考慮取步長(zhǎng)λ為有限值，以改變尋優(yōu)設(shè)計(jì)點(diǎn)的搜索方向，可得輔助點(diǎn)Hk+1，

(15)

式中λ為從Uk沿-G(Uk)方向移動(dòng)的步長(zhǎng)。

在平面G=ηG(Uk)上，Uk+1的搜索方向?yàn)樵c(diǎn)到輔助點(diǎn)Hk+1的方向，其方向余弦矢量

(16)

Uk+1=βk+1αk+1。

(17)

式中βk+1為第k+1次迭代的可靠度指標(biāo)。

將式(17)代入式(14)，求解得可靠度指標(biāo)

(18)

根據(jù)求得的可靠度指標(biāo)βk+1，由式(17)得第k+1次迭代的MPP的Uk+1。

在選擇參數(shù)η時(shí)，初始值設(shè)定為較小的值，可取η=0.1，在迭代過(guò)程中，設(shè)定ηk+1=0.1+ζηk以自適應(yīng)更新參數(shù)η，ζ為η的調(diào)整系數(shù)，一般取在0.01～0.05之間；在選擇參數(shù)λ時(shí)，初始值設(shè)定為較大的值，可取λ=50，在迭代過(guò)程中，設(shè)定λk+1=λk/ξ以自適應(yīng)更新參數(shù)λ，ξ為λ的參數(shù)調(diào)整系數(shù)，一般取在1.2～1.5之間。

3.2 區(qū)間分析

在第k+1次迭代中，含區(qū)間變量的單元函數(shù)為g(Y)，區(qū)間分析可表示為式(19)的區(qū)間優(yōu)化問(wèn)題：

s.t.

YL≤Y≤YR。

(19)

式中：Yk+1為第k+1次迭代的區(qū)間最優(yōu)點(diǎn)；YL為區(qū)間的下界矢量；YR為區(qū)間的上界矢量。

通過(guò)KKT條件判斷當(dāng)前迭代點(diǎn)是否為最優(yōu)點(diǎn)：

(20)

(21)

式中：j=1,2,…,m表示m維區(qū)間變量的第j維；Yk為第k次迭代的區(qū)間最優(yōu)點(diǎn)。

若滿足KKT條件，則Yk+1=Yk；若不滿足，則進(jìn)行區(qū)間優(yōu)化分析[23]。設(shè)當(dāng)前為第r次迭代，迭代點(diǎn)為Yr，將此區(qū)間優(yōu)化問(wèn)題轉(zhuǎn)化為二次規(guī)劃問(wèn)題：

s.t.

YL≤Y≤YR

。

(22)

式中Hs為Hessian矩陣。

式(22)為區(qū)間約束的二次規(guī)劃問(wèn)題，采用梯度投影法[24]可將迭代方向投影到區(qū)間可行域上，以高效計(jì)算含區(qū)間變量單元函數(shù)的最優(yōu)點(diǎn)。若新迭代點(diǎn)Yr+1滿足KKT條件，則區(qū)間迭代停止，Yk+1=Yr+1。

(23)

綜上所述，論文所提可靠性分析算法的流程圖如圖7所示。

4 案例分析

在兩個(gè)案例分析中，收斂精度ε1=ε2=1×10-6，在概率分析時(shí)，初始參數(shù)η=0.1，λ=30，參數(shù)調(diào)整系數(shù)ζ=0.05,ξ=1.2；N為可靠性分析算法的迭代循環(huán)次數(shù)，t為算法計(jì)算所耗時(shí)間。

4.1 懸臂梁案例

如圖8和圖9所示，懸臂梁長(zhǎng)度為L(zhǎng)，橫截面高度為h，寬度為w，材料的彈性模量為220 GPa，在其末端作用的水平和垂直載荷分別為Px，Py。懸臂梁末端的許用位移不超過(guò)δ=56 mm，考慮Px和Py為區(qū)間變量，L，w，h服從正態(tài)分布，其不確定特征參數(shù)如表1所示。

表1 懸臂梁的不確定特征

因?yàn)閼冶哿耗┒说脑S用位移不超過(guò)δ=56 mm，所以功能函數(shù)表示為

為驗(yàn)證所提算法的精度，將每個(gè)區(qū)間等分為10份，采用蒙特卡洛采樣法(Monte Carlo Sampling，MCS)在區(qū)間變量的組合值下對(duì)隨機(jī)變量抽取1×108次，調(diào)用功能函數(shù)1×1010次。根據(jù)表2可知，相比MCS，所提算法的相對(duì)誤差為4.2%，驗(yàn)證了該算法的精度。

表2 懸臂梁最壞情況下的失效概率

根據(jù)圖10不同可靠性分析算法收斂過(guò)程的對(duì)比和表2可知，基于CS+tHLRF的混合可靠性分析算法在迭代循環(huán)15次后達(dá)到收斂，可靠度計(jì)算時(shí)間為629 s；基于QIS+tHLRF的混合可靠性分析算法在迭代循環(huán)12次后達(dá)到收斂，可靠度計(jì)算時(shí)間為545 s。相比經(jīng)典抽樣法，驗(yàn)證了論文所提QIS具有更高的效率。

基于QIS+FORM分析算法的可靠度計(jì)算結(jié)果不收斂；基于QIS+iHLRF[25]分析算法迭代13次后收斂，可靠度計(jì)算時(shí)間為582 s。對(duì)比基于FORM的混合可靠性分析算法，當(dāng)功能函數(shù)的非線性程度較高時(shí)，所提基于tHLRF的混合可靠性分析算法同樣滿足收斂要求，相比基于iHLRF的混合可靠性分析算法，驗(yàn)證了所提基于tHLRF的混合可靠性分析算法具有更高的計(jì)算效率。文獻(xiàn)[15]為該類(lèi)方法典型的研究成果，主要用于解決功能函數(shù)為顯示表達(dá)式的問(wèn)題，可靠度計(jì)算時(shí)間為563 s；文獻(xiàn)[20]為該類(lèi)方法中較新的研究成果，可靠度計(jì)算時(shí)間為632 s。相比文獻(xiàn)[15,20]，同樣驗(yàn)證了本文所提算法具有較高的計(jì)算效率。

根據(jù)以上對(duì)比結(jié)果可知，在同樣收斂精度要求下，QIS較經(jīng)典抽樣法具有更高的計(jì)算效率；當(dāng)功能函數(shù)非線性程度較高時(shí)，tHLRF能解決收斂不穩(wěn)定的問(wèn)題，而且具有較高的計(jì)算效率。與該領(lǐng)域中的研究成果對(duì)比，說(shuō)明本文所提算法在保證計(jì)算精度的條件下具有較高的計(jì)算效率。

4.2 工程案例

表3 壓氣機(jī)輪盤(pán)的不確定特性

案例的工況條件為低循環(huán)疲勞試驗(yàn)，其轉(zhuǎn)速為9 550 r/min，30 s為一個(gè)周期，溫度為230℃恒溫。輪盤(pán)經(jīng)過(guò)105 300次試驗(yàn)循環(huán)，耗時(shí)約877 h，在銷(xiāo)釘孔處產(chǎn)生初始裂紋[26]，裂紋如圖11所示。

由于壓氣機(jī)輪盤(pán)沿圓周均勻?qū)ΨQ分布37個(gè)銷(xiāo)釘孔，現(xiàn)選取輪盤(pán)的1/37作為計(jì)算模型(其有限元模型如圖12)，并對(duì)輪盤(pán)銷(xiāo)釘孔周?chē)木W(wǎng)格進(jìn)行了加密。

通過(guò)上述模型對(duì)輪盤(pán)進(jìn)行有限元分析，應(yīng)變分析結(jié)果如圖13所示，同樣得出壓氣機(jī)輪盤(pán)的關(guān)鍵部位為銷(xiāo)釘孔?，F(xiàn)對(duì)該輪盤(pán)建立疲勞可靠性模型，計(jì)算該壓氣機(jī)輪盤(pán)在壽命為8 000周下的失效概率。

根據(jù)Coffin-Manson公式，結(jié)合Miner法則[27]建立輪盤(pán)的低循環(huán)疲勞可靠性模型。

裂紋形成壽命NS采用Coffin-Manson公式求解：

(24)

由Miner法則，若構(gòu)件在應(yīng)力載荷水平Sj作用下經(jīng)受nj次循環(huán)的損傷為T(mén)j=nj/Nj，則當(dāng)構(gòu)件承受k個(gè)載荷水平作用時(shí)，所受累計(jì)損傷

(25)

結(jié)構(gòu)危險(xiǎn)點(diǎn)的功能函數(shù)

g(Z)=g(X,Y)=a-T(Z)。

(26)

式中：Z=(X,Y)，X=[X1,X2,…,Xn]T為影響疲勞壽命的隨機(jī)變量，Y=[Y1,Y2,…,Ym]T為影響疲勞壽命的區(qū)間變量；a為損傷強(qiáng)度，a=1。

根據(jù)以上疲勞可靠性模型，采用不同可靠性分析算法求解輪盤(pán)的最大失效概率，收斂過(guò)程與計(jì)算結(jié)果分別如圖14和表4所示。

表4 輪盤(pán)最壞情況下的失效概率

在該工程對(duì)象中，無(wú)法得到壓氣機(jī)輪盤(pán)的顯式功能函數(shù)，首先根據(jù)輪盤(pán)失效機(jī)理建立疲勞可靠性模型；其次采用基于二次插值抽樣的響應(yīng)面法擬合功能函數(shù)；最后利用所提可靠性分析算法求解失效概率。根據(jù)表4可知，相比MCS，所提算法計(jì)算結(jié)果的相對(duì)誤差為4.8%，既驗(yàn)證了所提算法的精度，又說(shuō)明其可以解決工程中難以建立顯式功能函數(shù)的結(jié)構(gòu)可靠性問(wèn)題。

由圖14和表4可知MCS計(jì)算用時(shí)22 683 s，文獻(xiàn)[15]主要針對(duì)功能函數(shù)為顯式表達(dá)式的可靠度計(jì)算方法，其迭代37次收斂，可靠度計(jì)算用時(shí)13 268 s；文獻(xiàn)[18]為同類(lèi)方法中功能函數(shù)為隱式表達(dá)式時(shí)的典型研究成果，其迭代22次收斂，可靠度計(jì)算用時(shí)7 562 s；文獻(xiàn)[20]為同類(lèi)方法中較新的研究成果，其迭代19次收斂，可靠度計(jì)算用時(shí)6 923 s；本文所提算法迭代17次收斂，可靠度計(jì)算用時(shí)6 586 s。通過(guò)與文獻(xiàn)[15]對(duì)比，驗(yàn)證了結(jié)構(gòu)的功能函數(shù)為隱式表達(dá)式時(shí)，所提算法可以減少原始極限狀態(tài)函數(shù)的調(diào)用次數(shù)；通過(guò)與文獻(xiàn)[18,20]對(duì)比，驗(yàn)證了所提算法具有較高的計(jì)算效率。在大型工程實(shí)例中，功能函數(shù)多為隱式表達(dá)式，計(jì)算耗時(shí)較長(zhǎng)，而采用所提算法具有較高的計(jì)算效率，說(shuō)明所提算法具有較好的工程應(yīng)用價(jià)值。

5 結(jié)束語(yǔ)

針對(duì)概率—區(qū)間混合的可靠性問(wèn)題，本文提出一種高效的可靠性分析算法，為解決混合結(jié)構(gòu)可靠性問(wèn)題提供了新的參考。具體貢獻(xiàn)如下：

(1)提出QIS，使樣本點(diǎn)盡可能落在真實(shí)功能函數(shù)附近，降低了插值系數(shù)f對(duì)樣本點(diǎn)位置的影響，減少了對(duì)原始極限狀態(tài)函數(shù)的調(diào)用次數(shù)，提高了計(jì)算效率。

(2)提出基于二參數(shù)尋優(yōu)設(shè)計(jì)點(diǎn)的概率分析算法，在尋優(yōu)MPP時(shí)，引入兩個(gè)調(diào)整參數(shù)和兩個(gè)自適應(yīng)參數(shù)調(diào)整系數(shù)，分別控制搜索方向與搜索步長(zhǎng)，解決了功能函數(shù)非線性程度較高時(shí)計(jì)算結(jié)果收斂不穩(wěn)定的問(wèn)題，并提高了計(jì)算效率。

(3)本文提出的結(jié)構(gòu)可靠性分析算法適用于通用的隨機(jī)和區(qū)間變量共存的結(jié)構(gòu)可靠性問(wèn)題，對(duì)實(shí)際工程應(yīng)用具有一定的指導(dǎo)意義。

本文引入的兩個(gè)自適應(yīng)調(diào)整系數(shù)為在一定的范圍內(nèi)取值，后續(xù)研究將考慮采用優(yōu)化函數(shù)的形式，通過(guò)求解自適應(yīng)調(diào)整系數(shù)來(lái)提高可靠度的計(jì)算效率。