蒙特卡罗 (Monte－ Carlo)模拟，又称蒙特卡
罗方法、统计试验法等,
M－ C模拟是静态模拟，描述特定时间点上
的系统行为,
模拟过程中
不出现时间
参数。
基本思想,把随机事件
（变量）的概率特征与
数学分析的解联系起来,
概率特征,随机事件的概率和随机变量的
数学期望等, 用试验方法确定
一, 蒙特卡罗法计算定积分
例 7.3.1 用 M－ C 模拟求圆周率 π的估计值,
1
10
设二维随机变量
(X,Y)在正方形内
服从均匀分布,
(X,Y)落在圆内的概率为,
4}1{
22 ???? YXP
计算机上做 n次掷点试验：
产生 n 对二维随机点 (xi,yi), i＝ 1,2,…,n,
xi 和 yi是 RND 随机数对,
检查每对随机数是否满足,
122 ?? ii yx
相当于第 i
个随机点落
在 1/4圆内,
若有 k个点落在 l／ 4圆内
随机事件“点落入 1/4圆内”的
频率为 k/n
根据概率论中的大数定律,事件发生的频率
依概率收敛于事件发生的概率 p,即有
1}{lim ????
??
pP nk
n
得圆周率 π的估计值为
n
k4? ??
且当 试验次数足够大 时,其精度也随之提高,
分析,实际上概率值为
41
1
0
2 ??? ? dxx
恰为 1/4圆
的面积
频率法,利用随机变量落进指定区域内的频
率来计算定积分,
平均值法,利用随机变量的平均值 (数学期望 )
来计算定积分,
?? ba dxxfI )(
平均值法 的算法如下：
(1)产生 RND 随机数,r1,r2,…, rn；
（ 2） 令 ui=a＋ (b－ a)ri,i=1,2,…, n；
（ 3） 计算 作为 I 的估计值,??
?
n
i
iufn
ab
1
)(
原理分析：
设随机变量 ζ1,ζ2,…, ζn相互独立,
且 ζi～ U(0,1)
{f(ξi)},i=1,2,…, n 相互独立同分布
ab
Idxxf
abfE
b
ai ????? ? )(
1)]([
由 (强 )大数定律知
..)(1lim
1
ea
ab
I
f
n
n
i
i
n
?
??? ?
??
以概率为
1 成立
当 n 足够大时，得近似公式：
? ?????
?
n
i
in
b
a fabdxxfI
1
1 )()()(
注,平均值法本质上是用样本平均值作为总体教学期望的估计。
二, 蒙特卡罗模拟试验次数的确定
M－ C 模拟是一种试验近似方法,试验次数
如何确定？
？希望：模拟次数较少,模拟精度较高
频率法的讨论
用事件 A出现的频率作为概率 p的估计,
n
kp n??
问题,试验次数 n 多大时，对给定的置信度
1－ α（ 0<α<1）,估计精度达到 ε.
即问：取多大的 n使
? ? ???
?
?
?
?
?
? ??????? 1? p
n
kPppP n
成立？
证明
频率法是事件 A出现的频率作为概率 p的估计
n
kp n??
答案,2
2
)1(
??
?? zppn
其中,zα是正态分布的临界值,
n次独立试验中 A出现的次数 kn～ B(n,p).由中
心极限定理知
? ? ? ?)()(? ?????????? pnkpnPppP n
??
?
?
?
??
?
?
?
?
??
?
??
?
???
)1()1()1( pnp
n
pnp
npk
pnp
nP n
)
)1(
()
)1(
(
pnp
n
pnp
n
?
????
?
???
1)
)1(
(2 ?
?
???
pnp
n
????
?
?? 11)
)1(
(2
pnp
n令
平均值法 在给定 α和 ε下所需的试验次数
的估计式为
?
?
?
?
?
? ?
?
?
??
?
?
0
1
2
0
2
222
)(
1
1
/
n
i
i xx
n
S
Szn
查得正态分布的临界值 zα,可解得
2
2
)1(
??
?? zppn
试验次数估计式的分析
?
?
?
?
?
? ?
?
?
??
?
?
0
1
2
0
2
222
)(
1
1
/
n
i
i xx
n
S
Szn
为估计概率 p做
模拟,却又需要
用 p去估计模拟
次数 n.
如何计算 S2？
解决方法,先做 n0 次模拟（称为学习样本 )，根
据学习样本,
（ 1） 先求出 p的估计，再估计模拟次数 n:
2
2
)1(.1
??
?? zppn
2
2
)?1(?
?
?
?? zppn
（ 2） 计算出的样本方差 S2,用来估计 n.
2,M － C模拟的估计精度 ε与试验次数 n的平
方根成反比,若精度 ε提高 10倍,则试验次数 n
要增大 100倍,
P197表 8.2中列出了置信度为 0.95 时,在不同
精度 ε 及概率 p条件下频率法所需试验次数。
对该表进行分析，能得到什么结论？
1,精度提高，试验次数大幅提高；
2,事件发生概率越接近 0.5，试验次数越高；
核反应堆屏蔽层是用一定厚度的铅包围反应
堆，用以阻挡或减弱反应堆发出的各种射线,
在各种射线中,中子对人体伤害极大，因此，
在屏蔽层的设计中,了解中子穿透屏蔽层的概
率对反应堆的安全运行至关重要,
例 7.3.2 核反应堆屏蔽层设计问题
1.问题背景
假定屏蔽层是理想的均匀平板
一个中子进入屏蔽层后运动的物理过程,中
子以初速度 v0和方向角 α 射入屏蔽层,运动一
段距离后与铅核发生碰撞，中子获得新的速度
及方向 (v1,θ1),再游动一段距离后,与铅核发生
第二次碰撞,并获得新的状态 (v2,θ2),如此等等,
经过若干次碰撞后,出现下述情况之一时中子
终止运动过程
1） 中子被弹回反应堆；
2）中子穿透屏蔽层；
3） 第 n次碰撞后,中子被屏蔽层吸收,
D
返回
穿透
吸收
三种
状态
为使屏蔽层的厚度达到安全设计要求,在计
算机上对中子在屏蔽层的运动过程进行模拟
阐述中子的运动，
为模拟做理论准备
2,简化假设：
*1 假定屏蔽层平行板厚度为 D=3d,其中 d
为两次碰撞之间中子的平均游动距离；
*2 假设在第 10 次碰撞以后，中子速度下降
到为某一很小数值而终止运动（被引收）,
因每次碰撞后,中子因损失一部分能量而速
度下降,
*3 假定中子在屏蔽层内相继两次碰撞之间
游动的距离服从指数分布 ；
*4 中子经碰撞后的弹射角 θ～ U(0,2π).
思考,请仔细分析以上假设的合理性,
3,中子运动的数学描述
引进变量,
弹射角 θi— 第 i 次碰撞后中子的运动方向与
x 轴正向的夹角,
xi — 第 i 次碰撞后中子所处位置与屏蔽层内
壁的距离,
R i — 中子在第 i 次碰撞前后的游动距离,
D
0 x
xi
θi
三个变量
均为随机
变量
中子在屏蔽层里随机游动,第 i 次碰撞以后,
按照它的位置坐标 xi,可能有以下三种情况
发生：
（ 1） xi<0,中子返回反应堆 ；
（ 2） xi＞ D,中子穿透屏蔽层；
经过第 i 次碰撞，中子在屏蔽层内的位置是
xi=xi－ 1+Ricosθi, i=1,2,…, 10,
（ 3） 0＜ xi＜ D,若 i<10,中子
在屏蔽层内继续运动,若 i=10
中子被屏蔽层吸收,
中子三
状态判
别准则
4,模拟过程
(1) 产生 RND随机数对 (ri,ui )；
(2) 将 (ri,ui )代入公式计算
?
?
?
???
??
,2
,ln
ii
ii
u
rdR
第 i 次中子
的移动距离
和弹射角
(i=1,2,3,…, 10)
(3) 将 (Ri,θi) 代入公式
xi = xi－ 1+Ricosθi, i=1,2,…, 10
计算出第 i 次碰撞中子与内壁的距离 xi,
(4) 判断中子是否穿透屏蔽层,
5,模拟结果分析
要求穿透屏蔽层的概率数量级为 10－ 6～ 10－ 10,
按假设条件得到 一次模拟结果 如下,
中子数 (个 ) 穿透 (%) 吸收 (%) 返回 (%)
100
1000
3000
5000
30.0
26.0
26.5
26.3
28.0
23.4
21.8
22.0
42.0
50.6
51.7
51.7
中子穿透屏蔽层的百分比超过了 1/4,模拟
结果表明屏蔽层厚度 D=3d不合适,
多厚的屏蔽层才能使穿透的概率
W＜ 10－ 6？
问 题,
如何解决这个问题？ 思路？
1,计算机收索法
增大屏蔽层的厚度，如 D=6d,12d,24d、
36d,…, 交由计算机进行模拟,并搜索到所
求解,
2,分析法
1 2 m
……
D D D
设计屏蔽层的厚度,x=mD
将屏蔽层视为 m层厚度均为 D
的平行板,
由于碰撞的能量损失,中子穿过屏蔽层的平
均速度会逐层下降,
设 WD 是中子穿过厚度为 D 屏蔽层的概率,
则 穿过整个屏蔽层的概率 W满足
mDWW ?
利用模拟结果, 当 D=3d,WD≈0.25,令
(WD)m＜ 10－ 6,或 (m)4＞ 10 6
9 6 7.9
6 0 2.0
6
4l o g
10l o g 6 ???m
取屏蔽层的厚度 x=10D=30d，
可使穿透屏蔽层的概率 w＜ 10－ 6
注,模拟 5000个中子 的运动,用穿透屏蔽的
频率估计穿透概率,由表 8.2可知精度大约只有
1％,模拟精度太低,应适当增大模拟次数,
总结,模拟的意义？
1.模拟方法本质上是 试验性 的,模拟系统是
现实系统的仿真,
例中每模拟一次相当于对一个中子的运动做
一次, 试验, 或, 观察,,
2,是对思维结果的一种验证,
3,模拟本质上是一种求解问题的试验方法,
需要进行较多次数的重复模拟,并且对试验结
果还需进行统计分析,
4,上千万次的模拟计算工作可以借助计算机
完成,而且运算速度是非常快,