蒙特卡罗方法及应用.pdf

《蒙特卡罗方法及应用.pdf》由会员分享，可在线阅读，更多相关《蒙特卡罗方法及应用.pdf（5页珍藏版）》请在三一文库上搜索。

1、蒙特卡罗方法及应用 尹增谦1 管景峰2张晓宏1曹春梅1 ( 1 华北电力大学物理教学部,河北 保定 071003) ( 2 河北省保定市环境保护监测站 ,河北 保定071000) ( 收稿日期: 2002 -01-07) 摘要介绍蒙特卡罗方法的基本原理及其在计算物理中的应用. 关键词蒙特卡罗方法 THE MONTE CARLO METHOD AND ITS APPLICATION Yin Zengqian1 Guan Jingfeng2 Zhang Xiaohong1Cao Chunmei1 ( 1 Department of Physics, North China University o

2、f Electric Power, Baoding, Hebei071003) ( 2 Environmental Detecting Institute of Baoding City, Baoding, Hebei071000) AbstractThe principle of Monte Carlo method and its application are introduced . Key Words Monte Carlo method 1 引言 蒙特卡罗方法, 又称随机抽样方法 ,是一 种与一般数值计算方法有本质区别的计算方 法,属于试验数学的一个分支, 起源于早期的 用几率近似

3、概率的数学思想, 它利用随机数 进行统计试验 , 以求得的统计特征值( 如均 值、 概率等) 作为待解问题的数值解 . 随着现 代计算机技术的飞速发展, 蒙特卡罗方法已 经在原子弹工程的科学研究中发挥了极其重 要的作用 1, 2, 并正在日益广泛地应用于物理 工程的各个方面 , 如气体放电中的粒子输运 过程等 3 5 . 本文介绍蒙特卡罗方法的基本 原理及其在计算物理中的应用 . 2 蒙特卡罗方法的基本原理 就数学特性而言, 蒙特卡罗方法的发展 可以追溯到 18 世纪著名的蒲丰问题. 1777 年,法国科学家蒲丰( Buffon) 提出用投针试验 计算圆周率 值的问题. 这里我们用蒲丰问 题来

4、初步说明蒙特卡罗方法的基本原理和解 决问题的基本手续. 蒲丰问题是这样一个古典概率问题: 在 平面上有彼此相距为 2a 的平行线, 向此平面 任意投一长度为2l 的针, 假定 l a ,显然,所 投的针至多可与一条直线相交 ,那么 , 此针与 任意条平行线相交的概率可以求出 , 由下面 的分析可知,此概率与所取针长 2l 、 平行线间 距 2a 有关 ,并且包含有 值 . 在这里, 任投一 针的概率含义有以下三点 : ( 1)针的中点 Ml 在平行线之间等概率落入 ,即 Ml距平行线的 距离 x 均匀分布在区间 0, a 之内 ; ( 2)针与 线的夹角 均匀分布在区间- 2 , 2 之内;

5、( 3)x 与 互相独立 . 如图 1 所示,建立与平行线垂直且原点 45 物理与工程 Vol . 12No . 3 2002 本文工作获华北电力大学“大学物理教学研究”教改项目基金资助. 图 1蒲丰问题的概率分析 在某一条平行线上的 x 轴, 不失一般性 ,假定 针的中心处于图示中的 x 轴上 . 由于对称性, 我们只需分析针中心处在 x ( 0, a) 范围的 情况即可 . 令探针中心的坐标值为 x , 显然, 只有 x l 时才可能发生相交的事件. 我们来 分析在条件 x l 满足时 , 针与线相交的概 率: 只有当 0=arccos x l 时才能相交, 且 相交的概率为 P1= 2

6、arccos x l ( 1) 下面再来分析针中心位置在轴上的分布 ,显 然, 这是一个均匀分布 , 即针中心处于区间 ( x , x +dx) 内的概率为 dP2=dx a ( 2) 这样 ,一次投掷 ,针中心落入( x , x +dx) 且与 线相交的概率为 dP =P1dP2= 2 aarccos x l dx( 3) 则一次投掷,针与线相交的总概率为 P =dP = l 0 2 a arccos x l dx =2l a ( 4) 即: =2l Pa ( 5) 从( 5) 式可见, 可利用投针试验计算 值: 设 投针 N 次 , 其中 n 次针与线相交, 则可用频 率值 n/ N 作为

7、概率P 的估计值 ,从而求得 的估计值为 2l a N n ( 6) 这就是早期的用频率值作为概率近似值的方 法的应用实例, 表 1 是在历史上一些有名的 用投针试验计算值的结果 2 ,其中针长以 a 为单位. 表 1投针试验计算 值的结果 实验者时间( 年份)针长投针次数相交次数的估值 Wolf18500. 85 ,0002 ,5323. 1596 Smith18550. 63 ,2041, 218. 53. 1554 De Morgan, C18601. 0600382 . 53 . 137 Fox18840 . 751 ,0304893. 1595 Lezzerini19010 . 83

8、3 ,4081 ,8083 . 1415929 Reina19250. 54192 ,5208593. 1795 需要指出的是, 上述由投针试验求得 的近似值的方法 ,是进行真正的试验 ,并统计 试验结果, 要使获得的频率值与概率值偏差 小,就要进行大量的试验 , 这在实际中 , 往往 难以做到. 可以设想 ,对蒲丰问题这样一个简 单的概率问题, 若要进行 10 万次投针试验 , 以每次投针、作出是否相交判断并累加相交 次数用时 5 秒钟计算 ,则需用时 50 万秒 , 即 大约 139 个小时 . 那么, 可以设想, 对于象上 述确定条件下的核裂变、直流气体放电中粒 子的输运过程及粒子输运的

9、总效应, 若要用 多次掷骰子的方法近似求出就是不可能的了. 所以,在现代计算机技术出现之前 ,用频率近 似概率的方法 抑或称为雏形时代的蒙特 卡罗方法 并没有得到实质上的应用. 若用数值模拟方法代替上述的真正的投 针试验, 是利用均匀分布于( 0 ,1) 之间的随机 数序列, 并构造出随机投针的数学模型 ,然后 进行大量的随机统计并求得 的近似值. 46物理与工程Vol. 12 No . 3 2002 图 2 用数值模拟方法计算蒲丰问题 如图 2 建立坐标系, 平面上一根针的位 置可以用针中心 Ml的坐标 x 和针与平行线 的夹角 来决定 , 在 y 方向上的位置不影响 相交性质. 任意投针

10、,意味着 x 与 都是任意 取的. 但 的范围可限于 0, , x 的范围可 限于 0, a . 在这种情况下, 针与平行线相交 的数学条件是 x lsin , 0 x a( 7) 其次, 怎样模拟投针呢? 亦即如何产生任意 的 x , . x 在 0, a 任意取值, 意味着 x 在 0, a 上取哪一点的概率都一样, 即 x 的概 率密度函数为 f1( x)= 1 a ,当 0 x a 时 0,当 x 为其他值时 ( 8) 类似的, 的概率密度函数为 f2( )= 1 ,当 0 时 0,当 为其他值时 ( 9) 由此 ,产生任意( x , ) 的过程就变为由 f1( x) 抽样 x ,由

11、f2( ) 抽样 的过程 . 容易得到 x =a1 =2 ( 10) 式中 , 1, 2均为( 0,1) 上均匀分布的随机数 . 只要随机数的均匀性和独立性良好 , 如此构 造的数值模型就很好地模拟了实际试验中的 一次投针, 并用下式判断是否相交且记录统 计结果: s( xi, i)= 1,当 xi lsin i时 0,当为其他值时 如果投针 N 次 ,那么 sN= 1 N N i =1 s( xi, i)( 11) 是相交几率 P 的估计值 . 这样就实现了用数 值方法模拟真正的投针试验 . 用此方法计算 的 的近似值的情况如表 2 所示. 表 2用蒙特卡罗方法计算的的近似值 投针次数10

12、,00020, 000100, 000200, 000 的近似值3 . 1622333. 1379933. 1411793. 141354 表 2 中的计算结果表明, 随着模拟投针 次数的增大, 所计算的 的近似值越来越接 近于其真值. 而要进行这样的数值模拟 ,就需 要很大的计算量, 只有利用计算机才能实现 . 从蒲丰问题可以看出, 用蒙特卡罗方法 求解问题时, 应建立一个概率模型, 使待解问 题与此概率模型相联系 , 然后通过随机试验 求得某些统计特征值作为待解问题的近似 解 . 与此相似 ,在一些物理问题, 如核裂变 、 直 流气体放电等过程中 , 粒子的输运过程及粒 子输运的总效应 ,

13、 也是可以与某些概率过程 联系起来 ,例如,电子与原子 、 分子、离子的碰 撞过程, 实际上就是与碰撞截面有关的概率 过程 , 这样, 从数学物理特征来说, 类似于用 随机投针方法计算 的近似值, 确定条件下 的核裂变 、直流气体放电中粒子的输运过程 及粒子输运的总效应可以用多次掷骰子的方 法近似求出. 随着现代计算机技术的出现和飞速发 展 ,用计算机模拟概率过程, 实现多次模拟试 验并统计计算结果 , 进而可获得所求问题的 近似结果 . 计算机的大存储量 、 高运算速度使 得在短时间内 , 获得精度极高且内容丰富的 模拟结果 . 在历史上,也正是原子弹工程研究 初期阶段的工作 , 为模拟裂变

14、物质的中子随 机扩散, 提出了运用大存储量 、 高运算速度计 算机的要求 , 这也成为当时推动计算机技术 发展的重要动力 . 也就是在第二次世界大战 期间,冯诺依曼和乌拉姆两人把他们所从事 的与研制原子弹有关的秘密工作 对裂变 物质的中子随机扩散进行直接模拟 以摩 纳哥国的世界闻名赌城蒙特卡罗( Monte Car- lo) 作为秘密代号来称呼 . 用赌城名比喻随机 模拟 ,风趣又贴切 ,很快得到广泛接受, 此后, 人们便把这种计算机随机模拟方法称为蒙特 卡罗方法. 需要指出的是 ,正是由于广泛领域 47 物理与工程 Vol . 12No . 3 2002 的物理问题中存在着大量的随机过程,如

15、粒子 间的碰撞等,使得蒙特卡罗方法在计算机物理 和物理工程中得到日益广泛的应用, 并成为沟 通理论与实验研究的一个桥梁. 需要指出的是 , 蒙特卡罗方法不仅在处 理具有概率性质的问题方面获得广泛的应用, 对于具有确定性问题的计算也因其程序简单 等优点获得了广泛的应用. 这里以定积分的计 算简要说明其处理确定性问题的手续 . 对于定积分 s = b a f( x) dx 通过变量替换, 可以转换为下面的形式 s =k 1 0 g( x) dx 其中 g( x) ( 0, 1) ,当 x ( 0,1) 时( 12) 即转换为求积分 1 0 g( x) dx 亦即求边长为 1的正 方形中一个曲边梯形

16、的面积的问题,如图3所示. 图3用蒙特卡罗方法求定积分 我们可以设想这样一种随机投点求定积 分 1 0 g( x) dx 的方法 : 在一个边长为 1 的正方 形上并以其两边分别为坐标轴画出 曲线 g( x) ,实际上就是图 3 ,然后随机地正方形投 掷小球,那么, 小球击中 g( x)曲线下部分的 概率就等于所要求的积分 1 0 g( x) dx , 这样就 将确定性的定积分问题转化为一个概率问 题,同样可以通过数值模拟方法 蒙特卡 罗方法求得其近似解. 用此方法, 我们计算了 积分 1 0 x2dx , 当投球数为 1 万次时 , 得到的积 分近似值为 0. 332800 , 与其真值 1

17、 3 极为接 近 . 3蒙特卡罗方法的解题手续和特点 在用蒙特卡罗方法解算问题时 , 一般需 要这样几个过程: 构造或描述概率过程 ,对于 本身就具有随机性质的问题 , 如粒子输运问 题 ,主要是正确地描述和模拟这个概率过程. 对于本来不是随机性质的确定性问题, 比如 计算定积分、解线性方程组、偏微分方程边值 问题等, 要用求蒙特卡罗方法求解, 就必须事 先构造一个人为的概率过程 , 它的某些参量 正好是所要求问题的解 . 即要将不具有随机 性质的问题, 转化为随机性质的问题. 这构成 了蒙特卡罗方法研究与应用上的重要问题之 一 . 然后建立各种估计量, 使其期望值是所要 求解问题的解 . 最

18、后根据所构造的概率模型 编制计算程序并进行计算, 获得计算结果. 与其他的数值计算方法相比 , 蒙特卡罗 方法有这样几个优点: ( 1)收敛速度与问题 维数无关 , 换句话说, 要达到同一精度, 用蒙 特卡罗方法选取的点数与维数无关; 计算时 间仅与维数成比例 . 但一般数值方法, 比如在 计算多重积分时, 达到同样的误差, 点数与维 数的幂次成比例 , 即计算量要随维数的幂次 方而增加 . 这一特性,决定了对多维问题的适 用性 . ( 2)受问题的条件限制的影响小. ( 3) 程序结构简单 , 在电子计算机上实现蒙特卡 罗计算时 ,程序结构清晰简单 ,便于编制和调 试 .( 4)对于模拟象粒

19、子输运等物理问题具 有其他数值计算方法不能替代的作用. 蒙特 卡罗方法的弱点是收敛速度慢, 误差大的概 率性质. 这一情况在解粒子输运问题中仍然 存在 . 除此之外, 经验证明 ,只有当系统的大 小与粒子的平均自由程可以相比较时, 一般 在10 个平均自由程左右 ,这方法算出的结果 48物理与工程Vol. 12 No . 3 2002 较为满意. 而对于大系统深穿透问题, 算出的 结果往往偏低. 对于大系统, 其他数值方法往 往很适应,能算出较好的结果. 因此 ,已有人将 数值方法与蒙特卡罗方法联合起来使用,克服 这种局限性 ,取得了一定的效果. 随着计算机技术的飞速发展,蒙特卡罗方 法以其独

20、特的优点广泛应用于计算物理和物 理工程领域,对蒙特卡罗方法的推广必将使其 对物理学学科的发展发挥更大的作用 . 参考文献 1 裴鹿成, 张孝泽. 蒙特卡罗方法及其在粒子输运问题 中的应用. 北京: 科学出版社, 1980. 2 方再根. 计算机模拟和蒙特卡罗方法论. 北京: 北京工 业学院出版社, 1988. 3 B. Eliasson et al. M onte Carlo Simulation of Runaway ElectronsinO2/N2Mixture. Herbsttagungder SPG/SSP . 1987, Vol. 60 pp. 241 247. 4 陈俊英. H2/

21、CH4系统 EACVD 动力学过程研究. 河北 大学硕士学位论文, 2000. 5 张玉红. H2/C2H2系统电子助进热丝化学气相沉积动 力学过程研究. 河北大学硕士学位论文, 2001. 作者简介 尹增谦, 男, 河北省定州市人, 1991 年毕业于哈尔滨工 业大学并获工学硕士学位, 现任华北电力大学物理教学部 副教授, 河北大学光学工程专业在职博士研究生, 主要从事 光学工程领域及计算物理领域的科学研究工作和物理教学 工作, 已发表学术论文二十余篇. ( 上接第 44 页) 常用热管的工作温度范围与典型的工作介质及其相容壳体材料表 种类工作介质工作温度/ 相容壳体材料 低 温 热 管 氮

22、 氟里昂-21 ( CHCl2F) 氟里昂 -11 ( CCl3F) 氟里昂 -113 ( CCl2FCClF2) -60 100 -40 100 -40 120 -10 100 铝、 不锈钢、低碳钢 铝、铁 铝、 不锈钢、铜 铝、铜 常 温 热 管 己烷 丙酮 乙醇 甲醇 甲苯 水 0 100 0 120 0 130 10 130 0 290 30 250 黄铜、 不锈钢 铝、 铜、不锈钢 铜、不锈钢 铜、不锈钢、 碳钢 不锈钢、低碳钢、 低合金钢 铜、 碳钢( 内壁经化学处理) 中 温 热 管 萘 联苯 导热姆-A 导热姆-E 汞 147 350 147 300 150 395 147 300 250 650 铝、不锈钢、 碳钢 不锈钢、 碳钢 铜、不锈钢、 碳钢 不锈钢、碳钢、 镍 奥氏体不锈钢 高 温 热 管 钾 铯 钠 锂 银 400 1 ,000 400 1 ,100 500 1 ,200 1 ,000 1, 800 1 ,800 2, 300 不锈钢 钛、铌 不锈钢、因康镍合金 钨、 钽、钼、铌 钨、钽 参考文献 1 商政宋等译. 实用热管技术. 北京: 化学工业出版社, 1988. 2 胡亚范等. 热管技术应用于钻井余热回收的方案设 计. 应用能源技术, 2000, 5. 49 物理与工程 Vol . 12No . 3 2002