资源描述
实验一随机数的产生及蒙特卡洛随机模拟方法实验一1实验目的实验目的实验内容实验内容学习学习随机数的产生及蒙特卡洛随机模拟方法随机数的产生及蒙特卡洛随机模拟方法的基本过程与方法。的基本过程与方法。实验作业实验作业2 2、蒙特卡洛随机模拟实例。、蒙特卡洛随机模拟实例。1 1、产生随机数的计算机命令。、产生随机数的计算机命令。实验目的实验内容学习随机数的产生及蒙特卡洛随机模拟方法实验作2数学模拟的方法数学模拟的方法 在实际问题中，面对一些带随机因素的复杂系在实际问题中，面对一些带随机因素的复杂系统，用分析方法建模常常需要作许多简化假设，统，用分析方法建模常常需要作许多简化假设，与面临的实际问题可能相差甚远，以致解答根本与面临的实际问题可能相差甚远，以致解答根本无法应用。无法应用。这时，计算机模拟几乎成为唯一的选这时，计算机模拟几乎成为唯一的选择。择。在一定的假设条件下，运用数学运算模拟系统在一定的假设条件下，运用数学运算模拟系统的运行，称为数学模拟。现代的数学模拟都是在计的运行，称为数学模拟。现代的数学模拟都是在计算机上进行的，称为计算机模拟。算机上进行的，称为计算机模拟。计算机模拟可以反复进行，改变系统的结构和计算机模拟可以反复进行，改变系统的结构和系数都比较容易。系数都比较容易。数学模拟的方法 在实际问题中，面对一些带随机因素的复杂3一、随机数的产生一、随机数的产生一、随机数的产生4一）产生模拟随机数的计算机命令一）产生模拟随机数的计算机命令 在在MatlabMatlab软件中，可以直接产生满足各种分布的随软件中，可以直接产生满足各种分布的随机数，命令如下：机数，命令如下：1产生产生m*n阶阶(a，b)均匀分布均匀分布U(a，b)的随机数矩阵：的随机数矩阵：unifrnd(a,b,m,n)产生一个产生一个a，b均匀分布的随机数：均匀分布的随机数：unifrnd(a,b)当只知道一个随机变量取值在（当只知道一个随机变量取值在（a，b）内，但）内，但不知道（也没理由假设）它在何处取值的概率大，不知道（也没理由假设）它在何处取值的概率大，在何处取值的概率小，就只好用在何处取值的概率小，就只好用U（a，b）来模拟它。）来模拟它。一）产生模拟随机数的计算机命令 在Matlab软件中，52产生产生mm*nnmm*nn阶离散均匀分布的随机数矩阵：阶离散均匀分布的随机数矩阵：R=unidrnd(N)R=unidrnd(N,mm,nn)2产生mm*nn阶离散均匀分布的随机数矩阵：6当研究对象视为大量相互独立的随机变量之和，当研究对象视为大量相互独立的随机变量之和，且其中每一种变量对总和的影响都很小时，可以且其中每一种变量对总和的影响都很小时，可以认为该对象服从正态分布。认为该对象服从正态分布。当研究对象视为大量相互独立的随机变量之和，且其中每一种变量对7若连续型随机变量若连续型随机变量X的概率密度函数为的概率密度函数为其中其中 00为常数，则称为常数，则称X服从参数为服从参数为 的指数分布。的指数分布。指数分布的期望值为指数分布的期望值为 若连续型随机变量X的概率密度函数为指数分布的期望值为 8排队服务系统中顾客到达间隔、质量与可靠性排队服务系统中顾客到达间隔、质量与可靠性中电子元件的寿命通常服从指数分布。中电子元件的寿命通常服从指数分布。例例 顾客到达某商店的间隔时间服从参数为顾客到达某商店的间隔时间服从参数为1010(分钟分钟)的指数分布的指数分布(指数分布的均值为指数分布的均值为10)-指两个顾客到达商店的平均间隔时间指两个顾客到达商店的平均间隔时间是是1010分钟分钟.即平均即平均1010分钟到达分钟到达1 1个顾客个顾客.顾客到顾客到达的间隔时间可用达的间隔时间可用exprnd(10)exprnd(10)模拟。模拟。排队服务系统中顾客到达间隔、质量与可靠性中电子元件的寿命通常9设离散型随机变量设离散型随机变量X X的所有可能取值为的所有可能取值为0,1,2,0,1,2,且且取各个值的概率为取各个值的概率为其中其中 00为常数，则称为常数，则称X服从参数为服从参数为 的的泊泊松分布松分布。泊泊松分布在排队系统、产品检验、天文、物理等松分布在排队系统、产品检验、天文、物理等领域有广泛应用。领域有广泛应用。泊泊松分布的期望值为松分布的期望值为设离散型随机变量X的所有可能取值为0,1,2,且取各个值106产生产生1个参数为个参数为n,p的二项分布的随机数的二项分布的随机数binornd(n,p)，产生产生m n个参数为个参数为n,p的二项分布的的二项分布的随机数随机数binornd(n,p,m,n)。掷掷一一枚枚均均匀匀硬硬币币，正正面面朝朝上上的的次次数数X X服服从参数为，从参数为，p p的二项分布的二项分布,XB(1,p),XB(1,p)6产生1个参数为n,p的二项分布的随机数binornd(n11总结：常见分布的随机数产生语句总结：常见分布的随机数产生语句总结：常见分布的随机数产生语句12v补充：随机数的产生命令补充：随机数的产生命令vMATLAB可以直接产生满足各种分布的随机数可以直接产生满足各种分布的随机数v具体具体命令如下命令如下:v 产生产生mn阶阶0,1上均匀分布的随机数矩阵上均匀分布的随机数矩阵v rand(m,n)v产生一个产生一个0,1上均匀分布的随机数上均匀分布的随机数v randv 产生产生mn阶阶a,b上均匀分布的随机数矩阵上均匀分布的随机数矩阵v unifrnd(a,b,m,n)v产生一个产生一个a,b上上均匀分布的随机数均匀分布的随机数v unifrnd(a,b)v 产生一个产生一个1:n的随机排列的随机排列(元素均出现且不重复元素均出现且不重复)v p=randperm(n)v注意注意:randperm(6)与与unifrnd(1,6,1,6)的区别的区别补充：随机数的产生命令13 产生产生mn阶阶均值为均值为mu方差为方差为sigma的正态的正态分布的随机分布的随机数矩阵数矩阵 normrnd(mu,sigma,m,n)产生一个均值为产生一个均值为mu方差为方差为sigma的正态分布的随机数的正态分布的随机数 normrnd(mu,sigma)产生产生mn阶期望值为阶期望值为mu(mu=1/)的指数分布的随机的指数分布的随机数矩阵数矩阵 exprnd(mu,m,n)产生一个期望值为产生一个期望值为mu的指数分布的随机数的指数分布的随机数 exprnd(mu)注意注意:产生一个参数为产生一个参数为的指数分布的随机数应输入的指数分布的随机数应输入 exprnd(1/)产生mn阶均值为mu方差为sigma的正态分布的随机数14v 产生产生mn阶阶参数为参数为A1,A2,A3的指定分布的指定分布name的的随机数矩随机数矩阵阵 random(name,A1,A2,A3,m,n)v产生一个参数为为产生一个参数为为A1,A2,A3的指定分布的指定分布name的随机数的随机数 random(name,A1,A2,A3)v举例举例:产生产生24阶的均值为阶的均值为0方差为方差为1的正态分布的随机数矩阵的正态分布的随机数矩阵 random(Normal,0,1,2,4)vname的取值可以是的取值可以是(详情参见详情参见help random)：vnorm or Normal/unif or Uniformvpoiss or Poisson/beta or Betavexp or Exponential/gam or Gammavgeo or Geometric/unid or Discrete Uniformv 产生mn阶参数为A1,A2,A3的指定分布name15二、蒙特卡罗随机模拟二、蒙特卡罗随机模拟二、蒙特卡罗随机模拟16 蒙特卡洛（蒙特卡洛（Monte CarloMonte Carlo）方法是一种应）方法是一种应用随机数来进行计算机模拟的方法此方法对用随机数来进行计算机模拟的方法此方法对研究的系统进行随机观察抽样，通过对样本值研究的系统进行随机观察抽样，通过对样本值的统计分析，求得所研究系统的某些参数的统计分析，求得所研究系统的某些参数 蒙特卡洛（Monte Carlo）方法是一种应用17用蒙特卡洛方法进行计算机模拟的步骤用蒙特卡洛方法进行计算机模拟的步骤：1 1 设计一个逻辑框图，即模拟模型这个框设计一个逻辑框图，即模拟模型这个框图要正确反映系统各部分运行时的逻辑关系。图要正确反映系统各部分运行时的逻辑关系。2 2 模拟随机现象可通过具有各种概率分布模拟随机现象可通过具有各种概率分布的模拟随机数来模拟随机现象的模拟随机数来模拟随机现象用蒙特卡洛方法进行计算机模拟的步骤：1 设计一个逻辑框图18一）频率的稳定性模拟一）频率的稳定性模拟 1.1.事件的频率事件的频率在一组不变的条件下，重复作在一组不变的条件下，重复作n n次试验，记次试验，记m m是是n n次试验中次试验中事件事件A A发生的次数。发生的次数。频率频率 f=m/n f=m/n 2.2.频率的稳定性频率的稳定性 掷一枚均匀硬币，记录掷硬币试验中频率掷一枚均匀硬币，记录掷硬币试验中频率P*P*的波动情况。的波动情况。一）频率的稳定性模拟 1.事件的频率2.频率的稳定性 掷一枚19function liti1(p,mm)pro=zeros(1,mm);randnum=binornd(1,p,1,mm)a=0;for i=1:mm a=a+randnum(1,i);pro(i)=a/i;end pro=pronum=1:mm;plot(num,pro)在在MatlabMatlab中编辑中编辑.m.m文件如下：文件如下：function liti1(p,mm)在Matlab中编辑20在在MatlabMatlab命令行中输入以下命令：命令行中输入以下命令：liti1(0.5,1000)在Matlab命令行中输入以下命令：liti1(0.5,1021在在MatlabMatlab命令行中输入以下命令：命令行中输入以下命令：liti1(0.5,10000)在Matlab命令行中输入以下命令：liti1(0.5,1022练练习习掷掷一一枚枚不不均均匀匀硬硬币币，正正面面出出现现概概率率为为0.30.3，记记录录前前10001000次次掷掷硬硬币币试试验验中中正正面面频率的波动情况，并画图。频率的波动情况，并画图。练习掷一枚不均匀硬币，正面出现概率为0.3，记录前1000次23在在MatlabMatlab命令行中输入以下命令：命令行中输入以下命令：liti1(0.3,1000)在Matlab命令行中输入以下命令：liti1(0.3,1024 二）几何概率模拟二）几何概率模拟1.1.定定义 向向任任一一可可度度量量区区域域G G内内投投一一点点，如如果果所所投投的的点点落落在在G G中中任任意意可可度度量量区区域域g g内内的的可可能能性性与与g g的的度度量量成成正正比比，而而与与g g的的位位置置和和形形状状无无关关，则则称称这这个个随随机机试试验验为为几几何何型型随随机试验。或简称为几何概型。机试验。或简称为几何概型。二）几何概率模拟1.定义 向任一可度量区域G内投一点，如果252.概率计算概率计算 P(A)=A的度量的度量/S的度量的度量例例5 两两人人约约定定于于12点点到到1点点到到某某地地会会面面，先先到到者者等等20分钟后离去，试求两人能会面的概率？分钟后离去，试求两人能会面的概率？解：设解：设x,y分别为甲、乙到达时刻分别为甲、乙到达时刻(分钟分钟)令令A=两人能会面两人能会面=(x,y)|x-y|20，x60,60,y y6060P(A)=A的面积的面积/S的面积的面积=（602-402）/602=5/9=0.55562.概率计算 P(A)=A的度量/S的度量例526function proguji=liti5(mm)%mm 是随机实验次数是随机实验次数frq=0;randnum1=unifrnd(0,60,mm,1);randnum2=unifrnd(0,60,mm,1);randnum=randnum1-randnum2;proguji=0;for ii=1:mm if abs(randnum(ii,1)=20 frq=frq+1;endendproguji=frq/mmliti5(10000)proguji=0.5557在在MatlabMatlab中编辑中编辑.m.m文件如下：文件如下：function proguji=liti5(mm)%mm 27三）三）蒲丰投针实验蒲丰投针实验：v法国科学家蒲丰法国科学家蒲丰(Buffon)(Buffon)在在17771777年提出年提出的蒲丰投针实验是早期几何概率一个非常的蒲丰投针实验是早期几何概率一个非常著名的例子。著名的例子。蒲丰投针实验蒲丰投针实验的重要性并非的重要性并非是为了求得比其它方法更精确的是为了求得比其它方法更精确的值，而值，而是它开创了使用随机数处理确定性数学问题的先河，是它开创了使用随机数处理确定性数学问题的先河，是用偶然性方法去解决确定性计算的前导，是用偶然性方法去解决确定性计算的前导，由由此此可以可以领略到从领略到从“概率土壤概率土壤”上开出的一朵上开出的一朵瑰丽的鲜花瑰丽的鲜花蒙特卡罗蒙特卡罗方法方法(MC)MC)v蒲丰投针实验可归结为下面的数学问题：平面上画蒲丰投针实验可归结为下面的数学问题：平面上画有距离为有距离为a的一些平行线，向平面上任意投一根长的一些平行线，向平面上任意投一根长为为l(la)的针，假设针落在任意位置的可能性相同，的针，假设针落在任意位置的可能性相同，试求针与平行线相交的概率试求针与平行线相交的概率P(从而求从而求)三）蒲丰投针实验：28v蒲丰投针实验蒲丰投针实验：v 如右图所示，以如右图所示，以Mv表示针落下后的中点，表示针落下后的中点，v以以x表示表示M到最近一条到最近一条v平行线的距离，以平行线的距离，以表示针与此线的交角：表示针与此线的交角：v针落地的所有可能结果满足：针落地的所有可能结果满足：v其样本空间视作矩形区域其样本空间视作矩形区域,面积是面积是:v针与平行线相交的条件：针与平行线相交的条件：v它是样本空间它是样本空间子集子集A，面积是：，面积是：v积分计算积分计算 syms l phi;int(l/2*sin(phi),phi,0,pi);%ans=lv因此，针与平行线相交的概率为：因此，针与平行线相交的概率为：v从而有：从而有：特别当特别当 时时 p为统计频率为统计频率蒲丰投针实验：29蒲丰投针实验蒲丰投针实验的计算机模拟：的计算机模拟：format long;a=1;l=0.6;%显示精度显示精度,线线宽宽和和针长针长figure;axis(0,pi,0,a/2);%初始化绘图板初始化绘图板set(gca,nextplot,add);%初始化绘图方式为叠加初始化绘图方式为叠加counter=0;n=2010;%初始化计数器和设定初始化计数器和设定投针次数投针次数x=unifrnd(0,a/2,1,n);phi=unifrnd(0,pi,1,n);%样本空间样本空间for i=1:n if x(i)l*sin(phi(i)/2%满足此条件表示针与线的相交满足此条件表示针与线的相交 plot(phi(i),x(i),r.);frame(i)=getframe;%描点并取帧描点并取帧 title(Current Point,num2str(i),Total,num2str(n);counter=counter+1;%统计针与线相交的次数统计针与线相交的次数 endendfren=counter/n;pihat=2*l/(a*fren)%用用频率频率近似计算近似计算%movie(frame,1)%播放帧动画播放帧动画1次次蒲丰投针实验的计算机模拟：30蒲丰投针实验计算圆周率蒲丰投针实验计算圆周率v蒙特卡罗投点法是蒙特卡罗投点法是蒲丰投针实验蒲丰投针实验的推广：的推广：v在一个边长为在一个边长为a的正方形内随机投点，该点落在此的正方形内随机投点，该点落在此正方形的内切圆中的概率应为该内切圆与正方形正方形的内切圆中的概率应为该内切圆与正方形的面积比值，即的面积比值，即vn=10000;a=2;m=0;vfor i=1:nv x=rand(1)*a;y=rand(1)*a;v if(x-a/2)2+(y-a/2)2=(a/2)2)v m=m+1;v endvendvdisp(投点法近似计算的投点法近似计算的为:,num2str(4*m/n);xyo(a/2,a/2)xyo蒲丰投针实验计算圆周率蒙特卡罗投点法是蒲丰投针实验的推广：31作作业业：1.1.掷掷两两枚枚不不均均匀匀硬硬币币，每每枚枚正正面面出出现现概概率率为为0.40.4，记记录录前前10001000次次掷掷硬硬币币试试验验中中两两枚枚都都为为正正面面频频率的波动情况，并画图。率的波动情况，并画图。2:两船欲停靠同一个码头两船欲停靠同一个码头,设两船到达码头的时间设两船到达码头的时间各不相干，而且到达码头的时间在一昼夜内是等可各不相干，而且到达码头的时间在一昼夜内是等可能的能的.如果两船到达码头后需在码头停留的时间分如果两船到达码头后需在码头停留的时间分别是别是1 小时与小时与2 小小 时，试求在一昼夜内，任一船到时，试求在一昼夜内，任一船到达时，需要等待空出码头的概率达时，需要等待空出码头的概率.（频率估计概率）（频率估计概率）作业：1.掷两枚不均匀硬币，每枚正面出现概率为0.4，记录前32
展开阅读全文