粒子群優(yōu)化與Kriging 模型相結(jié)合的結(jié)構(gòu)非概率可靠性分析*

喬心州， 陳永婧， 劉 鵬， 方秀榮

（西安科技大學(xué) 機(jī)械工程學(xué)院，西安 710054）

引 言

工程結(jié)構(gòu)往往存在大量不確定因素，可靠性方法是處理上述不確定性因素的有效手段之一.傳統(tǒng)的概率可靠性方法已得到廣泛應(yīng)用，然而由于條件限制，許多情況下難以獲得用于確定概率分布的足夠樣本信息.近年來，能夠處理有限樣本信息的非概率可靠性方法受到廣泛關(guān)注.

Ben-Haim[1]首次提出了非概率可靠性的概念.Elishakoff 等[2]針對此概念提出了一種基于非概率安全系數(shù)的可靠性度量，該度量不是一特定數(shù)值而是一區(qū)間值.與概率可靠性指標(biāo)類似，郭書祥等[3]提出了一種基于區(qū)間模型的非概率可靠性指標(biāo)，并闡述了該指標(biāo)的物理和幾何意義.Jiang 等[4]發(fā)展了針對該可靠性指標(biāo)的半解析法.曹鴻鈞和段寶巖[5]提出了一種基于橢球模型的非概率可靠性指標(biāo)，該指標(biāo)采用標(biāo)準(zhǔn)空間中原點(diǎn)到失效面的最短歐式距離度量可靠性.Jiang 等[6]針對上述可靠性指標(biāo)提出中心點(diǎn)法和設(shè)計(jì)點(diǎn)法兩種求解方法.Meng 等[7]提出了一種基于p范數(shù)的超參數(shù)凸模型可靠性指標(biāo).Meng 等[8]進(jìn)一步提出了一種基于指數(shù)凸模型的可靠性指標(biāo).借鑒概率可靠度的概念，Wang 等[9]、Jiang 等[10]進(jìn)一步提出了區(qū)間模型和橢球模型的非概率可靠度概念，得到了不確定變量域落入可靠域的體積與整個不確定變量域的體積的比值度量可靠性.Jiang 等[10]根據(jù)此概念進(jìn)一步提出了橢球模型非概率可靠度的一階近似法和二階近似法.

上述可靠性分析方法均是針對顯式結(jié)構(gòu)功能函數(shù)，在復(fù)雜工程結(jié)構(gòu)中通常面臨的是隱式結(jié)構(gòu)功能函數(shù)，對其進(jìn)行可靠性分析通常需要采用計(jì)算成本較高、樣本需求量大的有限元分析或Monte-Carlo 模擬.鑒于此，采用計(jì)算量小、樣本需求少的代理模型來代替原分析模型對結(jié)構(gòu)進(jìn)行可靠性分析成為了一種有效手段.在基于代理模型非概率可靠性分析方面，陳江義等[11]利用響應(yīng)面法對橢球模型與區(qū)間模型進(jìn)行非概率可靠性分析.Bai 等[12]提出了一種基于不含交叉項(xiàng)二次響應(yīng)面法的非概率可靠性分析方法.馬超等[13]用支持向量機(jī)法擬合結(jié)構(gòu)功能函數(shù)，并將該方法應(yīng)用于實(shí)際的飛機(jī)機(jī)翼可靠性分析中.潘林鋒等[14]用Kriging 模型擬合回歸極限狀態(tài)方程，并用HL-RF 修正法計(jì)算了非概率可靠性指標(biāo).鄭嚴(yán)[15]利用支持向量機(jī)和三種改進(jìn)粒子群算法進(jìn)行了結(jié)構(gòu)可靠性分析及優(yōu)化設(shè)計(jì)，主要解決了隱式結(jié)構(gòu)概率及非概率可靠性分析問題.但Kriging 模型的搜索方式對初始值具有依賴性，初始值的設(shè)置誤差會導(dǎo)致后續(xù)模型建立時預(yù)測結(jié)果陷入局部最優(yōu)的情況，且Kriging 模型采用模式搜索法進(jìn)行單點(diǎn)搜索，搜索路徑單一.

針對Kriging 模型的上述不足，本文提出了一種粒子群優(yōu)化(particle swarm optimization，PSO)與Kriging 模型相結(jié)合的非概率可靠性分析方法.本文的結(jié)構(gòu)安排如下：第1 節(jié)對基于橢球模型的非概率可靠性指標(biāo)進(jìn)行概述；第2 節(jié)給出了本文所提方法及其實(shí)施流程；第3 節(jié)通過三個算例驗(yàn)證了本文方法的有效性；第4 節(jié)對本文進(jìn)行總結(jié).

1 基于橢球模型的非概率可靠性指標(biāo)

考慮某工程結(jié)構(gòu)功能函數(shù)：

為方便對結(jié)構(gòu)可靠性進(jìn)行定義，需首先將不確定參數(shù)X線性變換為其標(biāo)準(zhǔn)化形式δ，相應(yīng)的結(jié)構(gòu)功能函數(shù)和不確定域分別變換為G(δ)和 一個中心在原點(diǎn)的單位超球δTδ≤1.圖1 給出了二維空間的線性變換.通過上述處理，可有效避免由于工程實(shí)際問題中不確定參數(shù)的數(shù)量級不同導(dǎo)致病態(tài)矩陣的產(chǎn)生，同時變量的無量綱化也便于結(jié)構(gòu)可靠性的定義與比較.

圖1 線性變換與可靠性指標(biāo)Fig. 1 The linear transform and the reliability index

此時，結(jié)構(gòu)的非概率可靠性指標(biāo) η可定義為標(biāo)準(zhǔn)化空間內(nèi)原點(diǎn)到失效面G(δ)=0的最短距離，可通過如下優(yōu)化問題求解：

上述可靠性指標(biāo)與傳統(tǒng)的概率可靠性指標(biāo)幾何意義相同，可采用HL-RF 法進(jìn)行求解.

2 PSO-Kriging 方法的非概率可靠性分析

考慮到許多工程實(shí)際問題中結(jié)構(gòu)功能函數(shù)均為隱式的，以下發(fā)展一種將上述非概率可靠性指標(biāo)與PSOKriging 方法相結(jié)合的非概率可靠性分析方法.該方法主要由以下三個部分組成：① 構(gòu)建Kriging 模型；② 運(yùn)用PSO 方法對模型相關(guān)參數(shù)進(jìn)行尋優(yōu)；③ 基于代理模型求解可靠性指標(biāo).

對于代理模型的建立，首先需要選取合理分布的樣本點(diǎn).樣本點(diǎn)的選取對后續(xù)建模的精度有一定影響，因此對于樣本采樣的最佳策略應(yīng)在符合采樣精度的前提下，采樣時間最短、采樣點(diǎn)應(yīng)均勻分布在整個待測區(qū)域.Sobol 序列抽樣方法[17]無論樣本點(diǎn)N為多少，均可在面域內(nèi)取得均勻分布的點(diǎn)，且不會出現(xiàn)局部團(tuán)簇現(xiàn)象，故本文采用Sobol 序列抽樣方法進(jìn)行抽樣.

2.1 Kriging 模型的構(gòu)建

Kriging 模型將需代理求解的未知函數(shù)看作是一個隨機(jī)過程，其表達(dá)的隱函數(shù)為

上述最優(yōu)相關(guān)參數(shù)的確定為構(gòu)建最優(yōu)Kriging 模型奠定了基礎(chǔ)，進(jìn)而可求解相應(yīng)的可靠性指標(biāo).

2.2 PSO 方法求解模型相關(guān)參數(shù)

PSO 方法采用多點(diǎn)并行搜索，利用群體中的各個粒子對信息進(jìn)行共享，每個粒子根據(jù)個體極值和群體極值的變化方向判斷自身的搜索方向，從而使得整個群體的運(yùn)動在求解空間中產(chǎn)生從無序到有序的演變過程，無論初始值的位置如何，都可大概率尋找到最優(yōu)解.本文采用PSO 對相關(guān)參數(shù)θ 進(jìn)行尋優(yōu).

在尋優(yōu)過程中，以式(9)作為適應(yīng)度函數(shù)，以Kriging 模型的n個樣本點(diǎn)作為粒子，其中第i個粒子位置為xi，根據(jù)適應(yīng)度函數(shù)，將第i個粒子的個體極值pbest記 為pi， 群體最優(yōu)位置極值gi為 所有pbest的最優(yōu)解.每個粒子的移動速度記為vi，則粒子根據(jù)每次位置判斷下次更新的速度和位置為

式中，k為當(dāng)前迭代數(shù)， r and 為 ( 0,1)之 間的隨機(jī)數(shù)， ω為 慣性權(quán)重因子，更新的學(xué)習(xí)因子分別為c1，c2.當(dāng)粒子群搜索滿足式(9)時，最優(yōu) θ值即為該粒子的全局最優(yōu)位置gbest.此時構(gòu)建的Kriging 模型為最優(yōu)模型，可據(jù)此模型擬合出功能函數(shù)，利用式(4)計(jì)算設(shè)計(jì)點(diǎn)和可靠性指標(biāo)η.

由于代理模型是近似擬合，因此需要多次擬合結(jié)構(gòu)功能函數(shù)并求解可靠性指標(biāo)，直至滿足收斂條件|ηi?ηi?1|≤ε， ε為給定的允許值.

2.3 所提方法實(shí)施過程及流程

本文所提方法流程如圖2 所示.

圖2 算法流程圖Fig. 2 The flowchart for the proposed method

步驟1 利用Sobol 序列抽樣選取n個樣本點(diǎn)，根據(jù)隱式函數(shù)確定對應(yīng)的函數(shù)值.

步驟2 PSO 搜索最優(yōu)θ 值的過程：

1) 設(shè)置初始參數(shù).設(shè)置PSO 初始種群個數(shù)為N，設(shè)定最大迭代次數(shù)為Iger， 并將慣性權(quán)重設(shè)為 ω，學(xué)習(xí)因子為c1，c2，初始種群的位置以及粒子的位置參數(shù)及速度限制范圍.

2) 通過式(9)計(jì)算每個粒子的個體適應(yīng)度值，比較當(dāng)前粒子的適應(yīng)度與歷史最佳適應(yīng)度pbest，如果當(dāng)前位置的適應(yīng)值更高，更新粒子群的粒子的個體最優(yōu)適應(yīng)度pbest.

3) 將粒子當(dāng)前位置的適應(yīng)值與其全局最佳位置gbest對應(yīng)的適應(yīng)值進(jìn)行比對，如果當(dāng)前位置的適應(yīng)值更高，則用當(dāng)前位置更新全局最佳位置gbest.

4) 更新每個粒子的速度與位置.

5) 判斷是否滿足式(9)，若滿足，則得出最優(yōu)的θ；若不滿足，則繼續(xù)迭代.

步驟3 得出最優(yōu)的θ 值，構(gòu)建相應(yīng)的Kriging 模型.

步驟4 由預(yù)測出的Kriging 模型求解非概率可靠性指標(biāo)及設(shè)計(jì)點(diǎn)坐標(biāo).

步驟5 若可靠性指標(biāo)滿足收斂條件| ηi?ηi?1|≤ε，則停止；若不滿足，則在當(dāng)前設(shè)計(jì)點(diǎn)附近重新確定樣本范圍，返回步驟1.

3 算 例

式中，X1∈[0.0,1.0]，X2∈[0.0,1.0]. 在本算例中，樣本M=50， 參數(shù)c1，c2均 設(shè)置為0.3，慣性權(quán)重 ω均設(shè)為0.75，最大迭代次數(shù)為10 次，速度限制為V∈[?1,1]； 迭代收斂標(biāo)準(zhǔn)為| ηi?ηi?1|≤0.001.

表1 給出了相關(guān)系數(shù)為0.2 時本文方法的迭代過程；表2 列出了不同相關(guān)系數(shù)下，本文方法、Kriging 方法和HL-RF 法(參考值)的對比結(jié)果；圖3 給出了不同相關(guān)系數(shù)下兩種方法的相對誤差.

圖3 相對誤差對比圖（算例1）Fig. 3 Comparison of relative errors (example 1)

表1 相關(guān)系數(shù)為0.2 的迭代過程Table 1 The iterative process for a correlation coefficient of 0.2

表2 不同相關(guān)系數(shù)的可靠性指標(biāo)（算例1）Table 2 Reliability indexes for different correlation coefficients (example 1)

算例2 考慮圖4 所示的刮板輸送機(jī)GB/T 12718 鏈輪，在其鏈窩處受到一個120 kN 的切向力.在本算例中，齒弧半徑R1， 齒根圓弧半徑R2， 鏈窩部的平面圓弧半徑R3， 短齒根部的圓弧半徑R4以及平環(huán)中面至鏈輪中心距離H均為未知但有界的不確定變量，其區(qū)間分別為：RI1=[30,40] mm，RI2=[8,10] mm，RI3=[24,36] mm，RI4=[7,11] mm和HI=[125,151] mm， 齒頂位移最大許用值為dmax.

圖4 鏈輪Fig. 4 A chain wheel

鏈輪的功能函數(shù)為

本算例中樣本M=30， 參數(shù)c1，c2均 設(shè)置為0.2，慣性權(quán)重 ω均設(shè)為0.75，最大迭代次數(shù)為10 次，速度限制為V∈[?1,1]， 迭代收斂標(biāo)準(zhǔn)為| ηi?ηi?1|≤0.001.

表3 給出了相關(guān)系數(shù)為0.9 時本文方法的迭代過程；表4 列出了不同相關(guān)系數(shù)下，本文方法、Kriging 方法和參考值(MCS 法)的對比結(jié)果；表5 給出了MCS 求解的95%置信區(qū)間，表明參考值是有效的；圖5 給出了不同相關(guān)系數(shù)下兩種方法的相對誤差.

圖5 相對誤差對比圖（算例2）Fig. 5 Comparison of relative errors (example 2)

表3 相關(guān)系數(shù)為0.9 的迭代過程Table 3 The iterative process for a correlation coefficient of 0.9

表4 不同相關(guān)系數(shù)的可靠性指標(biāo)（算例2）Table 4 Reliability indexes for different correlation coefficients (example 2)

表5 MCS 求解的95%置信區(qū)間（算例2）Table 5 The 95% confidence intervals of MCS solutions (example 2)

算例3 如圖6 所示的25 桿桁架結(jié)構(gòu)[10]，桁架彈性模量E=199 949.2 MPa ， Poisson 比 μ=0.3，橫向、縱向桿件的長度為L，桿件① ～ ④、? ～ 25 、? ～ ?、⑤ ～ ⑩的截面積分別表示為A1，A2，A3和A4.作用在節(jié)點(diǎn)7、9 和11 的豎直載荷分別為F3=1 779.2 kN，F(xiàn)2=2 224 kN和F1=1 779.2 kN， 節(jié)點(diǎn)1 的水平載荷F4=1 334.4 kN，節(jié)點(diǎn)6 水平位移的最大許用值為dmax. 由于制造與測量過程中出現(xiàn)的誤差，桿件的截面積Ai(i=1,2,3,4)以及桿件的長度L為不確定變量，其不確定區(qū)間如表6 所示.

表6 不確定變量的分布參數(shù)Table 6 The distribution parameters of uncertain variables

圖6 25 桿桁架Fig. 6 A 25-bar truss structure

該桁架極限狀態(tài)函數(shù)為

本算例中樣本M=30， 參數(shù)c1，c2均 設(shè)置為0.2，慣性權(quán)重 ω均設(shè)為0.75，最大迭代次數(shù)為10 次，速度限制為V∈[?1,1]， 迭代收斂標(biāo)準(zhǔn)為| ηi?ηi?1|≤0.01.

表7 給出了相關(guān)系數(shù)為0 時本文方法的迭代過程；表8 列出了不同相關(guān)系數(shù)下，本文方法、Kriging 方法和參考值(MCS 法)的對比結(jié)果；表9 給出了MCS 求解的95%置信區(qū)間，表明參考值是有效的；圖7 給出了不同相關(guān)系數(shù)下兩種方法的相對誤差.

圖7 相對誤差對比圖（算例3）Fig. 7 Comparison of relative errors (example 3)

表7 相關(guān)系數(shù)為0 的迭代過程Table 7 The iterative process for a correlation coefficient of 0

表8 不同相關(guān)系數(shù)的可靠性指標(biāo)（算例3）Table 8 Reliability indexes for different correlation coefficients (example 3)

表9 MCS 求解的95%置信區(qū)間（算例3）Table 9 The 95% confidence intervals of MCS solutions (example 3)

由表1、3、7 可以看出，本文所提的PSO-Kriging 方法均能逐漸收斂于參考值附近，且與參考值誤差較小.

由表2、4、8 可知，不同相關(guān)系數(shù)下，相比Kriging 方法，PSO-Kriging 方法的計(jì)算時間較短，且給出的可靠性指標(biāo)更接近參考解，表明本文方法具備更高的計(jì)算效率與精度.

由圖3、5、7 可以看出，PSO-Kriging 方法的相對誤差對相關(guān)系數(shù)的變化不敏感，具有較強(qiáng)的魯棒性，而Kriging 方法則對相對系數(shù)的變化較敏感.

4 結(jié) 論

本文提出了一種PSO 與Kriging 模型相結(jié)合的非概率可靠性分析方法，該方法采用PSO 方法得到最優(yōu)相關(guān)參數(shù)，避免了Kriging 模型的單向搜索及對初始值敏感的缺點(diǎn).三個算例分析結(jié)果表明，本文方法有效可行，且精度和效率均高于Kriging 模型非概率可靠性分析方法.

參考文獻(xiàn)（ References ） ：

［1］BEN-HAIM Y. A non-probabilistic concept of reliability[J].Structural Safety, 1994, 14（4）: 227-245.

［2］ELISHAKOFF I, BEN-HAIM Y. Discussion on: a non-probabilistic concept of reliability[J].Structural Safety,1995, 17（3）: 195-199.

［3］ 郭書祥, 呂震宙, 馮元生. 基于區(qū)間分析的結(jié)構(gòu)非概率可靠性模型[J]. 計(jì)算力學(xué)學(xué)報(bào), 2001, 18（1）: 56-60. （GUO Shuxiang, Lü Zhengzhou, FENG Yuansheng. A non-probabilistic model of structural reliability based on interval analysis[J].Chinese Journal of Computational Mechanics, 2001, 18（1）: 56-60.（in Chinese））

［4］JIANG T, CHEN J J, XU Y L. A semi-analytic method for calculating non-probabilistic reliability index based on interval models[J].Applied Mathematical Modelling, 2007, 31（7）: 1362-1370.

［5］ 曹鴻鈞, 段寶巖. 基于凸集合模型的非概率可靠性研究[J]. 計(jì)算力學(xué)學(xué)報(bào), 2005, 22（5）: 546-549, 578. （CAO Hongjun, DUAN Baoyan. An approach on the non-probabilistic reliability of structures based on uncertainty convex models[J].Chinese Journal of Computational Mechanics, 2005, 22（5）: 546-549, 578.（in Chinese））

［6］JIANG C, LU G Y, HAN X, et al. Some important issues on first-order reliability analysis with non-probabilistic convex models[J].ASME Journal of Mechanical Design, 2014, 136（3）: 034501.

［7］MENG Z, HU H, ZHOU H L. Super parametric convex model and its application for non-probabilistic reliabilitybased design optimization[J].Applied Mathematical Modelling, 2018, 55: 354-370.

［8］MENG Z, ZHANG Z, ZHOU H L. A novel experimental data-driven exponential convex model for reliability assessment with uncertain-but-bounded parameters[J].Applied Mathematical Modelling, 2020, 77（1）: 773-787.

［9］WANG X, QIU Z P, ELISHAKOFF I. Non-probabilistic set-theoretic model for structural safety measure[J].Acta Mechanica, 2008, 198: 51-64.

［10］JIANG C, BI R G, LU G Y, et al. Structural reliability analysis using non-probabilistic convex model[J].Computer Methods in Applied Mechanics and Engineering, 2013, 254: 83-98.

［11］ 陳江義, 文尉超, 王迎佳. 基于響應(yīng)面方法的結(jié)構(gòu)非概率可靠性分析[J]. 鄭州大學(xué)學(xué)報(bào)(工學(xué)版), 2015, 36（5）: 121-124. （CHEN Jiangyi, WEN Weichao, WANG Yingjia. Analysis of structural non-probabilistic reliability based on response surface method[J].Journal of Zhengzhou University(Engineering Science), 2015, 36（5）: 121-124.（in Chinese））

［12］BAI Y C, HAN X, JIANG C, et al. A response-surface-based structural reliability analysis method by using nonprobability convex model[J].Applied Mathematical Modelling, 2014, 38（15/16）: 3834-3847.

［13］ 馬超, 呂震宙. 隱式極限狀態(tài)方程的非概率可靠性分析[J]. 機(jī)械強(qiáng)度, 2009, 31（1）: 45-50. （MA Chao, Lü Zhengzhou. Non-probablistic reliability analysis method for implicit limit state function[J].Journal of Mechanical Strength, 2009, 31（1）: 45-50.（in Chinese））

［14］ 潘林鋒, 周昌玉, 陳士誠. 基于Kriging模型的非概率可靠度計(jì)算[J]. 應(yīng)用力學(xué)學(xué)報(bào), 2010, 27（4）: 791-794. （PAN Linfeng, ZHOU Changyu, CHEN Shicheng. Non-probability reliability calculation based on Kriging model[J].Chinese Journal of Applied Mechanics, 2010, 27（4）: 791-794.（in Chinese））

［15］ 鄭嚴(yán). 基于智能算法的結(jié)構(gòu)可靠性分析及優(yōu)化設(shè)計(jì)研究[D]. 博士學(xué)位論文. 成都: 西南交通大學(xué), 2012. （ZHEN Yan.Research on structural reliability analysis and optimization based on intelligence algorithm[D]. PhD Thesis.Chengdu: Southwest Jiaotong University, 2012.（in Chinese））

［16］JIANG C, HAN X, LU G Y, et al. Correlation analysis of non-probabilistic convex model and corresponding structural reliability technique[J].Computer Methods in Applied Mechanics & Engineering, 2011, 200（33/36）:2528-2546.

［17］ 屈力剛, 劉洪俠, 李銘, 等. 基于Sobol序列采樣點(diǎn)分布策略的研究與應(yīng)用[J]. 鍛壓裝備與制造技術(shù), 2019, 54（6）: 101-105. （QU Ligang, LIU Hongxia, LI Ming, et al. Research and application of sampling point distribution strategy based on Sobol sequence[J].China Metalforming Equipment & Manufacturing Technology, 2019, 54（6）: 101-105.（in Chinese））