《数学建模蒙特卡罗方法学习教案.pptx》由会员分享，可在线阅读，更多相关《数学建模蒙特卡罗方法学习教案.pptx（30页珍藏版）》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、数学数学(shxu)建模蒙特卡罗方法建模蒙特卡罗方法第一页，共30页。什么(shn me)叫蒙特卡罗方法？蒙特卡罗方法又称统计模拟法、随机抽样技术，是一种随机模拟方法，以概率和统计理论方法为基础的一种计算方法，是使用随机数（或伪随机数）来解决很多计算问题的方法。将所求解的问题同一定的概率模型相联系(linx)，用电子计算机实现统计模拟或抽样，以获得问题的近似解。为象征性地表明这一方法的概率统计特征，故借用赌城蒙特卡罗命名。第1页/共30页第二页，共30页。基本(jbn)思想当所求问题的解是某个事件的概率，或者(huzh)是某个随机变量的数学期望，或者(huzh)是与概率，数学期望有关的量时，通
2、过某种试验的方法，得出该事件发生的概率，或者(huzh)该随机变量若干个具体观察值的算术平均值，通过它得到问题的解。当随机变量的取值仅为1或0时，它的数学期望就是某个事件的概率。或者(huzh)说，某种事件的概率也是随机变量（仅取值为1或0）的数学期望。第2页/共30页第三页，共30页。蒙特卡罗方法(fngf)的特点优点：优点：1、能够比较逼真地描述具有随机性质的事物的、能够比较逼真地描述具有随机性质的事物的 特点及物理实验过程特点及物理实验过程2、受几何条件、受几何条件(tiojin)限制小限制小3、收敛速度与问题的维数无关、收敛速度与问题的维数无关4、具有同时计算多个方案与多个未知量的能力
3、、具有同时计算多个方案与多个未知量的能力5、误差容易确定、误差容易确定6、程序结构简单，易于实现、程序结构简单，易于实现第3页/共30页第四页，共30页。缺点：缺点：缺点：缺点：1 1收敛速度慢收敛速度慢收敛速度慢收敛速度慢2 2误差具有概率性误差具有概率性误差具有概率性误差具有概率性3 3在粒子输运问题中，计算结果与系统大在粒子输运问题中，计算结果与系统大在粒子输运问题中，计算结果与系统大在粒子输运问题中，计算结果与系统大小有关小有关小有关小有关 所以在使用蒙特卡罗方法时，要所以在使用蒙特卡罗方法时，要所以在使用蒙特卡罗方法时，要所以在使用蒙特卡罗方法时，要“扬长扬长扬长扬长避短避短避短避短
4、”，只对问题中难以用解析（或数值），只对问题中难以用解析（或数值），只对问题中难以用解析（或数值），只对问题中难以用解析（或数值）方法处理的部分，使用蒙特卡罗方法计算，方法处理的部分，使用蒙特卡罗方法计算，方法处理的部分，使用蒙特卡罗方法计算，方法处理的部分，使用蒙特卡罗方法计算，对那些对那些对那些对那些(nxi)(nxi)能用解析（或数值）方法能用解析（或数值）方法能用解析（或数值）方法能用解析（或数值）方法处理的部分，应当尽量使用解析方法处理的部分，应当尽量使用解析方法处理的部分，应当尽量使用解析方法处理的部分，应当尽量使用解析方法第4页/共30页第五页，共30页。主要应用(yngyng)
5、范围粒子输运问题（实验物理，反应堆物理，高能物理）统计(tngj)物理典型数学问题真空技术激光技术以及医学生物探矿第5页/共30页第六页，共30页。什么(shn me)是随机数？在连续型随机变量的分布中，最简单而且最基本的分在连续型随机变量的分布中，最简单而且最基本的分布是单位均匀分布。由该分布抽取的简单子样称为随布是单位均匀分布。由该分布抽取的简单子样称为随机数序列，其中每一个体称为随机数机数序列，其中每一个体称为随机数符号符号(fho)(fho)：两个特点：独立性，均匀性两个特点：独立性，均匀性第6页/共30页第七页，共30页。产生(chnshng)随机数随机数表方法随机数表方法(fngf
6、)(fngf)物理方法物理方法(fngf)(fngf)第7页/共30页第八页，共30页。随机数表随机数表是由0,1,2,3,4,5,6,7,8,9十个数字组成，每个数字以0.1的等概率出现，数字之间相互独立，这些数字序列叫作随机数字序列。（如果(rgu)要得到n位有效数字的随机数，只需将表中每n个相邻的随机数字合并在一起，且在最高位的前边加上小数点即可。例如，某随机数表的第一行数字为7 6 3 4 2 5 8 9 1.，要想得到三位有效数字的随机数一次为0.763，0.425，0.891.）第8页/共30页第九页，共30页。物理(wl)方法利用某些物理现象，在计算机上增加些特殊设备，可以在计算
7、机利用某些物理现象，在计算机上增加些特殊设备，可以在计算机上直接产生随机数。上直接产生随机数。作为随机数发生器的物理源主要有两种：一种是根据放射性物质作为随机数发生器的物理源主要有两种：一种是根据放射性物质的放射性，另一种是利用计算机的固有的放射性，另一种是利用计算机的固有(gyu)(gyu)噪声。噪声。一般情况下，任意一个随机数在计算机内总是用二进制的数表示一般情况下，任意一个随机数在计算机内总是用二进制的数表示的：的：其中其中 或者为或者为0 0，或者为，或者为1 1。因此，利用物理方法。因此，利用物理方法在计算机产生随机数，就是要产生只取在计算机产生随机数，就是要产生只取0 0或或1 1
8、的随机数字序列，数的随机数字序列，数字之间相互独立，每个数字取字之间相互独立，每个数字取0 0或或1 1的概率均为的概率均为0.50.5第9页/共30页第十页，共30页。缺点(qudin)随机数表需在计算机中占有很大内存，而且也难以满足随机数表需在计算机中占有很大内存，而且也难以满足蒙特卡罗方法对随机数需要量非常大的要求，因此，该蒙特卡罗方法对随机数需要量非常大的要求，因此，该方法不适于在计算机上使用。方法不适于在计算机上使用。物理方法产生的随机数序列无法重复实现，不能进行程物理方法产生的随机数序列无法重复实现，不能进行程序复算。给验证结果带来很大困难。而且增加随机数发序复算。给验证结果带来很
9、大困难。而且增加随机数发生器和电路联接生器和电路联接(lin ji)(lin ji)等附加设备，费用昂贵。因此等附加设备，费用昂贵。因此该方法也不适合在计算机上使用。该方法也不适合在计算机上使用。第10页/共30页第十一页，共30页。伪随机数用递推公式用递推公式(gngsh)(gngsh)产生随机数序列。产生随机数序列。第11页/共30页第十二页，共30页。伪随机数存在的两个(lin)问题递推公式和初始值确定后，整个递推公式和初始值确定后，整个(zhngg)(zhngg)随机数序列便被唯一确定。不随机数序列便被唯一确定。不满足随机数相互独立的要求。满足随机数相互独立的要求。由于随机数序列是由递
10、推公式确定的，而在计算机上所能表示的由于随机数序列是由递推公式确定的，而在计算机上所能表示的00，11上的数又是有限的，因此，这种方法产生的随机数序列就不可能不出现上的数又是有限的，因此，这种方法产生的随机数序列就不可能不出现重复。随机数序列出现周期性的循环现象。重复。随机数序列出现周期性的循环现象。第12页/共30页第十三页，共30页。解决方案第一个问题：不能从本质上加以改变，但只要递推公式选的比较好，随第一个问题：不能从本质上加以改变，但只要递推公式选的比较好，随机数间的相互独立性是可以近似满足的。机数间的相互独立性是可以近似满足的。第二个问题：因为第二个问题：因为(yn wi)(yn w
11、i)用蒙特卡罗方法解任何具体问题时，所使用用蒙特卡罗方法解任何具体问题时，所使用的随机数的个数总是有限的，只要所用随机数的个数不超过伪随机数序的随机数的个数总是有限的，只要所用随机数的个数不超过伪随机数序列出现循环现象时的长度就可以了。列出现循环现象时的长度就可以了。第13页/共30页第十四页，共30页。应用(yngyng)：蒙特卡罗方法计算积分 可以通俗地说，蒙特卡罗方法是用随机试验的方法计算积分，即可以通俗地说，蒙特卡罗方法是用随机试验的方法计算积分，即将所要计算的积分看作服从某种分布密度函数将所要计算的积分看作服从某种分布密度函数f(r)f(r)的随机变量的随机变量(su j(su j
12、bin lin)bin lin)(r)(r)的数学期望的数学期望 通过某种试验，得到观察值通过某种试验，得到观察值r1r1，r2r2，rNrN（用概率语言来说，（用概率语言来说，从分布密度函数从分布密度函数f(r)f(r)中抽取个子样中抽取个子样r1r1，r2r2，rNrN，），将相应的，），将相应的个随机变量个随机变量(su j bin lin)(su j bin lin)的值的值g(r1)g(r1)，g(r2)g(r2)，g(rN)g(rN)的算术平均的算术平均值值 作为积分的估计值（近似值）。作为积分的估计值（近似值）。第14页/共30页第十五页，共30页。求积分求积分(jfn)(jfn
13、)（2.1）第15页/共30页第十六页，共30页。蒙特卡罗方法步骤(bzhu)如下：1 1、在区间【、在区间【a,ba,b】上利用计算机均匀产生】上利用计算机均匀产生n n个随机数个随机数x1x1，x2xnx2xn，这个可以，这个可以(ky)(ky)在在MATLABMATLAB软件中用软件中用unifrndunifrnd命令实现。命令实现。2 2、计算每一个随机数相应的被积函数值、计算每一个随机数相应的被积函数值f f（x1x1），），f f（x2x2）ff（xnxn）。）。3 3、计算被积函数值的平均值、计算被积函数值的平均值4 4、所以、所以2.12.1式的值式的值 第16页/共30页第十
14、七页，共30页。简单定积分简单定积分(jfn)(jfn)例子：用蒙特卡罗方法求例子：用蒙特卡罗方法求首先我们进行严格的数学计算，首先我们进行严格的数学计算，便于后面与蒙特卡洛计算方便于后面与蒙特卡洛计算方法所得结果形成对比：已知法所得结果形成对比：已知 的原函数是的原函数是 ，那么定积分，那么定积分值就是：值就是：我们可以在我们可以在MatlabMatlab中输入以下代码中输入以下代码(di m)(di m)进行精确计算：进行精确计算：expexp（2 2）-exp-exp（0 0），这个值是此定积分的真实值），这个值是此定积分的真实值.下面进行蒙特卡洛计算上述定积分，其下面进行蒙特卡洛计算上
15、述定积分，其MATLABMATLAB代码代码(di(di m)m)如下：如下：N=500;N=500;x=unifrnd(0,2,N,1);x=unifrnd(0,2,N,1);mean(2*exp(3*x.2)mean(2*exp(3*x.2)第17页/共30页第十八页，共30页。解释(jish)语句：上述三条语句完整实现了蒙特卡洛计算上述定积上述三条语句完整实现了蒙特卡洛计算上述定积分步骤分步骤.第一条语句是设定了停止条件第一条语句是设定了停止条件(tiojin)(tiojin)，共做共做N N次次Monte Carlo Monte Carlo 模拟模拟.第二条第二条语句实现了在积分区间上
16、均匀产生语句实现了在积分区间上均匀产生N N个随机数个随机数.第第三条语句实现蒙特卡洛计算方法的面积逼近三条语句实现蒙特卡洛计算方法的面积逼近.对对N N设置不同的值，观察所得蒙特卡洛计算方法定积设置不同的值，观察所得蒙特卡洛计算方法定积分值，如表分值，如表1 1所示，我们可以发现：当不断增大所示，我们可以发现：当不断增大N N值时，所得结果越来越接近真实值值时，所得结果越来越接近真实值.第18页/共30页第十九页，共30页。蒙特卡罗计算方法性质(xngzh)：蒙特卡洛计算方法依据概率统计理论，具有统计特性，蒙特卡洛计算方法依据概率统计理论，具有统计特性，主要表现在以下三个方面主要表现在以下三
17、个方面.1.1.蒙特卡洛计算方法具有随机性、不确定性蒙特卡洛计算方法具有随机性、不确定性.即每次运行即每次运行结果都会不一样，因为计算机产生随机数并不是可以重结果都会不一样，因为计算机产生随机数并不是可以重现的现的.2.2.蒙特卡洛计算方法具有统计稳定性。虽然每次运行产蒙特卡洛计算方法具有统计稳定性。虽然每次运行产生随机数是不一样的，但是随机数的概率分布是一样的，生随机数是不一样的，但是随机数的概率分布是一样的，所以蒙特卡洛计算方法可以满足概率统计的稳定性所以蒙特卡洛计算方法可以满足概率统计的稳定性.3.3.随着随机数数量增加，蒙特卡洛计算方法所得结果会随着随机数数量增加，蒙特卡洛计算方法所得
18、结果会更加逼近真实更加逼近真实(zhnsh)(zhnsh)值，这就是我们所讲的依概率收值，这就是我们所讲的依概率收敛到真实敛到真实(zhnsh)(zhnsh)值意义值意义.第19页/共30页第二十页，共30页。蒙特卡罗方法蒙特卡罗方法(fngf)(fngf)实验实验 面积、体积计算面积、体积计算(j sun)问问题题冰淇淋锥的体积计算冰淇淋锥的体积计算(j sun)思考题与练习题思考题与练习题第20页/共30页第二十一页，共30页。蒙特卡罗方法蒙特卡罗方法(fngf)随机投点试验求近似解随机投点试验求近似解引例引例.给定曲线给定曲线y=2 x2 和曲线和曲线y3=x2，曲线的交点，曲线的交点(
19、jiodin)为：为：P1(1，1)、P2(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/10000plot(x(II),y(II),g.)S=2.1136第21页/共30页第二十二页，共30页。例例5.14 5.14 计算计算(j sun)(j sun)其中其中D D为为y=x 2y=x 2与与y2=x y2=x 所围所围 D的的边边界界曲曲线线交交点点为为：(1，1)，(4，2)
20、，被被积积函函数数在在求求积积区区域域(qy)内内的的最最大大值值为为16。积积分分值值是是三三维维体体积积，该该三维图形位于立方体区域三维图形位于立方体区域(qy)0 x 4，1 y 2，0 z 16内，立方体区域内，立方体区域(qy)的体积为的体积为192。data=rand(10000,3);x=4*data(:,1);y=-1+3*data(:,2);z=16*data(:,3);II=find(x=y.2&x=y+2&zsqrt(x.2+y.2)&z=1&u=R&z=1+sqrt(1-R2);m=length(II);q(k)=8*m/N;enderror=q-pi;实验参考程序实验
21、参考程序(chngx)蒙特卡罗方法计算体积蒙特卡罗方法计算体积半球体积半球体积圆锥体积圆锥体积第25页/共30页第二十六页，共30页。实验实验(shyn)任务一：记录任务一：记录L次实验次实验(shyn)的实验的实验(shyn)数据及误差数据及误差实验任务二：修改实验程序实验任务二：修改实验程序MonteC计算计算L次实验数据均值次实验数据均值(jn zh)及均值及均值(jn zh)误差误差(mean 计算平均值计算平均值(jn zh)序号序号 1 2 3 4 5 6 7 数据数据误差误差 L 8 16 32 64 128 256均值均值误差误差第26页/共30页第二十七页，共30页。func
22、tion icecream(m,n)if nargin=0,m=20;n=100;endt=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冰淇淋锥体积冰淇淋锥体积冰淇淋锥图形冰淇淋锥图形(txng)绘制绘制程序程序第27页/共30页第二十八页，共30页。5.下面程序绘出二维图形填下面程序绘出二维图形填充图充图(右图右图)。分析。分析(fnx)每条每条语句功能语句功能,给程序中语句写注给程序中语句写注记记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第28页/共30页第二十九页，共30页。第29页/共30页第三十页，共30页。