第２３卷第３期 海南师范大学学报（自然科学版） Ｖｏ１．２３ Ｎｏ．３ ２０１０年９月 Ｊｏｕｒｎａｌ ｏｆ Ｈａｉｎａｎ Ｎｏｒｍａｌ Ｕｎｉｖｅｒｓｉｔｙ（Ｎａｔｕｒａｌ Ｓｃｉｅｎｃｅ） Ｓｅｐ．２０１０ 蒙特卡罗法在二重积分中的改进算法 李满枝，王洪涛，张广路 （海南师范大学数学与统计学院，海南海口５７１１５８） 摘要：应用蒙特卡罗方法计算二重积分，给出算法的详细计算原理和计算流程步骤，在现有的 算法基础上进行改进，并给出实例的具体计算结果．数值模拟结果表明：改进的新算法简化了计 算过程，且模拟的精度较高． 关键词：蒙特卡罗方法；二重积分；数值计算 中图分类号：０ ２４２ 文献标识码：Ａ 文章编号：１６７４—４９４２（２０１０）０３—０２４２—０３ Ｃａｌｃｕｌａｔｉｏｎ ｏｆ Ｄｏｕｂｌｅ Ｉｎｔｅｇｒａｌｓ Ｂａｓｅｄ ｏｎ Ｍｏｎｔｅ Ｃａｒｌｏ Ｍｅｔｈｏｄ ＬＩ Ｍａｎｚｈｉ，ＷＡＮＧ Ｈｏｎｇｔａｏ，ＺＨＡＮＧ Ｇｕａｎｇｌｕ （Ｃｏｌｌｅｇｅ ｏｆＭａｔｈｅｍａｔｉｃｓ ａｎｄ Ｓｔｔａｉｓｔｉｃｓ，Ｈａｉｎａｎ Ｎｏｒｍａｌ Ｕｎｉｖｅｒｓｉｔｙ，Ｈ￣ｋｏｕ ５７１ １５８，Ｃｈｉｎａ） Ａｂｓｔｒａｃｔ：Ｉｎ ｔｈｉｓ ｐａｐｅｒ Ｍｏｎｔｅ Ｃａｒｌｏ ｍｅｔｈｏｄ ｗａｓ ｕｓｅｄ ｔｏ ｃａｌｃｕｌａｔｅ ｄｏｕｂｌｅ ｉｎｔｅｇｒａ１．Ｔｈｅ ｒｅａｌｉｚａｔｉｏｎ ｓｔｅｐｓ，ａｌｇｏｒｉｔｈｍ ａｎｄ ｔｈｅ ｒｅＭｉ￣ｔｉｏｎ ｏｆ ａ ｃｏｎｃｒｅｔｅ ｅｘａｍｐｌｅ ｗｅｒｅ ｗｅｒｅ ｇｉｖｅｎ．Ｔｈｅ ｅｘａｍｐｌｅ ｓｈｏｗｓ ｔｈａｔ ｔｈｅ ｉｍｐｒｏｖｅｄ Ｍｏｎｔｅ－Ｃａｒｌｏ ｍｅｔｈｏｄ ｈａｓ ａ ｂｅｔｔｅｒ ｃｏｍｐｕｔｉｎｇ ａｌｇｏｒｉｔｈｍ ａｎｄ ａ ｂｅｔｔｅｒ ａｃｃｕｒａｃｙ． Ｋｅｙ ｗｏｒｄｓ：Ｍｏｎｔｅ Ｃａｒｌｏ ｍｅｔｈｏｄ；ｄｏｕｂｌｅ ｉｎｔｅｇｒａｌｓ；ｎｕｍｅｉｒｃａｌ ｃａｌｃｕｌａｔｅ 二重积分的传统计算方法要先化成累次积分， 计算过程，节省计算机时，提高计算效率．最后通过 再根据被积函数的原函数进行计算，但若原函数难 算例验证改进算法的可行性与优良性． 以求得，则该积分就无法直接计算，具有很大的局 限性．蒙特卡罗方法为求解二重积分提供了一个新 １二重积分的蒙特卡罗法的改进 的计算方法，该方法利用计算机的快速计算和高精 积分都可以看作是某个随机变量的数学期望． 度的特点来模拟随机投点实验，然后通过概率模 因此，在利用蒙特卡罗方法计算二重积分的时候， 型，由数学期望的计算来得到积分的近似值． 采用了这个随机变量的算术平均值来作为其近似 蒙特卡罗方法以随机模拟和统计试验为手段， 值［３－５］． 从随机变量的概率分布中，通过选择随机数的方法 定理（传统蒙特卡罗方法）设ｆ（ｘ，ｙ）为区域 产生一种符合该随机变量概率分布特性的随机数 Ｄ上的有界函数， 值序列，作为输入变量序列进行特定的模拟试验、 １）取一个包含Ｄ的矩形区域Ｑ：ａ≤ ≤ｂ， 求解的方法¨１ ．在应用蒙特卡罗方法求解二重积 分时，必须要在包含积分域的矩形区域中产生非均 匀的随机点，然后设法转换成我们需要的随机数序 ２）取任一个概率密度函数ｇ（ ，ｙ），满足』ｇ（ ， 列并以此作为数字模拟试验的输入变量序列进行 模拟求解［３ ３． ｙ）ｄｔｘｄｙ＝１； ３）（ ，Ｙ ），ｉ＝１，２，…，ｎ是以ｇ（ｘ，ｙ）为概率密 本文在传统蒙特卡罗方法的基础上做了改进， 直接在积分区域中产生随机点进行计算，由此简化 度的随机数列，设（ ，Ｙ ），ｉ＝１，２，…，ｋ为落在Ｄ 收稿日期：２０１０—０５－１７ 基金项目：海南省教育厅高等学校资助项目（Ｈｊｓｋ２０１０－４１）；海南师范大学青年教师科研资助项目（ＱＮ０９２１）；海南师范 大学应用数学重点学科项目 第３期 李满枝等：蒙特卡罗法在二重积分中的应用 ２４３ 中的ｋ个随机数，则ｒｔ充分大时，有 ，， 一｝ ｋ 定理的证明参见文【３—４】． 应用传统蒙特卡罗方法通常是计算形如ｌ ｌ ｆ（ｘ，ｙ）ｄｘｄｙ的二重积分，其数值计算公式为［ ］ ｒ ６ ｒ ６： ｆ Ｊ ，ｙ）ｄｘｄｙ一 Ｎ【ｉ （ 一 ・ （２） 定理 （改进蒙特卡罗方法）设ｆ（． ｘ，ｙ）为区域 上的有界函数， １）在区域ＪＤ上取任一个概率密度函数ｇ（ｘ，ｙ：｝， 满足ＪＪ ｇ（ ，ｙ）ｄｘｄｙ＝１，为了简单计算，取 占 ｇ（ ）＝｛ Ｊ—｝，。 （ ， ）∈Ｄ ， ｌ０，（ ， ）诺Ｄ 其中Ａ。是区域Ｄ的面积； ２）（ ，Ｙｉ），ｉ＝Ｉ，２，…，／１，是ｌ，Ｘｇ（ｘ，ｙ）为概率密 度的随机数列，设（ ，Ｙ ），ｉ＝１，２，…，，ｌ为落在Ｄ 中的ｎ个随机数，则ｌＤ ＂ｔ充分大时，有 Ｄ 』 其 中 一 ｎ袅 碧一 Ａ 是 ） ） ３）设（ ，Ｙ ），ｉ＝１，区 ２，…，，ｌ为落在Ｄ中的１１，个 随机数，则ｎ充分大时，有 ，６ ｒ ６２（ ） Ｊ。． ，ｙ）ｄｘｄｙ ［高地 （６２ ）一 ））］ ） 证明 根据二重积分的数学定义，形如』 ： ， 占 ｙ）ｄｏ＂的二重积分根据积分区间Ｄ的定义可以化为 类次积分Ｊ ａｌ Ｊ ，ｙ）ｄｘｄｙ，其中Ｄ ≤ ≤６ ， ａ２（．９Ｃ）≤ ≤ｂ２（ ），令 的面积，则满足条件ｊ ｇ（ ，ｙ）ａ￣ｄｙ＝１．则根据 公式（１）得 』 ，，，， １ 耋袅象参 由积分中值定理，曲边梯形ｏ（３ｃ￣边为Ｙ＝ ｂ （ ），下曲边为Ｙ＝ｎ２（ ））的面积 Ａ。＝（ｂｌ一口。）（６：（ ）一ｎ２（ ）），（０ ≤ ≤ｂ ）， 在区间ａ ，ｂ 】内任意插入ｎ一１个分点 ， 。：：ａ。， ＝ｂｌ，ａｌ＜ ＜ｂ１，（，ｉ＝１，２，…，凡一１），当凡很 大时’（６：（ ０２（ ｝ ；（ｂ２ ）一ａ２（ ））， 即 。一｝∑（，‘ｉ＝１ ６。一 ）（６ （ ）一ｎ （ ）），当ｎ充 分大时，Ａ。可近似看作为区域Ｄ的面积． 由上式展开得 ¨ ，）Ｊ Ｊ ， ｙ）ｄｘｄｙ ｝ ，（ ）Ｙｌ Ａ。 ）（６。 ）（６ （ ）一ｎ２㈤）卜 ［ｉ ６２（ ）］ 证毕． 与传统算法相比，改进算法省略了取包含Ｄ的 矩形区域ｎ的步骤，因此也省略－ｊ＇ｍｉｎ（ａ２（ ））＝ｃ， ｍａｘ（ｂ：（ ））＝ｄ的计算，省略了判断随机数落在积 分区域Ｄ的步骤，从而大大降低计算难度，简化了 计算步骤． ２算例 用蒙特卡罗法计算 』 ． 根据公式（４）可知 ２４４ 海南师范大学学报（自然科学版） ２０１０年 』 』 ｙｄ ｄｙ ［ 耋 ， ｃ２一・，（ 一 ）］一 加抽取样本次数，计算量增加，其精度也随之增加， 但是二者不成正比例关系．从图３可以看出：当随 ・）】 计算结果见图１，图２，图３． ． 图１ ９００个随机数的分布图 Ｆｉｇ．１ Ｄｉｓｔｒｉｂｕｔｉｏｎ ｏｆ ９００ ｒａｎｄｏｍ ｎｕｍｂｅｒ 图２ ２０００个随机数的分布图 Ｆｉｇ．２ Ｄｉｓｔｒｉｂｕｔｉｏｎ ｏｆ ２０００ ｒａｎｄｏｍ ｎｕｍｂｅｒ 注图１一图２中的曲线是函数ｙ＝ 在［１。２】区间的曲线， 随机数在曲线之下说明产生的是 Ｅ［１，２】，ｙ∈［１， ２】的 随机数 ３结论 在计算过程采用线性同余器产生均匀随机数， 应用以上步骤求解二重积分得到较好的结果．当增 机数的数目Ｎ＝９９０时，误差ｄ＝０．００９ ２５４，是最 佳数据．再产生更多的随机数，误差反而较大．所以 应用蒙特卡罗法求解二重积分存在最优数据问题． 总的看来，用蒙特卡罗方法求解的算法描述简单， 有较高的计算精度，对于求解二重积分是个简单有 效的方法． 椭 随机数的个数 图３随机数的数目与误差精度的分析圈 Ｆｉｇ．３ Ａｎａｌｙｓｉｓ ｃｈａｒｔ ｏｆ ｒａｎｄｏｍ ｎｕｍｂｅｒ ａｎｄ ｅｒｒｏｒ 参考文献 ［１］徐钟济．蒙特卡罗方法【Ｍ】．上海：上海科学技术出版社， １９８５． ［２］方再根．计算机模拟和蒙特卡罗方【Ｍ】．北京：北京工 业学院出版社，１９８８． ［３］尹增谦，管景峰，等．蒙特卡罗方法及应用【Ｊ】．物理与工 程，２００２，１２（３）：４６—５０． ［４］刘辉玲，叶锋．计算多重积分的均匀随机数蒙特卡罗法 的实现【Ｊ】．电脑知识与技术，２００８，４（８）：２２８９—２２９１． ［５］黎锁平．运用蒙特卡罗方法求解随机性问题【Ｊ】．甘肃工 业大学学报，２００１，２７（２）：９５－９７． ［６］Ｆｉｓｈｍａｎ Ｇ Ｓ．Ｍｏｎｔｅ Ｃａｄｏ—Ｃｏｎｃｅｐｔｓ，Ａｌｇｏｒｉｔｈｍｓ ａｎｄ Ａｐ— ｐｌｉｃａｄｏｎｓ［Ｍ］．Ｎｅｗ Ｙｏｒｋ：Ｓｐｒｉｎｇｅｒ，１９９６． ［７］Ｎｉｅｄｅｒｒｅｉｔｅｒ Ｈ．Ｋａｎｄｏｍ Ｎｕｍｂｅｒ Ｇｅｎｅｒａｔｉｏｎ ａｎｄ Ｑｕａｓｉ— Ｍｏｎｔｅ Ｃａｒｌｏ ｍｅｔｈｏｄｓ［Ｍ］．ＣＢＭＳ—ＮＳＦ Ｒｅｇｉｏｎａｌ Ｃｏｎｆｅｒ— ｅｎｃｅ Ｓｅｒｉｅｓ ｉｎ Ａｐｐｌｉｅｄ Ｍａｔｈｅｍａｔｉｃｓ，６３：ＳＩＡＭ，１９９２． ［８］Ｒｏｂｅｒｔ Ｃ Ｐ，Ｃａｓｅｌｌａ Ｇ．Ｍｏｎｔｅ Ｃａｒｌｏ Ｓｔａｔｉｓｉｔｃａｌ Ｍｅｔｈｏｄ （ｓｅｃｏｎｄ ｅｄｉｔｉｏｎ）【Ｍ】．Ｎｅｗ Ｙｏｒｋ：Ｓｐｒｉｎｇｅｒ，２００４． 责任编辑：毕和平
