蒙特卡罗方法及应用实验讲义2016资料

蒙特卡罗方法及应用

实验讲义

东华理工大学核工系

2016.8

实验一 蒙特卡罗方法基本思想

一、实验目的

1、了解蒙特卡罗方法方法的基本思想；

2、掌握蒙特卡罗方法计算面积、体积的方法；

3、掌握由已知分布的随机抽样方法。

二、实验原理

Monte Carlo 方法，又称统计模拟方法或计算机随机模拟方法，是一种基于“随机数”进行数值模拟的方法，一种采用统计抽样理论近似求解物理或数学问题的方法。

如待求量可以表述成某些特征量的期望值、某些事件出现的概率或两者的函数形式，那么可采用蒙特卡罗方法求解。在求解某些特征量的期望值或某些事件出现的概率时，必须构建合符实际的数学模型。例如采用蒙特卡罗方法计算某函数所围面积时，构建的数学模型是构造一已知面积的可均匀抽样区域，在该区域投点，由伯努利定理大数定理可知，进入待求区域投点的频率依概率1收敛于该事件出现的概率（面积之比）。

由已知分布的随机抽样方法指的是由已知分布的总体中抽取简单子样。具体方法很多，详见教材第三章。 三、实验内容

1、安装所需计算工具（MATLAB 、fortran 、C++等）；

2、学习使用rand(m,n)、unifrnd(a,b,m,n)函数

3、求解下列问题：

3.0、蒲丰氏投针求圆周率。

3.1、给定曲线y =2 – x 2 和曲线y 3 = x 2，曲线的交点为：P 1( – 1，1 )、P 2( 1，1 )。曲线围成平面有限区域，用蒙特卡罗方法计算区域面积；

3.2

、计算1z z ⎧≥⎪⎨≤⎪⎩所围体积

其中{(,,)|11,11,02}x y z x y z Ω=-≤≤-≤≤≤≤。

4、对以下已知分布进行随机抽样：

4.1、()()[]2

3

321,0,12

f x x x x =+

-∈； 4.2、()()

()[]11

,1,21E f x f x x E k E =

⋅∈+ 其中()()()()()2

123

221111211411ln 212221E x f x E x x x x E k E E E E E ⎧+-⎛⎫=+-+⎪ ⎪⋅⎝⎭⎪

⎨+⎡⎤⎪=-⋅+++-

⎢⎥⎪+⎣⎦⎩

。

四、实验报告编写

1、给出各题的抽样程序并解释语句的含义；

2、给出3.1和3.2抽样结果误差随抽样次数的关系图，并解释原因；

表1 实验记录表

3、给出4.1和4.2的抽样框图、试验累积频率与理论累积频率关系图，并给出抽样次数（>106）与抽样时间。

实验二 由已知分布的随机抽样方法

一、实验目的

1、掌握由已知分布的随机抽样方法。

2、用编程语言实现某具体随机抽样方法。 二、实验原理

由已知分布的随机抽样方法指的是由已知分布的总体中抽取简单子样。具体方法很多，本实验综合直接抽样方法、挑选抽样方法和替换抽样方法，以散射方位角余弦分布的抽样为例。实验原理详见教材对应章节。 1.连续型分布的直接抽样方法

对于连续型分布，如果分布函数F(x) 的反函数F －1(x)存在，则直接抽样方法是：

1()F X F ξ-=

2.挑选抽样方法

为了实现从己知分布密度函数f(x)抽样，选取与f(x)取值范围相同的分布密度函数h(x)，如果

()

sup

()

x f x M h x -∞<<∞

=<∞ 则挑选抽样方法为：

3.替换法抽样方法

为了实现某个复杂的随机变量 y 的抽样，将其表示成若干个简单的随机变量x 1，x 2，…，x n 的函数

12(,,

,)n y g x x x =

得到 x 1，x 2，…，x n 的抽样后，即可确定 y 的抽样，这种方法叫作替换法

抽样。

蒲丰氏问题的算法 如何产生任意的（x,θ）？

x 在［0，a ］上任意取值，表示x 在［0，a ］上是均匀分布的，其分布密度函数为：

11/,0()0,

a x a

f x ≤≤⎧=⎨

⎩其他 类似地，θ的分布密度函数为：

21/,0()0,f πθπ

θ≤≤⎧=⎨

⎩

其他 因此，产生任意的（x,θ）的过程就变成了由f 1(x)抽样x 及由f 2(θ)抽样θ的过程了。由此得到：

1

2

x a ξθπξ==

其中ξ1，ξ2均为（0,1）上均匀分布的随机变量。

每次投针试验，实际上变成在计算机上从两个均匀分布的随机变量中抽样得到（x,θ），然后定义描述针与平行线相交状况的随机变量s(x,θ)，为

1,sin (,)0x l s x θ

θ≤⋅⎧=⎨

⎩

,其他 如果投针Ｎ次，则

1

1(,)N

N i i i s s x N θ==∑

是针与平行线相交概率Ｐ的估计值。事实上，

12sin 0

0(,)()()d d d d 2ππl P s x f x f x x l

a a

πθθθθθ===

⎰⎰⎰

⎰

于是有

22πN

l l

aP as =

≈

1、给出源程序程序并解释语句的含义；

2、作出抽样框图、试验累积频率与理论累积频率关系图，并给出抽样次数（>106）与抽样时间。