Monte-Carlo算法在土木工程中的应用.docx

 PAGE \* MERGEFORMAT 3 Matlab 程序设计： Monte-Carlo 算法在土木工程中的应用 主讲人： 王佐才 引言 Monte-Carlo算法是一个非常有用的算法，在许多实际问题中，都有用武之地。目前，这个算法已经在金融学、经济学、工程学、物理学、计算科学及计算机科学等多个领域广泛应用。而且这个算法本身并不复杂，只要掌握概率论及数理统计的基本知识，就可以学会并加以应用。由于这种算法与传统的确定性算法在解决问题的思路方面截然不同，作为土木工程相关学生，掌握此算法，可以开阔思维，为解决问题增加一条新的思路。 Monte-Carlo算法引导 首先，我们来看一个有意思的问题：在一个1平方米的正方形木板上，随意画一个圈，求这个圈的面积。我们知道，如果圆圈是标准的，我们可以通过测量半径r，然后用 S = pi * r^2 来求出面积。可是，我们画的圈一般是不标准的，有时还特别不???则，如下图是我画的巨难看的圆圈。 显然，这个图形不太可能有面积公式可以套用，也不太可能用解析的方法给出准确解。不过，我们可以用如下方法求这个图形的面积： 假设我手里有很多均匀的豆子，我用n颗豆子铺满整个正方形，其中m颗豆子落在圈内，则圈的面积等于m/n*a^2。（解释“飞镖法”的问题） 以上这个过程，就是Monte-Carlo算法直观应用例子。非形式化地说，Monte-Carlo算法泛指一类算法。在这些算法中，要求解的问题是某随机事件的概率或某随机变量的期望。这时，通过“实验”方法，用频率代替概率或得到随机变量的某些数字特征，以此作为问题的解。 Monte-Carlo数理基础 总体来说，适用于Monte-Carlo算法的问题，比较常见的有两类。一类是问题的解等价于某事件概率，如上述求不规则图形面积的问题；另一类是判定问题，即判定某个命题是否为真。 我们主要考虑第一类。对于这类问题，通常的方法是通过大量重复性实验，用事件发生的频率估计概率。之所以能这样做的数学基础，是伯努利大数法则：事件发生的频率依概率收敛于事件的概率p。这个法则从数学生严格描述了频率的稳定性，直观意义就是当实验次数很大时，频率与概率偏差很大的概率非常小。此类问题的误差分析比较繁杂，此处从略。有兴趣的同学可以参考相关资料。 利用Matlab编程时需要用的主要命令 生成均匀随机数的matlab命令： rand(m,n) 生成正态分布随机数的matlab命令：randn(m,n) 寻找满足条件的数据的matlab命令：find 计算数组长度的matlab命令： length 应用实例一：求曲线y=x^2与曲线x=y^2所围成的面积。 %% Monte Carlo simulation x=rand(1,100000); y=rand(1,100000); Y=sqrt(x); Z=x.^2; k=find(y=YyZ); area=length(k)/length(y) 应用实例二： 有一长度为L的简支梁，荷载p1和p2随机的落在梁上，求p1到左端点的距离小于p1与p2之间距离的概率是多少。 % order L=10; L=10; x=rand(1,1000000)*L; y=rand(1,1000000)*L; k=find(abs(y-x)x); p=length(k)/length(x) 应用实例三： 在结构设计中，结构抵抗外部荷载的能力称为抗力，用R表示，外部荷载用L表示。一般而言，R与L均为随机变量。现假设R的均值为u_R=500，标准差s_R=50.; L的均值u_L=350, s_L=70. 且R和L均服从正态分布。利用Mente –Carlo 样本方法编写Matlab程序求解： （1）Ru_R+s_R的概率； （2）定义g=R-L，求g0的概率。 %Matlab solution: u_R=500; s_R=50; N=100000; R=randn(1,N)*s_R+u_R; K1=find(Ru_R+s_R); p1=length(k1)/N; u_L=350; s_L=70; L=randn(1,N)*s_L+u_L; g=R-L; k2=find(g0); p2=length(k2)/N; 作业 结构设计中的外荷载R往往由结构恒载D，活载L，以及其它荷载S，即R=D+L+S。 现有一结构的恒载和荷载均服从正态分布，且均值分别为500kN， 300 kN，变异系数（COV）分别为10%与20%。其它荷载S满足对数正态分布，且均值为100，变异系数（COV）为25%。 求R的均值与方差。