工程结构可靠度计算的Matlab实现

　 计算机应用 　

工程结构可靠度计算的Matlab实现

章慧健

(西南交通大学土木工程学院,四川成都610031)

　　【摘　要】　从基本原理和相关算例分析两方面,阐述利用功能极其强大的科学和工程计算数学软件系统Matlab的优化工具箱和随机数发生器实现工程结构可靠度的多种计算方法。　　【关键词】　结构可靠度;　优化工具箱;　蒙特卡罗法;　随机变量;　Matlab软件　　【中图分类号】　TP319　　　　　　　　　　　　【文献标识码】　B

　　工程结构可靠度计算除了如何确定极限状态功能函数

和各变量的概型分布以外,很大程度上是一种数学计算(本文不讨论如何确定功能函数和变量的概型分布,只对已知功能函数如何计算可靠度指标β的情况作一论述)。Matlab是一种功能极其强大的科学和工程计算数学软件系统,大量数学、统计、科学和工程所需的函数、C程语言相比,Matlab运算功能强等特点,可以大大提高编程效率。

2　基于Matlab优化工具箱的最优化方法

211　,1,,n是n,,,见下例),由这些随机变量表示的结构极限状:Z=g(X1,X2,…,Xn)=0,根据可靠度指标β的几何意义(在标准正态坐标系中,原点到极限状态曲面的最短距离)可知:

β=

1　　　可靠度计算常用方法有解析法和模拟法。

解析法中普遍采用的是一次二阶矩法,包括中心点法、验算点法(JC法)、映射变换法、实用分析法等。用Matlab实现这些方法的编程计算是相当方便快捷的,因为它有很多现成可用的计算程序,如求反函数、概率分布函数,概率密度函数,函数求导等等,而且采用矩阵运算可以避免循环结构。与同样方法的C语言实现,Matlab的实现语句要少的多。笔者从可靠度指标β的几何意义出发,建立求解可靠指标的优化模型,利用Matlab优化工具箱计算可靠度指标和验算点。

模拟法这里主要指蒙特卡罗法。蒙特卡罗法在目前可靠度计算中,被认为是一种相对精确法,它不受极限状态方程非线性的限制。对那些大型复杂结构,相应的极限状态功能函数往往是非线性的,如果采用解析法中求导的方法计算可靠度,当考虑因素较多时,极限状态功能函数将变得难以处理。蒙特卡罗法回避了可靠度分析中的数学困难,不需要考虑极限状态曲面的复杂性。此外,若基本变量相关,可利用条件概率密度,把多维问题化为一维问题来解决。因此从理论上来说,该方法的应用几乎没有什么限制。但在实际问题中,变量连续型分布是很复杂的。有的只能给分布函数的解析表达式,但给不出其反函数的解析表达式,如著名的β分布;有的则连分布函数的解析表达式都给不出。所以通常情况下对连续型分布采用直接抽样是有一定困难的。笔者利用Matlab的强大数值计算功能,实现了在Matlab中采用蒙特卡罗直接抽样计算结构可靠度,较好地解决了上述问题。将Matlab用于蒙特卡罗法的一个显著优点是它拥有功能强大的随机数发生器指令。

[1]

∑[(X

3

i

σXi]-μXi)/

2

开始时验算点未知,把β看成极限状态曲面上点P(X1,

X2,…,Xn)的函数,通过优化求解,找到β最小值,即可得可

333

靠度指标β和验算点P3(X1,X2,…,Xn)。求解可靠度指标

可以归结为以下约束优化模型:

β=min　

2

∑[(X

i=1

n

3

i

σXi]-μXi)/

2

333

S1t1　Z=g(X1,X2,…,Xn)=0

如果极限状态方程中某个变量(Xj)可用其他变量表示,即:

Xj=h(X1,X2,…,Xj-1,Xj+1,…,Xn),则上述模型可转化

为无约束优化模型:

minβ=

3

3

3

2

i=1,i≠j

∑[(X

n

3

i

33

σXi]2+{[h(X1-μ,X2,…,Xi)/

σXj}2Xj-1,Xj+1,…,Xn)-μXj]/

212　算　例

例:已知非线形极限状态方程Z=g(f,r,H)=567fr-2

015H=0,f服从正态分布,μ服从正f=016,σf=010786;r态分布,μr=2118,σr=010654;H服从对数正态分布,μH=

3218,σH=01984。f,r,H相互独立,求可靠度指标β及验算

[2]点(f3,r3,H3)。

μ解:先将H当量正态化:(变异系数为δH=σH/H=

0103)

[收稿日期]2006-05-25[作者简介]章慧健,硕士研究生。

154

　 计算

H服从对数正态分布,则h=lnH服从正态分布,且

机应用 　

μ′=H

μ+2

σH=3149,　′=2

ln(1+δH)=0103

他两个变量表示

CC=sqrt(((x(1)-016)/010786)^2+((x(2)-2118)/010654)^2+…((h-3149)/0103)^2);%β表达式

运行如下:

〉〉noncon

可靠度指标为bata=11964254

验算点为[01456167,21158971,331418877]

(1)采用约束优化法:

求解该问题可靠度指标可归结为如下约束优化数学模型:

23

σf]2+[(r3-μr)/σr]2+[(h3-minβ=[(f-μf)/

2

μ)/σH′′]H

s1t1　567fr-015HMatlab源程序如下:

3332

=567fr-015e

332h3

=0

3　在Matlab中实现蒙特卡罗法

311　基本原理

functionrelia　%定义主函数

x0=[016,2118,3149];%初始迭代点,这里取均值%调用优化工具箱求解

通常在用蒙特卡罗直接抽样法时,必须解X=F-1(r)以求得服从相应分布类型的随机变量X,再代入功能函数求

解。而Matlab(610版)提供了23种随机变量分布类型的随机数发生器(如正态分布、对数正态分布、泊松分布、威布尔分布等),可直接产生变量X,省去了可能会,极大地提高了效

2

率。r=norrnd(N(u,s)分布的m

);options=optimset(′LargeScale′,′off′