蒙特卡罗法用于结构抗震可靠度分析

2012年11月结构可靠度课程论文结构可靠度课程论文蒙特卡罗法用于结构抗震可靠度分析班 级： 土研12级1班 任课教师： 赵东拂 学 号： 1108521312023 学生姓名： 崔相东 日 期： 2012年11月10日 蒙特卡罗法用于结构抗震可靠度分析崔相东（北京建筑工程学院）摘要：地震作用及结构的地震反应存在非常大的不确定性，因此要用传统的解析方法来求出土木工程结构在地震作用的可靠度的过程是十分复杂，甚至无法实现的，而蒙特卡罗法（Monte Carlo）作为一种基于随机实验的模拟方法却能很好的实现这一过程，它可以模拟各种地震作用引起的破坏形式及结构反应，并将承载能力的不确定性纳入随机实验过程之中，很好地解决了常规地震危险性分析中存在的各种问题。下面就将（Monte Carlo）的基本思想做一简单介绍，并结合算例对其计算方法进行阐述。关键词： 地震作用、地震反应的不确定性；蒙特卡罗法Monte Carlo method for analysis of structural seismic reliabilityXiangDong cui (Beijing university of civil engineering and architecture)Abstract: there is extremely uncertainty in the earthquake function and the seismic reaction of structure, therefore, the process to find out reliability of civil engineering structure under seismic action via the traditional analytic method is very complex, and even be impossible, but the Monte Carlo method, as a kind of simulation method based on random experiment can fulfill the task, which can simulate various the failure form and structure response caused by earthquake effect, and evaluate the uncertainty of the bearing capacity into the experimental process, so solve all kinds of problems of the conventional seismic hazard analysis. The following content will introduce the basic thought of Monte Carlo method, combining with an example for its calculation methods.Keywords: uncertainty in the earthquake function and performance; Monte Carlo method引言：作为基础设施，工程结构不仅关系到国计民生，还会影响到一个国家的现代化进程，因此，保证结构在规定的使用期内能够承受各种作用，满足设计要求的各项使用功能，及具有不需过多维护而能保持其自身工作性能的能力是到至关重要的，即要保证结构的安全性、适用性和耐久性，这三个方面构成了工程结构可靠性的基本内容。为保证结构的可靠性，首先要研究建造结构所使用材料的各项力学性能，结构上各种作用的特性，结构的内力分析方具有明确的传力路径；精心施工，严格按照施工规程进行操作；正常使用，按设计要求使用结构并进行正常维护。然而，即便如此，也不能保证结构绝对的安全可靠，这是因为在结构的设计、建造和使用过程中，还存在着种种影响结构可靠性的不确定性，合理、正常的设计、施工和使用只是保证结构具有一定可靠性的前提和基本条件1。这些不确定性中最危险的因素就是地震作用，按Algermissen等人编的美国地震危险图估算, 地面运动峰加速度变异系数达1.38。而其他荷载的变异系数, 不论是我国还是美国的统计都在0.07 至0.37之间。此外,地震作用下结构往往处于非线性非弹性阶段, 破坏机制极复杂, 结构反应及承载能力的确定性比其它荷载作用时就更加突出2。综上所述，地震作用及结构的地震反应存在非常大的不确定性，因此要用传统的解析方法来求出土木工程结构在地震作用的可靠度的过程是十分复杂，甚至无法实现的，而蒙特卡罗法（Monte Carlo）作为一种基于随机实验的模拟方法却能很好的实现这一过程，它可以模拟各种地震作用引起的破坏形式及结构反应，并将承载能力的不确定性纳入随机实验过程之中，很好地解决了常规地震危险性分析中存在的各种问题3。下面将就蒙特卡罗法（Monte Carlo）的基本思想做一简单介绍，并结合简单算例对其计算方法进行阐述。1、 蒙特卡罗法介绍1.1 蒙特卡罗法基本思想蒙特卡罗的基本思想可概括为：为求研究问题的概率解, 构造一个表示所研究问题概率解的数学模型(计算模型)，记为:依据计算模型中各随机变量Xi所服从的分布进行随机抽样，并按计算模型计算Y的多个估计值，最终用频率统计法求出Y的概率解。方法的核心是随机抽样，而随机抽样的关键在于产生0,1区间上均匀分布的随机数。服从其它分布的随机数一般可通过该随机数变换得到4。1.2 蒙特卡罗法求解的基本过程由基本思想可知，求解过程大致分为四步：分析并拟定给定问题中的随机变量，构造表示给定问题概率解的数学模型；对模型中的随机变量X1, X2,Xn各进行L次随机抽样，获得L组抽样值：X1k, X2k,，Xnk (k=1,2, ,L)；把L组抽样值代入计算模型，求出随机变量Y的L个估计值Y1 , Y2 , , YL； 利用频率统计法，由Y的估计值Y1 , Y2 , , YL求出描述Y分布特征的分布曲线(概率解) 如图1-14。图1-1 随机变量Y的分布曲线1.3 随机数的概念随机数是随机变量的观测值，由其构成的数据序列叫做随机数序列，它是一个无周期的数据序列。 蒙特卡罗法需要对变量进行数以千计、万计、甚至是百万计的抽样，这在实现过程中几乎是不可能的。因此考虑用计算机模仿实际抽样过程，形成一个有周期的抽样数据序列。称这种数据序列为“伪随机数”序列，其中的元素叫做“伪随机数”。 伪随机数显然不是真正意义下的随机数，即这种抽样值并非是随机变量的真实观测值。尽管如此，只要对伪随机数序列进行一系列严格的统计检验，证明它可以满足统计的要求，则伪随机数就可以作为真随机数使用。 为了满足抽样问题的需要，在计算机上产生的伪随机数序列不仅要有足够长的周期，而且应当具有符合要求的概率统计性质。 从理论上讲，只要有一种连续分布的随机数，就可以采用数学变换的方法产生其它分布的随机数。0,1区间上均匀分布的随机变量的抽样值是最简单、最基本的一种连续分布的随机数，其它分布的随机数都可以借助它来产生，所以说：0,1均匀分布随机数技术是实现随机抽样的最基本工具4。 1.4 伪随机数的产生方法1.4.1. 乘同余法 该方法产生伪随机数序列的递推同余式为： Xn,Xn+1第n次和第n+1次产生的伪随机数；-乘子系数；M- 模；rn+1- 0,1区间上的伪随机数。Xn+1Xn (mod M)叫做以 M 为模的同余式，表示xn+1取值为：与Xn的积除以M的余数部分。1.4.2. 混合同余法该方法产生伪随机数序列的递推同余式如下： 混合同余法比乘同余法仅是增加了一个增量,其它含义与乘同余法相同。 如：M =219 = 524288 , = 55 = 3125 时，X0= 23 , 11 ,19 , 37;= 3,7,11,17分4套配合使用，混合同余法可产生周期为524288伪随机数序列。注:所获得的伪随机数是否能代表真正意义上的随机数,还需进行检验(检验方法略)。1.5 随机变量的抽样对随机变量的抽样，有经验分布函数抽样法、直接抽样法和变换抽样法等。在此仅介绍经验分布函数抽样法。1.5.1 随机变量的经验分布函数在概率论中，随机变量X的分布函数是随机变量X的取值不大于实数 x的概率。通常记为： F(x)=P(Xx) 经验分布函数是由X的 n 个观测值X1, X2,Xn,用统计方法得到的分布函数，记为Fn(x)。以统计所得的Fn(x)代替F (x), 并记4：1.5.2 随机变量经验分布函数的构造方法频率统计法构造经验分布函数的条件及步骤:(1)使用条件： 随机变量必须有观测值, 且个数足够多(如30)。 (2)构造分布函数的步骤确定频率统计区间数 要将n个观测值所在的大区间划分为k个小区间 ,一般n/k3 ，且k最好为奇数。 计算区间端值 将(xmin,xmax) k等分, 各小区间的k+1个端值为：其中xmin,xmax是观测值的最大和最小值。 经验分布函数 记ni 为观测值落入区间(xi , xi+1)内的频数(个数),fi为累加频率，则 其中n=n1+n2+nk 显然 f1=1.0,fk=nk/n,f1 f2fk 1.5.3 经验分布函数的抽样将坐标原点设为(xmin ,0)，若已知0,1区间上均匀分布的随机数ri ，则在图1-2纵轴上可确定点(xmin,ri), 过该点作横轴的平行线交分布曲线于点(xi,ri ), xi 则是对应于随机数 ri的一次随机抽样值。由此可得随机变量的一系列抽样值。称该抽样方法为经验分布函数抽样法4。图1-2 分布函数抽样过程示意图2、 蒙特卡罗法算例如图所示为一双杆受力体系，两杆的轴向承载能力分别为、，竖向荷载为P，、均服从对数正态分布，其平均值和标准差分别为；，P服从极值型分布，平均值和标准差分别为。试用直接重要抽样法估计该体系的失效概率。 图2-1 二杆受力体系 解 由于每一个杆失效都会导致体系失效，该体系为串联体系，由图可建立体系的两个功能函数 借助于一次二阶矩方法求得相应于的可靠指标，验算点坐标，;相应于的可靠指标，验算点坐标，;从而得 如果取总的抽样次数，以的验算点为抽样中心的抽样次数为 以的验算点为抽样中心的抽样次数为 、和P的概率密度函数分别为 在本例中分别用正态分布对、和P进行抽样，当以的验算点，为抽样中心时，抽样概率密度函数为 当以的验算点，为抽样中心时，抽样概率密度函数为 首先以的验算点，为抽样中心时，对、和P抽样。产生3个随机数0.0794、0.3157、0.7315，由反函数法求得、和P的样本值，。从而 由于，因而 ， 得 重复上述过程次，求得 同理。再以的验算点，为抽样中心时，对、和P抽样，模拟次，求得，最终求得体系失效概率的估计值为 此计算过程比较繁琐，人工手算太浪费时间，下面可以用Matlab进行计算，使计算大为简便。Matlab编写程序为：m =2;muX = 85;75;35;cvX = 0.15;0.15;0.3;sigmaX = cvX.* muX;sLn = sqrt(log(1+cvX(1:2).2);mLn = log(muX(1:2) - sLn.2/2;aEv = sqrt(6)* sigmaX(3) /pi;uEv = -psi(1)*aEv-muX(3);muX1= muX; sigmaX1 = sigmaX;a = 1,0,-1/sqrt(2);0,1,-1/sqrt(2);for k = 1 : m x = muX; normX = eps; while abs(norm(x) - normX)/normX 1e-6 normX = norm(x); g = a(k , :)* x; gX = a(k , :); cdfX = logncdf(x(1:2),mLn,sLn);1-evcdf(-x(3),uEv,aEv); pdfX = lognpdf(x(1:2),mLn,sLn);evpdf(-x(3),uEv,aEv); nc = norminv(cdfX); sigmaX1 = normpdf(nc)./pdfX; muX1 = x - nc.*sigmaX1; gs = gX.*sigmaX1;alphaX =-gs/norm(gs); bbeta = (g+gX*(muX1 - x)/norm(gs); x = muX1 + bbeta * sigmaX1.* alphaX; end xD( : , k) = x; pFL(k) = normcdf(-bbeta);endnS = 1e6 ;n = length(muX);nSi = ceil(pFL/sum(pFL)* nS);for l=1:m for k = 1:n, v(:,k)=normrnd(xD(k,l),sigmaX(k),nSi(l),1);end t = v(:,3)/sqrt(2);g1 = v(:,1)-t; g2 = v(:,2)-t; g = min(g1,g2,2); fp = lognpdf(v(:,1),mLn(1),sLn(1).*lognpdf(v(:,2),mLn(2).*evpdf(-v(:,3),uEv,aEv); pvi(1:nSi(1),1:m) = 1; for k1 =1:m for k=1:n pvi(:,k1) = pvi(:,k1).*normpdf(v(:,k),xD(k,k1),sigmaX(k); end endfp = fp./sum(pvi,2);clear v pvi;nFi = sum(fp(g0);pFi(l) = nFi/nSi(l);endpF = sum(pFi) 程序运行结果： ，和手算值相差级误差，可作为精确值。 参考文献1 贡金鑫. 工程结构可靠度计算方法. 大连理工大学出版社，2003. 2 董伟民，鲍布斌.地震作用下的结构可靠度. 北京：机械工业部设计总院、中国建筑科学研究院，2008.3 符圣聪，江静贝，黄世敏. 蒙特卡罗法用于地震危险性和液化势的估计. 北京: 中国建筑科学研究院工程抗震研究所，2009.4 张明.结构可靠度分析方法与程序.科学出版社，2009.5 GB50011-2001, 建筑结构抗震设计规范. 中国建筑工业出版社, 北京, 2002.9