实验十五：MATLAB的蒙特卡洛仿真

实验十五： MATLAB 的蒙特卡洛仿真
一、实验目的
1. 了解蒙特卡洛仿真的基本概念。

2. 了解蒙特卡洛仿真的某些应用
二．实验内容与步骤
1． 蒙特卡洛（Monte Carlo ）仿真的简介
随机模拟方法，也称为Monte Carlo 方法，是一种基于“随机数”的计算方法。

这一方法源于美国在第一次世界大战进行的研制原子弹的“曼哈顿计划”。

该计划的主持人之一、数学家冯·诺伊曼用驰名世界的赌城—摩纳哥的Monte Carlo 来命名这种方法，为它蒙上了一层神秘色彩。

冯·诺伊曼是公理化方法和计算机体系的领袖人物，Monte Carlo 方法也是他的功劳。

事实上，Monte Carlo 方法的基本思想很早以前就被人们所发现和利用。

早在17世纪，人们就知道用事件发生的“频率”来决定事件的“概率”。

18世纪下半叶的法国学者Buffon 提出用投点试验的方法来确定圆周率π的值。

这个著名的Buffon 试验是Monte Carlo 方法的最早的尝试！
历史上曾有几位学者相继做过这样的试验。

不过他们的试验是费时费力的，同时精度不够高，实施起来也很困难。

然而，随着计算机技术的飞速发展，人们不需要具体实施这些试验，而只要在计算机上进行大量的、快速的模拟试验就可以了。

Monte Carlo 方法是现代计算技术的最为杰出的成果之一，它在工程领域的作用是不可比拟的。

蒙特卡洛(Monte Carlo)模拟是一种通过设定随机过程，反复生成时间序列，计算参数估计量和统计量，进而研究其分布特征的方法。

具体的，当系统中各个单元的可靠性特征量已知，但系统的可靠性过于复杂，难以建立可靠性预计的精确数学模型或模型太复杂而不便应用时，可用随机模拟法近似计算出系统可靠性的预计值；随着模拟次数的增多，其预计精度也逐渐增高。

由于涉及到时间序列的反复生成，蒙特卡洛模拟法是以高容量和高速度的计算机为前提条件的，因此只是在近些年才得到广泛推广。

2． MC 的原理
针对实际问题建立一个简单且便于实现的概率统计模型，使问题的解对应于该模型中随机变量的概率分布或其某些数字特征，比如，均值和方差等。

所构造的模型在主要特征参量方面要与实际问题或系统相一致的。

根据模型中各个随机变量的分布，在计算机上产生随机数，实现一次模拟过程所需的足够数量的随机数。

通常先产生均匀分布的随机数，然后生成服从某一分布的随机数，再进行随机模拟试验。

收敛性: 由大数定律, Monte-Carlo 模拟的收敛是以概率而言的．
误差: 用频率估计概率时误差的估计，可由中心极限定理，给定置信水平 的条件下，有：
•
•
模拟次数:由误差公式得N U σεα2/1||-≤))((X g Var =σ2
2/1)(ε
σα-≥U N
3． 定积分的MC 计算原理
事实上，不少的统计问题，如计算概率、各阶距等，最后都归结为定积分的近似计算问
题。

设 a,b ，有限， , (){}M y b x a y x ≤≤≤≤=Ω0,:,并设()Y X ,是在Ω 上均匀分布的二维随机变量，其联合密度函数为(
)()M y b x a I a b M ≤≤≤≤-0,1 。

则易见 ()dx x f b a
⎰=θ 是 Ω 中 ()x f y = 曲线下方的面积。

假设我们向Ω 中进行随机投点，若点落在 ()x f y =下方，(即 ()x f y <称为中的，否则称为不中，则点中的概率为 ()a b M p -=θ。

若我们进行了 n 次投点，其中0n 次中的，则用频率来估计概率p 。

即n
n p 0= 。

4． 蒙特卡洛仿真应用举例
例1 计算定积分
事实上，其精确解为
用随机投点法求解如下：
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);
()dx x ⎰+=4
00.2cos θ7.2432
= 0.4sin 0.8+=θ()M
x f ≤≤0
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 的计算
（单位圆的面积等于 p ）
>> 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 计算原理，用随机投点法计算二重定积分
（注：可编辑下载，若有不当之处，请指正，谢谢!） ⎰⎰++=103
10.3sin cos dxdy y x θ。