实验十五 MATLAB的蒙特卡洛仿真.docx
实验十五 MATLAB的蒙特卡洛仿真实验十五: MATLAB的蒙特卡洛仿真 一、实验目的 1. 了解蒙特卡洛仿真的基本概念。 2. 了解蒙特卡洛仿真的某些应用 二实验内容与步骤 1 蒙特卡洛仿真的简介 随机模拟方法,也称为Monte Carlo方法,是一种基于“随机数”的计算方法。这一方法源于美国在第一次世界大战进行的研制原子弹的“曼哈顿计划”。该计划的主持人之一、数学家冯·诺伊曼用驰名世界的赌城摩纳哥的Monte Carlo来命名这种方法,为它蒙上了一层神秘色彩。冯·诺伊曼是公理化方法和计算机体系的领袖人物,Monte Carlo方法也是他的功劳。 事实上,Monte Carlo方法的基本思想很早以前就被人们所发现和利用。早在17世纪,人们就知道用事件发生的“频率”来决定事件的“概率”。18世纪下半叶的法国学者Buffon提出用投点试验的方法来确定圆周率的值。这个著名的Buffon试验是Monte Carlo方法的最早的尝试! 历史上曾有几位学者相继做过这样的试验。不过他们的试验是费时费力的,同时精度不够高,实施起来也很困难。然而,随着计算机技术的飞速发展,人们不需要具体实施这些试验,而只要在计算机上进行大量的、快速的模拟试验就可以了。Monte Carlo方法是现代计算技术的最为杰出的成果之一,它在工程领域的作用是不可比拟的。 蒙特卡洛(Monte Carlo)模拟是一种通过设定随机过程,反复生成时间序列,计算参数估计量和统计量,进而研究其分布特征的方法。具体的,当系统中各个单元的可靠性特征量已知,但系统的可靠性过于复杂,难以建立可靠性预计的精确数学模型或模型太复杂而不便应用时,可用随机模拟法近似计算出系统可靠性的预计值;随着模拟次数的增多,其预计精度也逐渐增高。由于涉及到时间序列的反复生成,蒙特卡洛模拟法是以高容量和高速度的计算机为前提条件的,因此只是在近些年才得到广泛推广。 2 MC 的原理 针对实际问题建立一个简单且便于实现的概率统计模型,使问题的解对应于该模型中随机变量的概率分布或其某些数字特征,比如,均值和方差等。所构造的模型在主要特征参量方面要与实际问题或系统相一致的。 根据模型中各个随机变量的分布,在计算机上产生随机数,实现一次模拟过程所需的足够数量的随机数。通常先产生均匀分布的随机数,然后生成服从某一分布的随机数,再进行随机模拟试验。 收敛性: 由大数定律, Monte-Carlo模拟的收敛是以概率而言的 误差: 用频率估计概率时误差的估计,可由中心极限定理,给定置信水平 的条件下,有: s=Var(g(X) |e|£U1-a/2sN Us模拟次数:由误差公式得 N³(1-a/2)2e 3 定积分的MC计算原理 事实上,不少的统计问题,如计算概率、各阶距等,最后都归结为定积分的近似计算问题。设 a,b,有限,0 ,£f(x)£M W=(x,y):a£x£b,0£y£M并设(X,Y)是在W 上均匀分布的二维随机变量,其联合密度函数为b1I(a£x£b,0£y£M) 。则易见 M(b-a)q=òf(x)dx 是 W 中 y=f(x) 曲线下方的面积。假设我们向W 中进行随机投a点,若点落在 y=f(x)下方,(即 y<f(x)称为中的,否则称为不中,则点中的概率为 p=qM(b-a) 。若我们进行了 n次投点,其中n0次中的,则用频率来估计概率p。即p=n0 。 n4 蒙特卡洛仿真应用举例 例1 计算定积分 q=ò(cosx+2.0)dx04事实上,其精确解为 q=8.0+sin4.0 = 7.2432用随机投点法求解如下: sjtdf(0,4,4,1000000) result = 7.2336 function result=sjtdf(a,b,m,mm) %a是积分的下限 %b是积分的上限 %m是函数的上界 %mm 是随机实验次数 frq=0; xrandnum = unifrnd(a,b,1,mm); yrandnum = unifrnd(0,m,1,mm); for ii=1:mm if (cos(xrandnum(1,ii)+2>=yrandnum(1,ii) frq=frq+1; end end result=frq*m*(b-a)/mm 例2 p的计算 10.80.60.40.20.20.40.60.8>> sjtdf_pi1(1000) piguji = 3.0520 >> sjtdf_pi1(10000) piguji = 3.1204 >> sjtdf_pi1(100000) piguji = 3.1296 function piguji=sjtdf_pi1(mm) %mm 是随机实验次数 frq=0;xrandnum = unifrnd(0,1,1,mm); yrandnum = unifrnd(0,1,1,mm); for ii=1:mm if xrandnum(1,ii)2+yrandnum(1,ii)2<=1 frq=frq+1; end end piguji=4*frq/mm 思考题 运用定积分的MC计算原理,用随机投点法计算二重定积分 q=òòcosx+siny+3.0dxdy0113