0 引言
采用代理模型技术来近似逼近结构功能函数是目前机械结构可靠度计算领域的研究热点，其中经典响应面法是一种较为常用的可靠度计算方法，它通过二次多项式模型逼近结构功能函数，并结合一次二阶矩法来获得可靠度指标及验算点[1]，计算过程简单且易于实现。BUCHER等[2]的研究发现，经典响应面法适用于一般的线性结构功能函数求解可靠度问题。而一些学者认为经典响应面法只能在一定程度上反映结构功能函数的非线性，当结构功能函数的非线性程度较高时，该方法就很难真实地反映函数的非线性程度，这就会造成计算精度上的误差[3]。鉴于上述问题，一些其他形式的代理模型逐渐应用到响应面法中来代替多项式模型，如神经网络响应面模型[4]、支持向量机模型[5]和Kriging模型等。与其他模型相比，Kriging模型有如下两个优点：①Kriging模型可利用样本点处的有效信息进行构造，而非全部的信息；②Kriging模型可预测已知样本信息的不确定性。基于上述特性，Kriging模型已逐渐应用于国内外的结构可靠度计算领域。RANJAN等[6]、BICHON等[7]分别提出了两种不同的自主选点改善函数，来提高Kriging模型对结构功能函数的拟合精度。ECHARD等[8]提出了一种新型的选点策略，并将其应用在可靠度计算中。谢延敏等[9]通过Kriging模型预测结构响应来进行可靠度计算。张崎等[10]结合Kriging模型与重要抽样法，提高了可靠度计算精度。
本文在Kriging响应面法基础上，提出了如下三点改进措施：①提出双点加点策略；②采用遗传算法优化Kriging模型的重要参数；③结合重要抽样法提高计算精度。最后，用算例验证了所提出方法的可行性。
1 Kriging代理模型的双点加点策略
考虑到验算点附近是失效概率较大的区域，对可靠度计算结果的影响极大，本文只需要在验算点处对结构功能函数进行局部精确拟合，即在迭代过程中采用某种加点策略来提高样本质量，并充分使用较少的样本信息使Kriging模型精确逼近真实验算点附近的真实模型。
双点加点是指在每次迭代只增加两个样本点，其中每次迭代产生的验算点是必须增加的，此外还需增加一个对可靠度计算具有重要价值的样本点。当样本点越靠近极限状态曲面(或称为失效面) 并具有越大的概率密度时，其重要性就越大，即越靠近极限状态曲线处的区域，其失效概率越大，对可靠度计算结果的影响也越大。为此，本文提出了一种基于Kriging预测的评价函数来判定样本点的价值程度，具体方法如下：①采用均匀设计方法生成构造Kriging模型的初始样本点；②通过蒙特卡罗(Monte Carlo，MC)法在抽样中心处生成一定数量的候选样本集；③利用初始Kriging模型预测候选样本集中的响应值；④使用评价函数从候选样本集中选出价值程度最大的样本点作为最佳样本点。
当样本点具有较大的概率密度且越靠近极限状态曲面时，该样本点越重要，故定义评价函数的表达式如下：
(1)
其中，g(x)为结构功能函数，f(x)为概率密度函数。从候选样本集中找到C(x)最小的点作为最佳样本点。
对于候选样本集，抽样中心应尽量靠近极限状态曲线(即结构功能函数g(x)=0)，这样生成的候选样本集效果更好。在经典响应面法中，由线性插值得到新的展开点往往会比验算点更加靠近极限状态曲线，其表达式如下：
(2)
式中，
为插值得到的第k个展开点；(x*)(k)为利用响应面法得到的第k个验算点；
为样本点的均值。
以二维结构功能函数为例，分别将验算点和展开点作为抽样中心生成候选样本集，图1给出了候选样本集的分布情况，其中x1、x2表示随机变量。
(a)验算点为抽样中心
(b)展开点为抽样中心
图1 候选样本集分布
Fig.1 Distribution of candidate samples
从图1中可以看出，以验算点x*为抽样中心时，生成的候选样本集分布在距离极限状态曲线较远的区域内，只有少量的样本点落入极限状态曲线附近，因此获得最佳样本点的概率会偏低。而以插值得到的展开点
为抽样中心时，其候选样本集的分布效果明显更加理想，只需生成少量的候选样本点，就有足够多的点落入极限状态曲线附近，该方法不仅可以提高获得最佳样本点的概率，同时减少了Kriging预测的工作量。但当结构功能函数呈现较高非线性程度时，依据式(2)得到的展开点到极限状态曲线的距离不一定比验算点到极限状态曲线的距离更近[11]，图2给出了一维高度非线性结构功能函数，其展开点处的结构功能函数值
大于验算点处的结构功能函数值g((x*)(k))，即展开点距离极限状态曲线更远。针对这种情况，本文定义了新的内插公式来确定下一次展开点，以保证与验算点相比,展开点距极限状态曲线更近，其表达式如下：
(3)
式中，
为插值得到的第k-1个展开点。
图2 非线性结构功能函数中的验算点与展开点
Fig.2 The checking point and the expansion point in the nonlinear structural performance function
综上所述，在Kriging模型迭代重构过程中可根据新的内插公式(式(3))来确定抽样中心，并采用双点加点策略来更新模型。
2 Kriging代理模型核函数的参数优化
由Kriging模型的相关理论可知，相关函数的选择对Kriging模型拟合能力的影响较为显著。从图3中可以看出，采用高斯核函数拟合出的模型曲面平滑度要比线性核函数拟合出的模型曲面平滑度高得多。曲面光滑可保证在迭代求解过程中的每一点都是可导可微的，这对可靠度指标和验算点采用一次二阶矩法(first order reliability method，FORM)求解时涉及到泰勒级数展开和求导是非常有利的[12]，故本文使用高斯核函数作为Kriging模型的相关函数。
核函数中的参数θ对Kriging模型的建立有着极其重要的影响，可通过求解如下优化问题得到[13]：
(a)高斯核函数拟合图形
(b)线性核函数拟合图形
图3 Kriging模型拟合的函数曲面
Fig.3 Function surface fitted by Kriging model
(4)
式中，m为设计变量个数；σ2为方差；R(θ)为参数θ的相关函数。
针对高斯核函数，式(4)可转化为一个最小化问题，即
(5)
姚拴宝等[14]利用DACE工具箱中编制的模式搜索方法对参数θ进行寻优，但该搜索方式对起始点的依赖性很强，若采用该方法确定参数θ，则会影响后续可靠性求解的精度。为解决上述问题，本文采用一种全局性搜索算法，利用Gatbx工具箱中的遗传算法(GA)对式(5)中的参数θ进行优化，进而得到优化Kriging的模型。
3 重要抽样法修正
为提高MC法的抽样效率，通过改变随机抽样的“重心”，增加结构功能函数g(x)<0的机会，使得抽取的样本点有较多的机会落入失效域内，这就是重要抽样法的基本思想。
结构失效概率的表达式如下：
(6)
式中，px(x)为重要抽样概率密度函数;I(g(x))为示性函数。
为了避免FORM法处理高度非线性结构功能函数所产生的误差，本文根据重要抽样思想对FORM法的可靠度计算结果进行了修正。
重要抽样法的关键在于如何获得真实验算点的位置[15]，而本文所提出的改进Kriging响应面法能够找到较为精确的验算点，故可以采用重要抽样法对FORM法得到的可靠度指标或失效概率进行修正。迭代过程中得到的可靠度指标和验算点是不同的，因此没有必要对每次迭代的计算结果均进行修正，只需修正最后一次迭代的计算结果。
4 算法流程
结合上述提出的改进措施，基于双点加点策略的改进Kriging响应面法的可靠度计算流程见图4，具体步骤如下。
图4 可靠度算法流程图
Fig.4 Reliability algorithm flow chart
(1)初始设计点可以选择在均值点处，即
(2)采用均匀设计法在
范围内生成ns个样本点，得到初始样本点集合X0={x1,x2，…,xns}，并通过有限元分析得到相应样本点的真实结构功能函数值集合Y0={y1,y2，…,yns}。其中，f为大于零的任意数，σx为标准差。
(3)优化参数θ，并由初始样本点集合X0和相应样本点的真实结构功能函数值集合Y0来构造结构功能函数g(x)的优化Kriging模型。
(4)采用FORM法求解当前验算点(x*)(k)及其可靠度指标β(k)，k(k=0,1,2,…)为序列号。
(5)计算当前验算点(x*)(k)的真实结构功能函数值(y*)(k)，并将该点加入样本库Sk中。
(6)依据式(2)生成候选样本集抽样中心
若
则使用新的内插公式(式(3))来确定抽样中心
(7)采用MC法以为抽样中心，生成10 000点候选样本集，根据Kriging模型预测候选样本集中的响应值，然后将其代入评价函数C(x)进行计算，并以此选取评价函数值C(x)最小的点作为最佳样本点(x′)(k)。
(8)计算当前最佳样本点(x′)(k)的真实结构功能函数值(y′)(k)，并将该点也加入到样本库Sk中。
(9)返回步骤(3)，直到前后两次失效概率满足如下收敛条件：
(7)
其中，ε一般取0.001。
(10)以验算点(x*)(k)为重要抽样中心，采用重要抽样法对计算结果进行修正。
5 算例
5.1 算例1
简支梁结构如图5所示，集中力偶Me作用在截面C处，已知b=1.2 m，L=2 m，该结构功能函数可表示为
(8)
式中，Me为施加的力偶，kN·m；D为梁截面的直径，m。
图5 简支梁示意图
Fig.5 Schematic diagram of simple supported beam
本算例将式(8)中的Me、D作为随机变量，其分布参数见表1。
表1 随机变量的分布参数(算例1)
Tab.1 Distribution parameters of random variables(case one)
随机变量均值分布类型相应变异系数力偶Me(kN·m)1930对数正态0.1直径D(mm)80正态0.1
图6 遗传算法对参数θ的优化过程
Fig.6 Optimization process of θ by genetic algorithm
本文采用均匀设计生成40个初始样本，图6所示为参数θ的迭代优化过程，遗传代数为50，种群数量为100，可以看出，迭代到第30代后处于收敛状态，目标函数最优值为0.06×10-3。
为了与优化Kriging模型相比较，笔者构造了基于本文所提方法的标准Kriging模型进行可靠度计算，图7给出了两种模型的可靠度指标迭代过程，可以看出，标准Kriging模型共迭代了7次，优化Kriging模型迭代了6次。两种模型迭代过程中的新增样本点分别见图8和图9。对比图8和图9中的验算点和最佳样本点可知，与标准Kriging模型相比，优化Kriging模型的验算点和最佳样本点更加贴近极限状态曲线，这表明优化Kriging模型产生的验算点和最佳样本点的重要性较高，能够在有限的样本数量下更好地反映可靠度的计算结果。由图10可以看出，在验算点处两种模型拟合的曲线均与真实极限状态曲线几乎完全重合，其中优化Kriging模型的全局拟合效果更好。
图7 可靠度指标迭代过程
Fig.7 Iterative process of reliability index
图8 优化Kriging模型的新增样本点
Fig.8 New sample points in the optimized Kriging model
图9 标准Kriging模型的新增样本点
Fig.9 New sample points in the standard Kriging model
图10 Kriging模型拟合的极限状态曲线与真实极限状态曲线
Fig.10 Limit state curve fitted by Kriging model and the real limit state curve
表2列出了本文所提算法、经典响应面法(classical RSM)、Kriging响应面法(Kriging RSM)、重要抽样+Kriging模型(IS+Kriging)、重要抽样+BP神经网络模型(IS+BPANN)和MC法的可靠度计算结果，其中由MC法得到的结果为精确解并用于精度检验。通过比较各方法的结果可以发现，本文提出的基于双点加点策略的改进Kriging响应面法的计算结果精度和计算效率均有良好表现，且优化Kriging模型可以获得更高的计算精度。与经典响应面法和Kriging响应面法相比，改进Kriging响应面法的计算量大幅减少，这是因为该方法在迭代计算过程中，每次迭代只增加两个高质量的样本点，同时在下一次迭代时，并不舍弃之前的样本点，而是重复利用这些点来构造Kriging模型，这与经典响应面法在展开点处重构模型(舍弃前次迭代的样本点)有着本质的区别。由于IS+Kriging法和IS+BPANN法均为一步可靠度求解，不需要迭代计算，这样对初始样本的要求较高，其中IS+BPANN法的计算结果精度虽然很高，但由于BP神经网络对初始样本的依赖性较强，因此需要大量的样本点来构造该模型；IS+Kriging法的计算量虽不大，但精度达不到要求。
表2 不同方法的计算结果比较
Tab.2 Comparison of calculation results with different methods
方法本文方法优化Kriging本文方法标准KrigingclassicalRSMKrigingRSMIS+KrigingIS+BPANNMC法初始样本404054050100106迭代次数67157结构功能函数计算次数57608928650100106随机变量Me(kN·m)2059.152059.192046.832057.742056.262059.17随机变量D(mm)55.352.150.155.851.952.3β3.18673.18033.23053.17963.17433.18553.1879pf(10-3)0.71950.73560.61790.73740.75100.72250.7165相对误差(%)0.4182.66513.7612.9164.8150.837
最后分析不同的加点方式对Kriging响应面法的影响，除本文所提出的双点加点策略，笔者规定：单点加点策略Ⅰ，每次迭代只增加一个验算点；单点加点策略Ⅱ，每次迭代只增加一个最佳样本点。
每种加点方式的验算点迭代过程见图11，单点加点策略Ⅰ迭代11次收敛到验算点坐标(2 058.147 5 kN·m,55.1 mm)，计算了57次结构功能函数，可靠度指标为3.181 7；单点加点策略Ⅱ迭代8次收敛到验算点坐标(2 060.765 6 kN·m，55.2 mm)，计算了55次结构功能函数，可靠度指标为3.183 2；双点加点策略迭代6次收敛到验算点坐标(2 059.093 1 kN·m，55.1 mm)，计算了55次结构功能函数，可靠度指标为3.186 7。
1.单点加点策略Ⅰ坐标(2 058.147 5 kN·m,55.1 mm) 2.单点加点策略Ⅱ坐标(2 060.765 6 kN·m，55.2 mm) 3.双点加点策略坐标(2 059.093 1 kN·m，55.1 mm)
图11 验算点迭代过程
Fig.11 Iterative process of checking point
综合上述分析得出，在考虑精度和效率时，与单点加点策略相比，同时增加验算点和最佳样本点能够加速验算点迭代收敛，其计算结果也更加精确。
5.2 算例2
某集装箱驮背车上的活动鞍座结构见图12，材料为Q450EW高强度耐候钢，有限元模型及加载情况见图13。结构功能函数要求结构最大应力小于材料屈服极限，即
g=450-σs(Fv,E,ρ,t1,t2,t3)
(9)
式中，σs为应力函数；Fv为垂直向下的集中力；E为材料弹性模量；ρ为材料密度；t1、t2、t3分别为顶板、支撑臂及支撑垫板的厚度。
图12 活动鞍座结构
Fig.12 Structure of active saddle
图13 有限元模型及边界条件
Fig.13 Finite element model and boundary conditions
本算例将上述参数作为随机变量，其具体参数分布见表3。
表3 随机变量的分布参数(算例2)
Tab.3 Distribution parameters of random variables (case two)
随机变量均值分布类型相应变异系数Fv(kN)79.168对数正态0.10E(MPa)206000正态0.05ρ(kg/m3)7.8×103正态0.05t1(mm)12正态0.10t2(mm)45正态0.10t3(mm)30正态0.10
采用本文方法计算该结构可靠度，其计算结果和500次MC法校验结果见表4。从表4中可以看出，本文所提方法具有良好的计算精度，仅计算了34次结构功能函数，与500次抽样相比，显著减少了计算量。
表4 计算结果
Tab.4 Calculation results
方法本文方法优化Kriging本文方法标准KrigingMC法初始样本2020500迭代次数55结构功能函数计算次数3434500β3.05723.04353.0537pf(10-2)0.11170.11690.1130相对误差(%)1.15042.5663
6 结论
(1)在Kriging响应面法的基础上，通过定义评价函数和新的内插公式，提出了双点加点策略来更新Kriging模型，并利用遗传算法优化参数θ，以充分发挥有效样本点的信息，使Kriging模型在真实功能函数的验算点处精确拟合。并结合重要抽样法进行修正，提高了可靠度计算结果精度。
(2)算例分析结果表明，优化Kriging模型相比标准Kriging模型能够更加准确地模拟结构功能函数。与其他方法相比，本文方法减少了功能函数的计算次数，且得到的可靠度计算结果与MC法产生的精确值更加接近。
[1] 智鹏鹏, 李永华, 陈秉智. 考虑参数不确定性的转向架构架结构强度分析[J]. 中国机械工程, 2019, 30(1)：22-29.
ZHI Pengpeng, LI Yonghua, CHEN Bingzhi. Structural Strength Analysis of Bogie Frames Considering Parameter Uncertainty[J]. China Mechanical Engineering, 2019, 30(1)：22-29.
[2] BUCHER C G, BOURGAND U. A Fast and Efficient Response Surface Approach for Structural Reliability Problems[J]. Structural Safety, 1990, 7(1)：57-66.
[3] ELHEW A H, MESBAHI E, PU Y. Reliability Analysis of Structures Using Neural Network Method[J]. Probabilistic Engineering Mechanics, 2006, 21(1)：44-53.
[4] 郭姣姣, 刘伟, 翟玮昊,等. 随机故障注入结合神经网络法的机电系统可靠性计算方法[J]. 机械工程学报, 2017, 53(6)：195-202.
GUO Jiaojiao, LIU Wei, ZHAI Weihao, et al. Reliability Calculation Method of Electromechanical System Based on Random Fault Injection Combined with Artificial Neural Network[J]. Journal of Mechanical Engineering, 2017, 53(6)：195-202.
[5] 赵安龙, 李洪双. 基于概率支持向量机的可靠性分析与设计方法[J]. 应用力学学报, 2017, 34(1) :50-56.
ZHAO Anlong, LI Hongshuang. The Analysis of Stochastic Response with Equivalent Linearization Method in Weak Nonlinear Systems[J]. Chinese Journal of Applied Mechanics, 2017, 34(1)：50-56.
[6] RANJAN P, BINGHAM D, MICHAILIDIS G. Sequential Experiment Design for Contour Estimation from Complex Computer Codes[J]. Technometrics, 2008, 50(4):527-541.
[7] BICHON B J, ELDRED M S, SWILER L P, et al. Efficient Global Reliability Analysis for Nonlinear Implicit Performance Functions[J]. AIAA Journal, 2008, 46(10):2459-2468.
[8] ECHARD B, GAYTON N, LEMAIRE M. AK-MCS：an Active Learning Reliability Method Combining Kriging and Monte Carlo Simulation[J]. Structural Safety, 2011, 33(2)：145-154.
[9] 谢延敏, 于沪平, 陈军,等. 基于Kriging模型的可靠度计算[J]. 上海交通大学学报, 2007,41(2)：177-180.
XIE Yanmin, YU Huping, CHEN Jun, et al. The Reliability Estimation Based on Kriging Model[J]. Journal of Shanghai Jiao Tong University, 2007, 41(2)：177-180.
[10] 张崎, 李兴斯. 结构可靠性分析的模拟重要抽样方法[J]. 工程力学, 2007, 24(1)：33-36.
ZHANG Qi, LI Xingsi. Analysis of Structural Reliability Based on Kriging Model[J]. Engineering Mechanics, 2007, 24(1)：33-36.
[11] 张哲, 李生勇, 石磊, 等. 结构可靠度分析中的改进响应面法及其应用[J]. 工程力学, 2007, 24(8)：111-115.
ZHANG Zhe, LI Shengyong, SHI Lei, et al. An Improved Response Surface Method for Structural Reliability Analysis and Application[J]. Engineering Mechanics, 2007, 24(8)：111-115.
[12] 杨文军, 张锴锋, 王磊, 等. 基于Kriging插值的失谐叶盘气动结构耦合振动特性分析[J]. 中国机械工程, 2018, 29(12):1405-1413.
YANG Wenjun, ZHANG Kaifeng, WANG Lei, et al. Aerodynamics-structure Coupling Virbrations of Mistuning Blade-disk Based on Kriging Interpolation[J]. China Mechanical Engineering, 2018, 29(12):1405-1413.
[13] 陈志英, 任远, 白广忱,等. 粒子群优化的Kriging近似模型及其在可靠性分析中的应用[J]. 航空动力学报, 2011, 26(7)：1522-1530.
CHEN Zhiying, REN Yuan, BAI Guangchen, et al. Particle Swarm Optimized Kriging Approximate Model and Its Application to Reliability Analysis[J]. Journal of Aerospace Power, 2011, 26(7)：1522-1530.
[14] 姚拴宝, 郭迪龙, 孙振旭,等. 基于Kriging代理模型的高速列车头型多目标优化设计[J]. 中国科学:技术科学, 2013, 43(2)：186-200.
YAO Shuanbao, GUO Dilong, SUN Zhenxu, et al. Multi-objective Optimization of the Streamlined Head of High-speed Trains Based on the Kriging Model[J]. Sci. China:Tech. Sci., 2013, 43(2)：186-200.
[15] 董文俊, 夏艳君, 唐樟春,等. 基于Kriging模型的迭代重要抽样方法研究及应用[J]. 机械设计与研究, 2016(2):6-10.
DONG Wenjun, XIA Yanjun, TANG Zhangchun, et al. Kriging Model Based Iterative Importance Sampling Method and Its Application to an Inside Flap of an Aircraft[J]. Machine Design and Research, 2016(2)：6-10.
