数学建模-蒙特卡罗方法.pptx

单击此处编辑母版标题样式,单击此处编辑母版文本样式,第二级,第三级,第四级,第五级,*,Click to edit Master title style,Click to edit Master text styles,Second level,Third level,Fourth level,Fifth level,*,美赛培训：蒙特卡罗方法,主讲人：老教练,汇报人：陈蓉蓉 陈雨路 杨园园,第1页,什么叫蒙特卡罗方法？,蒙特卡罗方法又称统计模拟法、随机抽样技术，是一个随机模拟方法，,以概率和统计理论方法为基础一个计算方法，是使用随机数（或伪随机数）来处理很多计算问题方法。,将所求解问题同一定概率模型相联络，用电子计算机实现统计模拟或抽样，以取得问题近似解。为象征性地表明这一方法概率统计特征，故借用赌城蒙特卡罗命名。,第2页,基本思想,当所求问题解是某个事件,概率,，或者是某个随机变量,数学期望,，或者是,与概率，数学期望相关量时,，经过某种试验方法，得出该事件发生,概率,，或者该随机变量若干个详细观察值,算术平均值,，经过它得到问题解。,当随机变量取值仅为,1,或,0,时，它数学期望就是某个事件概率。或者说，某种事件概率也是随机变量（仅取值为,1,或,0,）数学期望,。,第3页,蒙特卡罗方法特点,优点：,1、能够比较逼真地描述含有随机性质事物 特点及物理试验过程,2、受几何条件限制小,3、收敛速度与问题维数无关,4、含有同时计算多个方案与多个未知量能力,5、误差轻易确定,6、程序结构简单，易于实现,第4页,缺点：,1收敛速度慢,2误差含有概率性,3在粒子输运问题中，计算结果与系统大小相关,所以在使用蒙特卡罗方法时，要“扬长避短”，只对问题中难以用解析（或数值）方法处理部分，使用蒙特卡罗方法计算，对那些能用解析（或数值）方法处理部分，应该尽可能使用解析方法,第5页,主要应用范围,粒子输运问题（试验物理，反应堆物理，高能物理）,统计物理,经典数学问题,真空技术,激光技术以及医学,生物,探矿,第6页,什么是随机数？,在连续型随机变量分布中，最简单而且最基本分布是单位均匀分布由该分布抽取简单子样称为随机数序列，其中每一个体称为随机数,符号：,两个特点：独立性，均匀性,第7页,产生随机数,随机数表方法,物理方法,第8页,随机数表,随机数表是由,0,1,2,3,4,5,6,7,8,9,十个数字组成，每个数字以,0.1,等概率出现，数字之间相互独立，这些数字序列叫作随机数字序列。,（假如要得到,n,位有效数字随机数，只需将表中每,n,个相邻随机数字合并在一起，且在最高位前边加上小数点即可。比如，某随机数表第一行数字为,7 6 3 4 2 5 8 9 1.,，要想得到三位有效数字随机数一次为,0.763,，,0.425,，,0.891.,）,第9页,物理方法,利用一些物理现象，在计算机上增加些特殊设备，能够在计算机上直接产生随机数。,作为随机数发生器物理源主要有两种：一个是依据放射性物质放射性，另一个是利用计算机固有噪声。,普通情况下，任意一个随机数在计算机内总是用二进制数表示：,其中 或者为0，或者为1。所以，利用物理方法在计算机产生随机数，就是要产生只取0或1随机数字序列，数字之间相互独立，每个数字取0或1概率均为0.5,第10页,缺点,随机数表需在计算机中占有很大内存，而且也难以满足蒙特卡罗方法对随机数需要量非常大要求，所以，该方法不适于在计算机上使用。,物理方法产生随机数序列无法重复实现，不能进行程序复算。给验证结果带来很大困难。而且增加随机数发生器和电路联接等附加设备，费用昂贵。所以该方法也不适合在计算机上使用。,第11页,伪随机数,用递推公式,产生随机数序列。,第12页,伪随机数存在两个问题,递推公式和初始值确定后，整个随机数序列便被唯一确定。不满足随机数相互独立要求。,因为随机数序列是由递推公式确定，而在计算机上所能表示,0,，,1,上数又是有限，所以，这种方法产生随机数序列就不可能不出现重复。随机数序列出现周期性循环现象。,第13页,处理方案,第一个问题：不能从本质上加以改变，但只要递推公式选比很好，随机数间相互独立性是能够近似满足。,第二个问题：因为用蒙特卡罗方法解任何详细问题时，所使用随机数个数总是有限，只要所用随机数个数不超出伪随机数序列出现循环现象时长度就能够了。,第14页,应用：蒙特卡罗方法计算积分,能够通俗地说，蒙特卡罗方法是用随机试验方法计算积分，即将所要计算积分看作服从某种分布密度函数,f,(,r,),随机变量,(,r,),数学期望,经过某种试验，得到,观察值,r,1,，,r,2,，,，,r,N,（用概率语言来说，从分布密度函数,f,(,r,),中抽取,个子样,r,1,，,r,2,，,，,r,N,，），将对应,个随机变量值,g,(,r,1,),，,g,(,r,2,),，,，,g,(,r,N,),算术平均值,作为积分预计值（近似值）。,第15页,求积分,（,2.1,）,第16页,蒙特卡罗方法步骤以下：,1在区间【a,b】上利用计算机均匀产生n个随机数x1，x2xn，这个能够在MATLAB软件中用unifrnd命令实现。,2、计算每一个随机数对应被积函数值f（x1），f（x2）f（xn）。,3、计算被积函数值平均值,4、所以2.1式值,第17页,简单定积分例子：用蒙特卡罗方法求,首先我们进行严格数学计算，便于后面与蒙特卡洛计算方法所得结果形成对比：已知 原函数是 ，那么定积分值就是：,我们能够在Matlab中输入以下代码进行准确计算：exp（2）-exp（0），这个值是此定积分真实值.,下面进行蒙特卡洛计算上述定积分，其MATLAB代码以下：,N=500,;,x=unifrnd,(,0,2,N1,);,mean,(,2*exp,(3*x.2),第18页,解释语句：,上述三条语句完整实现了蒙特卡洛计算上述定积分步骤,.,第一条语句是设定了停顿条件，共做,N,次,Monte Carlo,模拟,.,第二条语句实现了在积分区间上均匀产生,N,个随机数,.,第三条语句实现蒙特卡洛计算方法面积迫近,.,对,N,设置不一样值，观察所得蒙特卡洛计算方法定积分值，如表,1,所表示，我们能够发觉：当不停增大,N,值时，所得结果越来越靠近真实值,.,第19页,蒙特卡罗计算方法性质：,蒙特卡洛计算方法依据概率统计理论，含有统计特征，主要表现在以下三个方面,.,1.,蒙特卡洛计算方法含有随机性、不确定性,即每次运行结果都会不一样，因为计算机产生随机数并不是能够重现,.,2.,蒙特卡洛计算方法含有统计稳定性。即使每次运行产生随机数是不一样，不过随机数概率分布是一样，所以蒙特卡洛计算方法能够满足概率统计稳定性,.,3.,伴随随机数数量增加，蒙特卡洛计算方法所得结果会愈加迫近真实值，这就是我们所讲依概率收敛到真实值意义,.,第20页,蒙特卡罗方法试验,面积、体积计算问题,冰淇淋锥体积计算,思索题与练习题,第21页,蒙特卡罗方法,随机投点试验求近似解,引例,.,给定曲线,y,=2,x,2,和曲线,y,3,=,x,2,，,曲线交点为,：,P,1,(1,，,1),、,P,2,(1,，,1),。,曲线围成平面有限区域，用蒙特卡罗方法计算区域面积。,P=rand(10000,2);,x=2*P(:,1)-1;,y=2*P(:,2);,II=find(y=x.2);,M=length(II);,S=4*M/10000,plot(x(II),y(II),g.),S=2.1136,第22页,例,5.14,计算,其中,D,为,y=x,2,与,y,2,=,x,所围,D,边界曲线交点为,：,(1,，,1),，,(4,，,2),，,被积函数在求积区域内最大值为,16,。,积分值是三维体积，该三维图形位于立方体区域,0,x,4,，,1,y,2,，,0,z,16,内，立方体区域体积为,1,92,。,data=rand(10000,3);,x=4*data(:,1);,y=-1+3*data(:,2);,z=16*data(:,3);,II=find(x=y.2,M=length(II);,V=192*M/10000,第23页,例,5.15,用蒙特卡罗方法计算,其中,积分区域,是由 和,z,=1,所围成。,被积函数在求积区域上最大值为,2,。所以有四维超立方体,1,x,1,，,1,y,1,，,0,z,1,，,0,u,2,P=rand(10000,4);,x=-1+2*P(:,1);y=-1+2*P(:,2);,z=P(:,3);u=2*P(:,4);,II=find(zsqrt(x.2+y.2),M=length(II);,V=8*M/10000,第24页,试验,:,蒙特卡罗方法计算体积,&,x=2*rand-1,产生,1,到,1,之间随机数,y=2*rand-1,产生,1,到,1,之间随机数,z=2*rand;,产生,0,到,2,之间随机数,冰淇淋锥含于体积,=8,六面体,2,2,因为,rand,产生,0,到,1,之间随机数,所以,N,个点均匀分布于六面体中,锥体中占有,m,个,则锥体与六面体体积之比近似为,m:N,第25页,function q,error=MonteC(L),if nargin=0,L=7;end,N=10000;,for k=1:L,P=rand(N,3);,x=2P(:,1)-1;,y=2*P(:,2)-1;,z=2*P(:,3);,R2=x.2+y.2;R=sqrt(R2);,II=find(z=R,m=length(II);q(k)=8*m/N;,end,error=q-pi;,试验参考程序,蒙特卡罗方法计算体积,半球体积,圆锥体积,第26页,试验任务一,：统计,L,次试验试验数据及误差,试验任务二,：,修改试验程序,MonteC,计算,L,次试验数据均值及均值误差,(mean,计算平均值,),序号,1 2 3 4 5 6,7,数据,误差,L 8 16 32 64 128 256,均值,误差,第27页,function icecream(m,n)if nargin=0,m=20;n=100;end,t=linspace(0,2*pi,n);,r=linspace(0,1,m);,x=r*cos(t);y=r*sin(t);,z1=sqrt(x.2+y.2);,z2=1+sqrt(1+eps-x.2-y.2);,X=x;x;Y=y;y;,Z=z1;z2;,mesh(X,Y,Z),view(0,-18),colormap(0 0 1),axis off,冰淇淋锥体积,冰淇淋锥图形绘制程序,第28页,5.,下面程序绘出二维图形填充图,(,右图,),。分析每条语句功效,给程序中语句写注记,x1=-1:0.1:1;,y1=x1.2.(1/3);,x2=1:-0.1:-1;,y2=2-x2.2;,fill(x1,x2,y1,y2,g),axis off,第29页,谢谢,第30页,