机械可靠性分析的高精度响应面法

文章编号 :1000-0887 (2007) 01-0017-08 ν 应用数学和力学编委会 ,ISSN 1000-0887 机械可靠性分析的高精度响应面法 Ξ 吕震宙 , 　赵　洁 , 　岳珠峰 (西北工业大学 航空学院 ,西安 710072) (我刊编委岳珠峰来稿) 摘要 : 　通过对已有可靠性分析中的响应面法的研究 ,提出了一种高精度的响应面法 ,该方法通过 迭代线性插值的策略 ,来保证确定响应面的抽样点比经典的响应面法更接近真实的极限状态方 程 ,并且该方法通过序列线性插值的方法来控制抽样点与插值中心点的距离 ,保证随着插值中心 点收敛于真实 设计 领导形象设计圆作业设计ao工艺污水处理厂设计附属工程施工组织设计清扫机器人结构设计 点 ,抽样点提供更多的关于设计点附近真实极限状态方程的信息 ,进而保证了 收敛的响应面能够在设计点附近更好地拟合真实的极限状态方程 ,并得到高精度的失效概率计算 结果·　算例充分说明了所提方法的合理性与适用性·　 关 　键 　词 : 　响应面法 ; 　隐式极限状态 ; 　失效概率 中图分类号 : 　TB114. 3 　　　文献标识码 : 　A 引　　言 对于具有隐式极限状态方程的结构系统 ,传统的一次可靠性方法和二次可靠性方法[1-8 ]难 以实施 ,在此情况下响应面法被广泛地推荐使用[9-24 ] ,响应面法的基本思想是 :采用响应面函 数来构造隐式极限状态方程的近似显式 表 关于同志近三年现实表现材料材料类招标技术评分表图表与交易pdf视力表打印pdf用图表说话 pdf 达式[9-10 ]·　已有的工作表明 ,如果响应面函数能够 很好的近似实际的隐式极限状态方程 ,它将可以得到精度相当高的失效概率估计值·　为了提 高响应面法的计算精度 ,必须考虑两个方面的因素 :响应面函数的确定和试验点的确定·　最常 用的响应面函数是含待定常数的多项式 ,通过设计试验点采用回归分析或拟合的方法来确定 多项式中的待定常数[9 ,12 ,14 ]·　通常多项式次数的提高可以得到更高精度的计算结果 ,但这是 以付出更多的计算工作量为代价的[12 ,16 ] ,考虑计算工作量以及数学概念方面的因素 ,一般采 用二次多项式 (更多的是采用不含交叉项的二次多项式)作为响应面函数的形式·　 本文主要研究试验点的选取方法 ,响应面法分析隐式极限状态方程可靠性的精度依赖于 试验点[9-16 ]·　由于设计点附近的区域 ,或者说失效域中基本变量联合概率密度较大的区域对 失效概率的贡献大 ,因此响应面函数应该在此区域对真实极限状态函数有较高的近似精度 ,由 此得到了试验点选择的一般原则 ,即试验点应落在设计点附近·　大部分文献中采用序列响应 面的方法来实现这个原则以提高计算精度·　文献[18 ]采用加权响应面法来近似隐式极限状态 方程 ,以提高可靠性分析的精度 ,并指出极限状态方程才是响应面法应着重近似的·　文献 [ 14 ] 71 　应用数学和力学 ,第 28 卷 第 1 期 　2007 年 1 月 15 日出版 　　　　　　　　　　　　　　　　 Applied Mathematics and Mechanics 　　　Vol. 28 ,No. 1 ,Jan. 15 ,2007 　 Ξ 收稿日期 : 　2005-11-15 ; 修订日期 : 　2006-10-30 基金项目 : 　国家自然科学基金资助项目 (10572117) ;新世纪优秀人才支持计划资助项目 (NCET-05- 0868) 作者简介 : 　吕震宙 (1966 —) ,女 ,湖北黄石人 ,教授 ,博士生导师 ( E-mail :zhenzhoulu @nwpu. edu. cn)·　 则采用梯度投影的方法来选择响应面的试验点 ,该方法的核心是希望确定响应面的试验点落 在离真实的极限状态方程更近的区域上·　梯度投影方法的主要局限性在于 :其所选择的试验 点需要通过摄动来保证求解响应面的矩阵非病态 ,而试验点的摄动很难控制·　本文在已有方 法的基础上提出了一种直接的将试验点选在更接近实际的失效面上的方法 ,这种方法不需要 摄动 ,并且所提方法还通过序列线性插值来强调了设计点的重要性·　 1 　确定响应面试验点的改进方法 1. 1 　传统试验设计方法 为了讨论简单起见 ,假设基本随机变量 x i ( i = 1 ,2 , ⋯, n) 是不相关的正态变量 ,其均值和 标准差分别为μi 和σi ,也即 xi ～ N (μi ,σ2i)·　以式 (1) 所示的二次不含交叉项的响应面多项式…g ( €x) 替代隐式极限状态函数 g ( €x) [9 ]·　 　　…g ( €x) = a0 + 6n i = 1 bixi + 6n i =1 cix 2 i , (1) 其中 a0、bi 、ci ( i = 1 ,2 , ⋯, n) 是 2 n + 1 待定常数 ,它们可以通过试验设计来确定·　 为了确定 2 n + 1 待定常数 ,至少需要 2 n + 1 试验点 ,传统响应面法所选择的 2 n + 1 个试 验点如 (2)式至 (4)式所示·　 　　€x1 = (μ1 ,μ2 , ⋯,μn) , (2) 　　€x j = (μ1 ,μ2 , ⋯,μi + fσi , ⋯,μn) , 　　j = 2 ,3 , ⋯, n + 1 ; i = j - 1 , (3) 　　€x j = (μ1 ,μ2 , ⋯,μi - fσi , ⋯,μn) , 　　j = n + 2 , n + 3 , ⋯,2 n + 1 ; i = j - ( n + 1) , (4) 其中试验点 €x1 是基本随机向量 …x 的均值向量 …xμ = (μ1 ,μ2 , ⋯,μn) 的试验点 ,它是中心试验 点 ,其它的周围试验点 €x j ( j = 2 ,3 , ⋯,2 n + 1) 是在中心试验点附近沿每个坐标轴的正、负方向 偏离 fσi 而形成的 , f 被称为偏离系数 ,常取经验值 1～3 ,它确定了试验点的取值范围·　 通过在试验点处进行拟合 ,可以唯一地确定 (1) 式中的待定常数·　以 a (0)0 、b (0)i 、c (0)i 表示 确定的示于式 (5) 中的响应面函数 …g (0) ( €x) 中的常数 ,其中上标 (0) 表示第一次确定响应面函 数·　 　　…g (0) ( €x) = a (0)0 + 6n i = 1 b (0)i xi + 6n i = 1 c (0) i x 2 i = 0·　 (5) 对于具有显式表达的极限状态方程 …g (0) ( …x) = 0 ,可以采用一次二阶矩方法来确定其设计 点 €x (0)D = ( x (0)D1 , x (0)D2 , ⋯, x (0)D n )·　通过均值点 ( …xμ , g ( …xμ) ) 和设计点 ( €x (0)D , g ( €x (0)D ) ) 的线性插值 , 可以得到一个接近 g ( €x) = 0的新的试验中心点 €x (1)M = ( x (1)M1 , x (1)M2 , ⋯, x (1)M n ) ,其第 i 个坐标€x (1)M i 如下所示[9 ]·　 　　 x (1)M i = μi + ( x (0)D i - μi) g ( …xμ) g ( …xμ) - g ( €x (0)D ) ·　 (6) 以 €x (1)M 替代 (2)式中的中心试验点 ,可得到新的试验点如下 　　€x (1)1 = €x (1)M = ( x (1)M1 , x (1)M2 , ⋯, x (1)M n) , (7) 　　€x (1)j = ( x (1)M1 , x (1)M2 , ⋯, x (1)M i + fσi , ⋯, x (1)M n ) , 　　j = 2 ,3 , ⋯, n + 1 ; i = j - 1 , (8) 　　€x (1)j = ( x (1)M1 , x (1)M2 , ⋯, x (1)M i - fσi , ⋯, x (1)M n ) , 　　j = n + 2 , n + 3 , ⋯,2 n + 1 ; i = j - ( n + 1)·　 (9) 81 吕 　震 　宙 　　　赵 　　洁 　　　岳 　珠 　峰 在第一次修正的试验点处 ,又可以得到新的响应面函数 …g (1) ( €x) ,此修正过程可以一直进 行下去 ,直至收敛准则被满足为止·　一般的收敛准则可表示为不等式 | €x ( k)D - €x ( k - 1)D | ≤ε(ε 是预先给定的小数)成立·　 1. 2 　改进的试验设计方法 为了提高响应面法的精度 ,选择试验点应满足以下两个准则 :1) 试验点应该位于真实的 极限状态方程附近 ,以使得 g ( €x) = 0能够被 …g ( €x) = 0高精度近似 ;2) 试验点应落在设计点附 近 ,以便使得对失效概率贡献大的区域能够很好的被近似·　 为满足上述准则 ,本文提出了具有两个创新点的改进试验设计方法·　其一是通过均值点 与传统周围试验点的线性插值来得到改进的周围试验点 ,这可以保证改进的试验点落在真实 的失效面附近 ;其二是采用序列线性插值来控制周围试验点与中心试验点的距离 ,以保证中心 试验点收敛于真实设计点时 ,周围试验点收敛于对失效概率贡献大的区域·　通过这两点改进 , 可以提高响应面方程对真实极限状态方程重要区域的近似 ,进而提高可靠性分析的精度·　以 下给出了改进试验设计方法的步骤·　 1) 选择二次不含交叉项的多项式作为响应面函数的形式 ,第一次确定响应面函数的试验 点与传统方法相同 ,第一次修正响应面的传统试验点采用如下 (10) 式至 (12) 式所示的形式表 示 : 　　€x (1)1c = ( x (1)1c1 , x (1)1c2 , ⋯, x (1)1c n) = ( x (1)M1 , x (1)M2 , ⋯, x (1)M n ) , (10) 　　€x (1)j = ( x (1)jc1 , x (1)jc2 , ⋯, x (1)jc n) = ( x (1)M1 , x (1)M2 , ⋯, x (1)M i + fσi , ⋯, x (1)M n ) , 　　j = 2 ,3 , ⋯, n + 1 ; i = j - 1 , (11) 　　€x (1)j = ( x (1)jc1 , x (1)jc2 , ⋯, x (1)jc n) = ( x (1)M1 , x (1)M2 , ⋯, x (1)M i - fσi , ⋯, x (1)M n ) , 　　j = n + 2 , n + 3 , ⋯,2 n + 1 ; i = j - ( n + 1) , (12) 其中下标中的字母“c”表示传统方法 ,上标中的“(1)”表示第一次修正响应面函数 ; 2) 以 ( €x (1)jc , g ( €x (1)jc ) ) ( j = 2 ,3 , ⋯,2 n + 1) 和 ( …xμ , g ( …xμ) ) 进行线性插值 ,可得到 g ( €x T(1)jA ) ≈ 0 的改进的周围试验点 €x T(1)jA , €x T(1)jA 的第 i 个坐标 xT(1)jA i 如下所示 　　 xT(1)jA i = μi + ( x (1)jc i - μi) g ( …xμ) g ( …xμ) - g ( €x (1)jc ) , i = 1 ,2 , ⋯, n ; j = 2 ,3 , ⋯,2 n + 1 , (13) x T(1) jA i 用来暂时存放改进的周围试验点的第 i 个坐标 ; 3) 如果下列不等式成立 　　 k0σj - 1 < | €x T(1)jA - €x (1)1c | < k1σj - 1 , 　　j = 2 ,3 , ⋯,2 n + 1 , (14) 　　 k0σj - ( n +1) < | €x T(1)jA - €x (1)1c | < k1σj - ( n +1) , 　　j = n + 2 , n + 3 , ⋯,2 n + 1 , (15) 则最终的改进的周围设计点 €x (1)jA = ( x (1)jA1 , x (1)jA2 , ⋯, x (1)jA n) 就被确定了 ,它就等于暂时周围试验 点 ,即 　　€x (1)jA i = €x T(1)jA1 , 　　( i = 1 ,2 , ⋯, n ; j = 2 ,3 , ⋯,2 n + 1) , (16) 在不等式 (14)和式 (15)中 , k0 和 k1 是用来控制周围试验点与中心试验点距离的经验参数 ,太 远会降低设计点的重要性 ,而太近的距离又会使得求解响应面系数的矩阵病态 ,建议的经验取 值为 k0 取 1 ～ 1. 5 ,而 k1 取 2～3 ; 4) 如果不等式 (14) 和式 (15) 不成立 ,则需要采用序列线性插值来选取合适的周围试验 点 ,此时先计算 €x T(1)jA ( j = 2 ,3 , ⋯,2 n + 1) 和 €x (1)1c 的中点 , 然后以中点和均值点进行线性插值 91机械可靠性分析的高精度响应面法 来得到新的周围试验点 , 此线性插值的过程可以一直进行下去直至不等式 (14) 和式 (15) 被满 足 ; 5) 在确定的改进试验点处修正响应面函数 ,响应面函数的修正一直进行下去 ,直至响应 面的收敛准则被满足·　 通过本文提出的改进 ,周围试验点被选在了设计点附近 ,并且更直接的接近真实的失效 面 ,这种改进使得设计点的重要性被充分考虑 ,从而使得设计点附近的区域能够被高精度的近 似 ,进而使得可靠性分析的精度也得到了提高·　 2 　算　　例 2. 1 　算例 1 非线性极限状态函数 g ( €x) = exp (0. 2 x1 + 6. 2) - exp (0. 47 x2 + 5. 0) 取自文献[14 ] ,其中 变量 x1 和 x2 为相互独立的标准正态变量·　图 1 给出了利用本文方法与经典响应面法对算例 1 进行可靠性分析的结果对照 ,图1中可靠度指标β和失效概率 Pf 的相对误差均表示响应面法 相对于精确解的相对误差·　 　　　　　(a) 可靠度指标的对照 　　　　　　　　　　(b) 可靠度指标相对误差的对照 　　　　　　　(c) 失效概率的对照 　　　　　　　　　　(d) 失效概率相对误差的对照 图 1 　算例 1 中本文方法与经典响应面法结果随偏离系数变化曲线的对照 2. 2 　算例 2 非线性极限状态方程 g ( €x) = exp [0. 2 x1 + 1. 4 ] - x2 也取自文献[14 ] ,其中的变量 x1 和 x2 是相互独立的正态变量·　图 2 给出了利用本文方法与经典响应面法对算例 2 进行可靠性分析 的结果对照·　 从图 1 和图 2 的结果可以看出 ,当 f 取较大值时传统响应面法计算的失效概率结果与本 文所提方法的结果有较大区别 ,这种区别随 f 值的减小而减小·　本文所提方法计算的可靠度 02 吕 　震 　宙 　　　赵 　　洁 　　　岳 　珠 　峰 　　　　　(a) 可靠度指标的对照 　　　　　　　　　　(b) 可靠度指标相对误差的对照 　　　　　　(c) 失效概率的对照 　　　　　　　　　　　　(d) 失效概率相对误差的对照 图 2 　算例 2 中本文方法与经典响应面法结果随偏离系数变化曲线的对照 指标和失效概率的误差明显小于传统响应面法 ,而且这种优越性随 f 值增大而增大·　所提方 法的一个重要优点是其可靠性分析结果对 f 的取值不敏感 ,在较大的 f 取值范围内 ,本文所提 方法都可以得到稳健的解·　 2. 3 　算例 3 考虑一悬臂梁的自振频率 ,假设弹性模量 x1 ,材料密度 x2 ,梁的厚度 x3 和长度 x4 均为独 立的正态分布 ,其均值和变异系数如表 1 所示 ,对于此悬臂梁进行可靠性分析所建立的非线性 极限状态方程为 g ( x1 , x2 , x3 , x4) = 1. 016 x1 x23/ ( x2 x44) - 360·　图 3 给出了利用本文方法与 经典响应面法对算例 3 进行可靠性分析的结果对照·　 表 1 算例 3 中基本变量的分布形式和分布参数 基本变量 分布 均值 变异系数 x1 正态 107 N/ cm2 0. 1 x2 正态 2. 5 ×10 - 4 kg/ cm3 0. 1 x3 正态 0. 98 cm 0. 1 x4 正态 20. 0 cm 0. 4 　　图 3 的计算结果表明 f 在 1. 0 ～ 9. 0 之间取值时 ,所提方法均可以算出精确的结果 ,而传 统方法只在 f 取 1. 0 ～ 4. 0 之间的值时才可以算出有一定精度的近似解 ,这充分说明本文所提 方法对 f 取值不敏感 ,从而使得所提方法具有广泛适用性·　 本文亦对其它大量的算例进行了所提算法的验证 ,为节省篇幅 ,未将这些算例的结果列在 本文中·　这些未列入文中的大量算例结果也充分表明 ,所提算法的结果对偏离系数 f 的变化 12机械可靠性分析的高精度响应面法 是不敏感的·　对于经典响应面法得不到收敛解的情况 ,本文方法均可以得到高精度的结果 ,并 且在这种情况下 ,本文算法的计算工作量将大大小于传统的经典算法·　 　　　　　(a) 可靠度指标的对照 　　　　　　　　　　(b) 可靠度指标相对误差的对 照 　　　　　(c) 失效概率的对照 　　　　　　　　　　　　(d) 失效概率相对误差的对照 图 3 　算例 3 中本文方法与经典响应面法结果随偏离系数变化曲线的对照 3 　结　　论 1) 本文所提出的响应面的改进试验设计方法可以显著提高隐式极限状态方程可靠性分 析的精度·　在传统响应面法中 ,能够保证试验点落在实际失效面附近的线性插值只被用来选 择试验中心点·　而在所提出的方法中 ,确定周围试验点也采用了相同的线性插值的策略 ,以保 证周围试验点也落在真实失效面附近 ,进而使得真实失效面能够很好地被近似·　为了更好的 拟合对失效概率贡献大的区域中的极限状态方程 ,本文还提出了序列线性插值的方法 ,通过控 制周围试验点与中心试验点的距离 ,来保证周围试验点在中心试验点收敛于设计点时落在感 兴趣的区域·　 2) 数值算例充分说明了本文所提方法的优点 ,与传统的响应面法相比 ,所提方法不仅在 可靠性分析的精度方面有显著提高 ,而且在计算的稳定性方面也有较大提高·　本文算法的计 算工作量在每步迭代中较传统方法有所提高 ,但由于所提算法的收敛性较传统响应面方法好 , 因此总的计算工作量并不一定总是大于传统响应面方法·　 3) 本文所提方法的概念与梯度投影方法类似 ,但本文方法的实现更为直接和简便 ,因此 本文方法可以看作是梯度投影方法的进一步发展 ,另外 ,本文方法可以保证得到的确定响应面 函数的线性系统是非病态的·　 22 吕 　震 　宙 　　　赵 　　洁 　　　岳 　珠 　峰 [参 　考 　文 　献 ] [1 ] 　 Nowak A S , Collins K R. Reli ability of St r uct ures [M ] . Boston :McGraw- Hill ,2000. [2 ] 　 Ibrahim Y. Observations on applications of importance sampling in structural reliability analysis [J ] . St r uct ural Safety ,1991 ,9 (4) :269-281. [3 ] 　 Olsson A ,Sandberg G ,Dahlblom O. On latin hypercube sampling for structural reliability analysis [J ] . St r uct ural Safety ,2003 ,25 (1) :47-68. [4 ] 　 Melchers R E. Radial importance sampling for structural reliability [J ] . J Engng Mech , ASCE ,1990 , 116 (1) :189-203. [5 ] 　 Moarefzadeh M R , Melchers R E. Directional importance sampling for ill-proportioned spaces [J ] . St r uct ural Safety ,1999 ,21 (1) :1-22. [6 ] 　 Carlier J , Li Y, Lutton J L . Reliability evaluation of large telecommunication networks [J ] . Discrete Applied Mathematics ,1997 ,76 (1/ 3) :61-80. [7 ] 　 Melchers R E , Ahammed M , Middleton C. FORM for discontinuous and truncated probability density functions [J ] . St r uct ural Safety ,2003 ,25 (3) :305-313. [8 ] 　 Fiessler B , Neumann H J , Rackwitz R. Quadratic limit states in structural reliability [ J ] . J Engng Mech ASCE ,1979 ,105 (4) :661-676. [9 ] 　 Bucher C G. A fast and efficient response surface approach for structural reliability problems [J ] . St r uct ural Safety ,1990 ,7 (1) :57-66. [10 ] 　Faravelli . Response surface approach for reliability analysis [J ] . J Engng Mech ASCE ,1989 ,115 (12) : 2763-2781. [ 11 ] 　Liu Y W , Moses F. A sequential response surface method and its application in the reliability analysis of aircraft structural system[J ] . St r uct ural Safety ,1994 ,16 (1) :39-46. [12 ] 　Rajashekhar M R , Ellingwood B R. A new look at the response safe approach for reliability analysis [J ] . St r uct ural Safety ,1993 ,12 (3) :205-220. [13 ] 　Bayer V , Schueller G I. Discussion on : A new look at the response surface approach for reliability [J ] . St r uct ural Safety ,1994 ,16 (3) :227-234. [14 ] 　Kim S - H , Na S - W. Response surface method using vector projected sampling points [J ] . St r uct ural Safety ,1997 ,19 (1) :3-19. [15 ] 　Bohm F , Bruckner- Foit A. On criterion for accepting a response surface model [ J ] . Probabilistic Engi neeri ng Mecha nics ,1992 ,7 (3) :183-190. [16 ] 　Das P K, Zheng Y. Cumulative formation of response surface and its use in reliability analysis [ J ] . Probabilistic Engi neeri ng Mecha nics ,2000 ,15 (4) :309-315. [17 ] 　Maymon G. Probability of failure of structural without a closed form failure function [ J ] . Comp u St r uct ,1993 ,49 (2) :301-313. [18 ] 　Kaymaz I , McMahon C A. A response surface method based on weighted regression for structural re2 liability analysis [J ] . Probabilistic Engi neeri ng Mecha nics ,2005 ,20 (1) :11-17. [19 ] 　Guan X L , Melchers R E. Effect of response surface parameter variation on structural reliability esti2 mates [J ] . St r uct ural Safety ,2001 ,23 (4) :429-444. [20 ] 　Gupta S , Manohar C S. An improved response surface method for the determination of failure proba2 bility and importance measures [J ] . St r uct ural Safety ,2004 ,26 (2) :123-139. [21 ] 　Guan X L , Melchers R E. Multitangent-plane surface method for reliability calculation [ J ] . J Eng Mech ,1997 ,123 (10) :996-1002. [22 ] 　Mahadevan S , Shi P. Multiple linearization method for nonlinear reliability analysis [J ] . J Eng Mech 32机械可靠性分析的高精度响应面法 ASCE ,2000 ,127 (11) :1165-1173. [23 ] 　杨成永 ,张弥 ,白小亮. 用于结构可靠性分析的多响应面法 [J ] . 北方交通大学学报 , 2001 ,25 (1) :1- 4. [24 ] 　Deng J , Gu D S , Li X B , et al . Structural reliability analysis for implicit performance function using artificial neural neural network[J ] . St r uct ural Safety ,2005 ,27 (1) :25-48. Adva nce d Resp o ns e S urf ace Me t hod f or Mec ha nical Relia bili t y Analysis L Β Zhen-zhou , 　ZHAO J ie , 　YU E Zhu-feng ( School of Aerona utics , Northwester n Polytechnical Univers ity , Xi’a n 710072 , P. R. Chi na) Abs t ract : Based on the classical response surface method (RSM) , a novel RSM using improved ex2 perimental points ( EPs) is presented for reliability analysis . Two novel points are included in the pre2 sented method. One is the use of linear interpolation , from which the total EPs for determining the RS are selected to be closer to the actual failure surface . The other is the application of sequential linear interpolation to control the distance between the surrounding EPs and the center EP , by which the presented method can ensure that the RS fits the actual failure surface in the region of maximum like2 lihood as the center EPs converging to the actual most probable point (MPP) . Since the fitting preci2 sion of the RS to the actual failure surface in the vicinity of the MPP , which has significant contribu2 tion to the probability of the failure surface being exceeded , is increased by the presented method , the precision of the failure probability calculated by RS is increased as well . Numerical examples illustrate the accuracy and efficiency of the presented method. Key words : response surface method ; implicit limit state ; failure probability 42 吕 　震 　宙 　　　赵 　　洁 　　　岳 　珠 　峰