基于MATLAB的可靠度分析算法的性能比较.pdf

2014年4月 ? 第40卷 总第178期 2014年 第2期 Sichuan Building Materials DOI:10. 3969/ j. issn. 1672 -4011. 2014. 02. 031 基于 MATLAB 的可靠度分析算法的性能比较 殷 杰 (中国建筑西南设计研究院有限公司, 四川 成都 610081) 作者简介:殷杰(1984 - ),男,江西九江人,硕士研究生,工程师,主要 研究方向:结构工程。 摘 要: 简要介绍和讨论可靠度分析算法原理; 基于 结构可靠度分析方法: 一次二阶矩法、 JC 法、 蒙特卡罗法、 重要抽样法、 响应面法的基本原理, 利用 MATLAB 软件对 两个算例编制了相应计算程序, 并对计算结果进行了比较。 计算结果的比较表明, 利用 MATLAB 编制的程序进行结构 可靠度计算, 可节约大量运算时间, 并能得到较为有效的、 合理的结果。 关键词: 结构可靠度; MATLAB; 一次二阶矩法; JC 法; 蒙特卡罗法; 重要抽样法; 响应面法 中图分类号: TU317 + 1 文献标志码: B 文章编号:1672 -4011(2014)02 -0069 -04 0 前 言 实际工程中结构可靠度的计算方法比较多。 其中, 一 次二阶矩法应用相对较为广泛, 是结构可靠度分析和计算 的基本方法。 该类方法是一种近似的计算方法, 但具有很 强的适用性, 计算精度能够满足工程需求。 以一次二阶矩 法为基础, 发展出了均值一次二阶矩法、 改进的一次二阶 矩法、 JC 法、 几何法。 而蒙特卡洛法一类方法则是对各随机变量进行大量抽 样, 结构失效次数占抽样数的频率即为其失效概率。 不需 考虑功能函数的非线性和极限状态曲面的复杂性, 直观、 精确、 通用性强; 缺点是计算量大、 效率低。 蒙特卡洛法 又包括大量子类方法, 其中重要抽样法运用较多, 重要抽 样法在满足同样精度的条件下, 减少了模拟次数, 大大提 高蒙特卡洛法的效率。 在结构系统中, 结构的随机输人量与输出量之间经常 不存在明确的解析关系式, 结构的功能函数为随机变量的 隐式函数, 这个时候常常利用响应面法, 其思想是选用一 个适当的明确表达的函数或曲面来近似替代一个不能明确 表达的函数或极限状态面, 当响应面在一系列取样点上拟 合之后, 再用构件可靠度分析的方法进行分析1。 MATLAB 是一款基于矩阵的数学软件, 其语法简单、 用户界面友善、 矩阵运算功能强大, 其已经成为数值分析 计算、 数值仿真、 信号处理等高级课程的基本运算工具。 在可靠度指标优化计算过程中需要建立大量与概率、 统计 和最优化方法相关的数值计算方法在 MATLAB 环境中均可 很容易实现2 3 1 结构可靠度的计算方法与原理 可靠度的计算方法从研究的对象来说, 可分为点可靠 度计算方法和体系可靠度计算方法。 方法主要有一次二阶 矩法, 均值一次二阶矩法, 改进一次二阶矩法, JC 法, Monte Carlo 法, Monte Carlo 重要抽样法、 响应面法等等。 而可靠度计算主要是计算结构可靠指标。 1 1 一次二阶矩法 一次二阶矩法4 5在实际工程中应用相当广泛, 其要 点是非正态随机变量的正态变换及非线性功能函数的线性 化由于将非线性功能函数作了线性化处理。 均值一次二阶 矩法、 改进的一次二阶矩法、 JC 法都是以一次二阶矩法为 基础的可靠度计算方法。 1 1 1 均值一次二阶矩法 设结构构件的功能函数为: Z = g(X) = g(X1,Xn)(1) 式中,Xi(i = 1,2,3,n) 均服从正态分布,且相互独 立。 将功能函数在均值点 xi(i = 1,n) 展开 Taylor 级数, 仅保留线性项, 得到: Z g(X1- Xn) + n i =1 g Xi X (Xi- Xi)(2) 并得到 Z 的均值、 均方差为: Z= g(X1,Xn)(3) Z= n i =1 g Xi X Xi () 2 1 2 (4) 结构构件可靠度指标表示为: = Z Z = g(X1,Xn) n i =1 g Xi X Xi () 2 1 2 (5) 该方法对于非线性功能函数, 因略去二阶及更高阶项, 误差将随着线性化点到失效边界距离的增大而增大, 而均 值法中所选用的均值点一般在可靠区而不在失效边界上, 结果往往带来相当大的误差, 同时选用不同的极限状态方 程不能得到相同的可靠指标, 此为该方法的严重问题。 1 1 2 改进一次二阶矩法 针对均值一次二阶矩法将功能函数线性化点去做基本 随机变量均值点带来的问题, 改进一次二阶矩法把线性化 点选在失效边界上, 且选在与结构最大可能失效概率对应 的设计验算点上, 以克服均值一次二阶矩法存在的问题。 该方法无疑优于均值一次二阶矩法, 为工程实际可靠度计 算中求解的基础。 但该方法只是在随机变量统计独立、 正 态分布和线性极限状态方程才是精确的, 否则只能得到近 似的结果。 设结构构件的功能函数为: Z = g(X) = g(X1,Xn)(6) 式中, Xi(i = 1,2,3,n) 均服从正态分布, 且相互独 96 2014年4月 ? 第40卷 总第178期 2014年 第2期 Sichuan Building Materials 立。 将功能函数在设计验算点 p(x 1 ,x n) 展开 Taylor 级数, 仅保留线性项, 得到: Z g(x 1 ,x n) + n i =1 g Xi p (Xi- x i )(7) 可得到 Z 的均值和均方差为: Z= g(x 1 ,x n) + n i =1 g Xi p (Xi- x i )(8) Z= n i =1 g Xi p Xi () 2 1 2 (9) 结构构件的可靠度指标可表示为: = Z Z = Z= g(x 1 ,x n) + n i =1 g Xi p (Xi- x i ) Z= n i =1 g Xi p Xi () 2 1 2 (10) 设计验算点坐标: x i = Xi+ XicosXi(11) 定义变量 Xi灵明度系数: cosXi= g Xi p Xi n i =1 g Xi p Xi () 2 1 2 (12) 考虑到设计验算点 p应位于极限状态曲面上, 故 g(x 1 ,x n) = 0 , 所以可靠度指标可进一步表示为: = Z Z = n i =1 g Xi p (Xi- x i ) Z= n i =1 g Xi p Xi () 2 1 2 (13) 计算步骤: 1) 假定 x i (一般可采用 x i = Xi); 2) 采用式 (12), 求 cosXi; 3) 采用式 (10), 求 ; 4) 采用式 (11), 求新的 x i ; 5) 用新的 x i 重复步骤(2) - (4),直到最后两次的 值之差小于允许值。 1 1 3 JC 法 由于工程结构各随机变量的非正态性, RudigerRackwitz 和 Bernd Fiessler 提出了 JC 法。 其原理为: 首先把随机变量 原来的非正态分布用正态分布代替, 但对所代替的正态分 布要求, 在设计验算点处的概率分布 CDF 和概率密度函数 PDF 与原分布的 CDF 和 PDF 相同, 然后根据这两个条件求 得等效正态分布的均值和标准差, 这个做法即为当量正态 化转换, 在等效正态分布的均值和标准差确定之后, 求解 结构可靠度指标的过程与改进一次二阶矩法大致相同。 该 方法克服上述两方法的不足, 适用于随机变量为任意分布 下结构可靠指标的求解, 运算简捷, 对线性程度不高的结 构功能函数, 其精度能满足工程实际需要, 并已为国际联 合委员会 (JCSS) 所采用4。 上文提到的当量正态化转换所需条件可表述如下: 设计验算点坐标 xi (i = 1,n) 处非正态变量与其等 效正态变量的分布函数值及分布密度函数值相等, 即: FXi(x i ) = x i - Xi Xi () (14) fXi(x i ) = x i - Xi Xi () Xi (15) 并由此求的等效变量的均值和均方差为: Xi= x i - -1FXi(x i )Xi(16) Xi= -1FXi(x i ) fXi(x i ) (17) JC 法的计算步骤4: 1) 假定初始验算点 x i (一般可采用 x i = Xi); 2) 对非正态变量 Xi, 根据 (16) 和 (17) 求等效正 态变量的 Xi, Xi, 并用 Xi, Xi替代 Xi, Xi; 3) 采用式 (12), 求 cosXi; 4) 采用式 (10), 求 ; 5) 采用式 (11), 求新的 x i ; 6) 用新的 x i 重复步骤(2) - (5),直到最后两次的 值之差小于允许值。 1 2 Monte Carlo 法 蒙特卡罗法5 7, 即统计试验法。 蒙特卡罗法求解结 构失效概率的基本思路是: 先对影响其可靠度的随机变量 进行大量随机抽样, 然后把这些抽样值分别代入功能函数 式, 确定结构的失效与否, 最后求得结构的失效概率, 失效 概率即结构失效次数占总抽样数的频率。 蒙特卡罗法的优 点在于其精度随着 N 的次数增加而渐次提高。 若 N 值选取 足够大时, 则可得到相对精确值。 设结构构件的功能函数为: Z = g(X1,Xn)(18) 式中, Xi(i = 1,2,3,n) 具有任意分布的随机变量。 对 X 进行 N 次随机抽样, 得到 N 组 (j = 1, , N), 将第 j 组 (j =1, , N) 的值 Xji(i = 1,n) 代入功能函 数, 得到 N 个 Zj(j = 1,N) 。 设在 N 个 Zj值中存在 Nf 个 Zj 0 , 则结构的失效概率为: Pf Nf N (19) 1 3 重要抽样的 Monte Carlo 法 将重要抽样方法运用到蒙特卡洛法中, 在满足同样精 度的条件下, 减少了模拟次数, 提高了效率。 当选择了合 适的重要抽样概率密度函数, 便可以有效减小方差, 大大 提高蒙特卡洛法的效率5。 设: Ig(X1,Xn) = 0,g(X1,Xn) 0 1,g(X1,Xn) 0 (20) 结构构件失效概率为: Pf= Pg(X) 0 (23) 所以有: 07 2014年4月 ? 第40卷 总第178期 2014年 第2期 Sichuan Building Materials Pf= Ig(x) fX(x) hX(x)hX(x)dx (24) 令 Y = Ig(x) fX(x) hX(x) (25) Pf可视为随机变量 Y 的数学期望 Y, 其无偏估计表示 为: Y= Pf P f = 1 N n i =1 Ig(xi) fX(xi) hX(xi) (26) P f的方差为: P f 1 N - 1 1 N n i =1 Ig(xi) fX(xi) hX(xi) - P f 2 (27) 1 4 响应面法 响应面法8 10的思想是选用一个适当的明确表达的函 数或曲面来近似替代一个不能明确表达的函数或极限状态 面, 对结构可靠度分析来说, 就是通过尽可能少的一系列有 限元计算来拟合一个响应面以替代未知的极限状态曲面。 当响应面在一系列取样点上拟合之后, 就可以利用构件可 靠度分析的方法对其进行分析。 响应面的形式一般取二次多项式: Z = g(X) Zr= a + n i =1 biXi+ n i =1 ciX2 i + 1i jnd ijXiXj (28) 采用时可以忽略交叉项得到非完全二次多项式: Zr= g (X) = a + n i =1 biXi+ n i =1 ciX2 i (29) 首先确定设计 X 的一系列试验点 xi, 并计算结构相应 的一系列功能函数值 zi。 可采用二水平因子设计在沿均值 X附近选取试验点, 在得到要展开的点 x 之后, 用得到的 响应面函数找到验算点的估计值 x, 下一步新的展开点 x 通过线性插值得到5: x = X+ g(X) g(X) - g(x)(x - X)(30) 在迭代中利用 x 形成求解一下线性方程组: A = g (31) 其中: = (a1,b1,b2,bn,c1,c2,cn)T(32) g = (g 1,g 2,g 2n+1) T (33) A = 1 1 xi1n xi1n xi21n xi21n 1xi1nxi 2 1n + 000 0- fdiagXidiagf2Xi2- 2fxiXi 0fdiagXidiagf2Xi2+ 2fxiXi (34) 利用上述原理进行结构可靠度分析的方法可归纳为如 下的算法步骤 5 : 1) 假定初始迭代点 x = (x1,x2,xn)T, 可取平均值 X; 2) 选取 f 值; 3) 在各个站开点处计算功能函数估计值 g i(i = 1,2, ,2n + 1) , 求出相应的系数矩阵 A; 4) 求解 a,bi,ci(i = 1,2,n) , 解方程 (31); 5) 由 (29) 计算可靠指标 ; 6) 计算在 x处功能函数的估计值; 7) 由 (30) 插值得到新的 x; 8) 重复步骤(3)到(7),直到最后两次的 值之差小于 允许值。 2 算 例 例 14 5已知结构的功能函数为 Z = R - S , R 、 S 均服 从正态分布, R N(20,4) , S N(15,3 6) , 试分别采用 一次二阶距法、 蒙特卡洛法、 重要抽样法、 响应面法估算 结构的失效概率。 例24 5已知结构的功能函数为 Z = R - G - Q , R 、 G 、 Q 均服从对数正态分布、 正态分布、 极值 I 型分布, R= 298 2 、 R= 51 3 、 G= 49 、 G= 3 2 、 Q= 68 、 Q= 19 5 , 试分别采用 JC 法、 蒙特卡洛法、 重要抽样法估算结 构的失效概率。 3 计算结果 表 1算例 1 计算结果 计算方法蒙特卡洛法 一次二阶距法 重要抽样法响应面法 可靠度指标0 92970 92910 92860 9291 验算点 x 值17 237617 2376 失效概率0 17630 176420 17650 17642 抽样次数1e61e6 表 2算例 2 计算结果 计算方法蒙特卡洛法重要抽样法JC 法 可靠度指标3 60183 59393 6014 失效概率1 5700e 0041 6013e 0041 5823e 004 抽样次数1e61e6 为进行算法比较, 这里的标准解以 Monte Carlo 计算值 表示, 由表 1、 表 2 计算结果比较可得, 上述几种方法的计 算结果与 Monte Carlo 计算值都较为接近, 但任然需要注意 的是采用 Monte Carlo 计算出的可靠度指标仍显得很不稳定 (如例 1 中采用 Monte Carlo 法就算出来的可靠度指标就在 0 960 到 0 930 左右变动)。 在例 1 中一次二阶距法、 响应 面法均有很稳定的结果, 同样在例 2, JC 法有很稳定的结 果。 基于 MATLAB 的蒙特卡罗方法, 在分析结构可靠度时 可大大简化计算, 节省大量时间, 且迭代次数越多, 精度越 高。 同时在 MATLAB 中运行时, 重要抽样法、 响应面法的 计算时间均快于蒙特卡罗法, 这说明了这两种方法相对于 蒙特卡罗方法的优势 4 结 论 本文通过一次二阶矩法、 JC 法、 蒙特卡罗法、 重要抽 样法、 响应面法的比较, 在各种基本原理的基础上利用 MATLAB 编制的可靠度分析的程序, 比较两组算例的结果, 评价了各种可靠度计算方法在编程计算方面的优劣, 最后 得出, 利用 MATLAB 编制的程序进行结构可靠度计算, 可 节约大量运算时间, 并能得到较为有效的、 合理的结果。 ID: 001106 (下转第 74 页) 17 2014年4月 ? 第40卷 总第178期 2014年 第2期 Sichuan Building Materials 1 2 2 混合振动控制优化 主动振动控制和被动振动控制同时应用于被控制对象 时, 二者形成的混合振动控制系统要达到良好的控制效果 需对整体进行优化。 这是因为: 单纯的被动最优控制与 单纯的主动最优控制方案简单叠加后组成的混合振动控制 方案不一定是整体最优, 因而要达到整体优化必须寻找一 种新的优化算法; 混合振动控制中的被动耗能器等往往 受材料性能的影响, 达到最优振动抑制效果的混合振动控 制方案, 使得相同技术参数的耗能材料或耗能器应用范围 更加广阔, 将大大推动被动耗能器件的发展; 主动振动 控制设计中, 输入控制能量最小, 是主动振动控制优化目 标的重要组成因素, 被动振动控制结构参数适当微小的变 化将会显著降低原控制结构对输入能量的要求, 并能产生 同样甚至更好的控制效果。 振动控制的优化包括阻尼器和作动器个数、 位置及二 者的比例, 目标函数和约束条件视不同要求而定, 可以是 模态阻尼比、 附加质量、 特定自由度上的响应等。 在阻尼 器和作动器的优化配置中, 有 3 种优化算法得到普遍重视: 遗传算法; 模拟退化算法; 试凑法。 2 结论与展望 2 1 结 论 本文首先介绍了结构振动控制的概念和工作原理; 指 出了高层建筑风振响应的危害及产生机理, 以及对其施加 控制的意义。 然后通过对目前常用的几种控制方法进行对 比, 阐述了混合控制在建筑振动控制中的巨大优势。 并以 目前常用的 AMD - TMD 混合控制为例, 简要介绍了混合控 制中主动控制系统对于风振反应所需施加作用力的计算方 法。 主要得出以下结论: 1) 混合控制方法较主、 被动控制方法有效果好、 造价 低、 易于实现等优点, 将其推广并应用在高层建筑中具有 良好的实际意义。 2) 混合控制对于高层建筑风振响应控制有完整的设计 理论和设计方法, 其中的一些方法(如 AMD - TMD 系统)已 经可以作为一种较成熟的方法应用于高层建筑风振控制中。 2 2 展 望 在高层建筑如雨后春笋般拔地而起的今天, 无论是高 度还是造型上都不断地冲击着人类想象力的极限, 建筑物 抵抗振动的能力也时刻在接受着挑战。 传统提高建筑结构 抵抗振动能力的方法, 实际上是略带盲目性质的。 而振动 控制技术却能够有计划、 有目标的去进行建筑结构的振动 控制, 做到有的放矢。 这样就可以根据控制的对象和控制 方法进行优化, 既能提高建筑的安全水平又能提高经济效 益。 混合控制技术在高层建筑风振控制中的应用还存在许 多难点需要克服, 例如: 1) 由于目前主动控制技术还不成熟, 混合控制的发展 也受到很大的限制。 因此, 在今后的研究当中应加大对主 动控制技术的力度, 克服主动控制技术存在的种种缺陷。 2) 混合控制优化的研究。 在考虑安全性、 经济性、 适 用性的前提下, 如何最佳控制效果? 在风荷载较小时、 较 大时和强烈时如何区别对待等问题也非常值得研究。 3) 混合控制的概念还在进一步的拓展, 智能材料的利 用还处在起步阶段, 如何有效地利用新兴的控制技术、 控 制方法进行组合叠加, 从而创造出更加合理可靠的混合控 制技术, 也是一个非常有价值的研究领域。 由于高层建筑遭受风荷载作用的频率远远大于地震作 用, 因此对于高层建筑的风振反应, 采用振动控制方法进 行削弱具有十分重要的意义。 一些初步研究表明混合控制 的性能大大优于被动控制, 甚至可达到或超过主动控制的 性能, 并在稳定性与适用性方面要优于后者, 因此成为当 前研究的一个热点。 混合控制可以利用各控制方法的优点, 拓宽控制系统的应用范围, 既保证控制效果又降低控制力, 特别是在强烈地震作用下, 混合控制更具有优越性, 因此 混合控制将会得到更广泛的研究与应用。 ID: 001105 参 考 文 献: 1 瞿伟廉,等. 高层建筑和高耸结构的风振控制设计M. 武汉: 武汉测绘科技大学出版社,1991. 2 周福霖. 工程结构减振控制M. 北京:地震出版社,1997. 3 李宏男,李忠献,祁皑,等. 结构振动与控制M. 北京:中国建 筑工业出版社,2005. 4