蒙特卡罗计算三重积分

蒙特卡罗方法计算三重积分平均值

摘要：本文是我读过文献[1]后学习到的一些知识和技巧。通过阅读文献我对蒙特卡罗方法有了一些新的认识，并对利用蒙特卡罗方法求三重积分有了较好的掌握。

1、蒙特卡罗方法

人类通过算法使计算机完成复杂的工作。人们通过算法设计和程序编制是计算机具有了“智能”。一般的算法都是确定性算法，即每一工作步骤都是确定的。还有一些算法叫概率算法，比如说舍伍德算法、拉斯维加斯算法，以及本文将要介绍的蒙特卡洛算法。这种算法允许在执行过程中随机选择下一个计算步骤。在某些问题中，随机性算法比确定性算法更有效率。这种算法的一个特征时，对于一个例子，用该算法计算几次得出的结果可能是不想同的，但是如果我们多计算几次，反复求解，就回降低误差，得到满意的结果。

本文介绍蒙特卡罗方法（随机抽样法，统计实验法）。它使用概率模型进行近似计算。蒙特卡罗方法最早起源于18世纪法国学者蒲丰提出的头针问

题：…….。后来，从事原子弹研制工作的冯·诺依曼和乌拉姆将其以摩纳哥的城市名称“蒙特卡罗”作为秘密代号称呼。

2、蒙特卡罗方法的应用

蒙特卡罗方法在数学、机械、金融、医学等领域的有着广泛的应用，成功解决了许多经典的数学和物理问题（前文提到的诺依曼从事的就是对裂变物质的中子随机扩散进行直接模拟的工作）。蒙特卡罗方法还可以用来求解多重积分，解线性方程组和非线性方程组[2]。下面将主要阐述利用蒙特卡罗方法解三重积分，并给出几个例子。

利用蒙特卡罗方法计算实际问题的基本过程[2]:

1）构造并描述随机概率过程。定积分计算问题本身并不是随机性质的问题，在用蒙特卡罗方法时要首先构造认为的概率过程，将非随机问题转化为随机问题。

2）建立估量，根据随机过程建立某些估计值，即为所求问题的目标解。

3）获取结果，对构造的概率模型进行程序设计，获得预期结果。

第一个过程是很关键的，这个转化过程

3、用蒙特卡罗法解三重积分[1]

均匀随机数计算三重积分可分为平均值法和掷点法两种，在这里我们只考虑三重积分的平均值法。为了将非随机问题转化为随机问题，下面介绍几个定理并证明。

定理1：设f(x,y,z)是区域Ω上的有界连续函数，对于∭f(x,y,z)dxdydz

Ω

1)取包含积分区域Ω的长方体V，V决定于a

2)取任一概率密度函数g(x,y,z)，满足∭g(x,y,z)dxdydz

=1;

Ω

3)(x i,y i,z i)，i=1, … ,N是以g(x,y,z)为概率密度的随机数列,设(x i,y i,z i)，

i=1，…,n为落在Ω中的n个随机数，当N充分大时，有

∭g(x,y,z)dxdydz Ω≈(1

N

)∑f(x i,y i,z i)

g(x i,y i,z i)

n

i=1

。

证明 :

设（X,Y,Z）是以g(x,y,z)为概率密度的随机向量，（X i,Y i,Z i）,i=1,…,N为（X,Y,Z）的一组样本，（x i,y i,z i）,i=1,…,N为其相应的值。

∭g(x,y,z)dxdydz=∭f(x,y,z)

g(x,y,z)

g(x,y,z)dxdydz

V =E(

f(x,y,z)

g(x,y,z)

)

Ω

由大数定理，在给定的概率密度下，随机变量的统计平均值的极限是其概率平均，即

E(f(x,y,z)

g(x,y,z)

)=∭f(x,y,z)

Ωdxdydz≈

1

N

∑

f(x i,y i,z i)

g(x i,y i,z i)

n

i=1

证毕。

接下来需要建立一定的估量，也就是建立蒙特卡罗估计值的标准误差。这里直接给出结论，即标准误差的值。

D(1

N

∑

f(x i,y i,z i)

g(x i,y i,z i)

n

i=1

)=

1

N

[∭

f(x,y,z)2

g(x,y,z)

Ω

dxdydz−(∭f(x,y,z)

Ω

dxdydz)2]

我们要选择一个概率密度函数g(x,y,z),使得蒙特卡罗估计值的标准误差最小。为此我们介绍定理2。

设f(x,y,z)是区域Ω上的有界连续函数，0≤f(x,y,z)≤M,按如下步骤得到计

算∭f(x,y,z)

Ω

dxdydz值的最优的蒙特卡罗算法。

1）去包含积分区域Ω的长方体V，V决定于a

2）取有利概率密度函数g(x,y,z)=f(x,y,z)/c,其中c=∭f(x,y,z)dxdydz

Ω

3）(x i,y i,z i),i=1,…,N是以g(x,y,z)为概率密度的随机数列，设(x i,y i,z i)，i=1,…，

n为落在Ω中的n个随机数，当N充分大时，有∭f(x,y,z)dxdydz≈

Ω

1 N ∑f(x i,y i,z i)

g(x i,y i,z i) n

i=1

4、算例

算例1.用蒙特卡罗方法计算积分∭ln(1+x+y+z)dxdydz

Ω

Ω:0

（注：这个算例是文献原有算例，我重复文献的过程，验证程序的正确性）

用mathematica计算得到的结果如下。

使用蒙特卡罗方法计算得到的计算值为0.8995，精确值为0.895129.误差为0.49%。

算例2：计算积分∭еx+y+z

dxdydz

ΩΩ:0