经典和量子耗散系统的随机模拟方法

作者简介

内容简介

第一章 随机方法概述
　　随机模拟方法是利用计算机进行数值计算的一类方法，又称计算机模拟方法.该方法是对构造的模型作统计实验，故也称统计实验法或统计模拟法.这个名字表明，这些模拟方法要使用计算机内的“随机数发生器”所产生的“伪随机数”.
　　这种方法同确定性方法并不是截然相反的，布朗(Brown)动力学就是把两种方法结合起来形成一种混合方法的例子，其中系统的漂移(平均）由牛顿(Newton)方程来控制，而它的扩散行为取决于涨落以及噪声.
　　随机模拟方法的应用广泛，可以用来研究各种类型的问题，既能解决随机问题，也可求解确定性问题.若问题中不包含时间因素，称为静态模型，研究此类问题的模拟方法习惯上称为蒙特卡罗方法；而将求解随机性的动态模型的模拟方法称为随机模拟方法.但现在人们已不再将两者严格地区分开来.本书使用随机方法的名称，表明涉及的问题更广泛.
　　1.1 预备知识
　　定义1.1 分子动力学
　　分子动力学方法是计算一组分子的相空间轨道，其中每个分子各自服从经典牛顿定律.
　　定义1.2 布朗动力学
　　布朗动力学方法是模拟构成一个统计系统的N个粒子的轨道，其中每一个粒子都遵守朗之万方程.
　　定义1.3 蒙特卡罗方法
　　该方法使问题的解等于一个假设的统计模型的参数，用随机数序列建立这个统计模型的一个样本，从它可以算出这个参数的估计值.
　　定义1.4 马尔可夫(Markov)链
　　序列，如果对任何的n，条件概率具有如下的性质：
　　(1.1)
　　则称此序列是一个马尔可夫链.
　　定义1.5 弱平稳过程
　　设是一个随机过程，如果对于任意的t和，
　　(1.2)
　　不依赖于t，则是一个弱平稳过程.
　　定理1.1 大数定理
　　若x1， ，xn是按照分布密度函数f(x)抽取的随机变量，设存在，有
　　(1.3)
　　这个定理表明，只要样本足够大，求和平均将趋近想要计算的积分值.即使随机变量是相关的，上面的结论也仍然成立.
　　定理1.2 中心极限定理
　　(1.4)
　　其中，σ2是统计量g(x)在分布密度函数f(x)下的方均差，定义为
　　(1.5)
　　称为置信水平.这个定理告诉人们，如果样本为有限个，蒙特卡罗计算所引起的估计误差为.并且，统计量的估计值偏离真值的程度不大于ε的事件发生的概率为.
　　1.2 蒙特卡罗方法的发展
　　1946年，物理学家冯 诺伊曼(Von Neumann)等在电子计算机上用随机抽样方法模拟了裂变物质的中子链式反应，因为此项研究是与研制原子弹有关的秘密工作，他们把此方法称为蒙特卡罗方法.以赌城的名字作为随机模拟方法的代号，既风趣又贴切，是具有象征意义的.虽然随机模拟方法始于20世纪40年代，即系统地用蒙特卡罗方法解决科学问题，利用早期的电子计算机模拟中子在裂变物质中的输运.但一个最早包含了随机模拟方法思想雏形的实验，可追溯到著名的蒲丰“投针求π”.
　　1777年，法国学者蒲丰用随机投针方法求在地面上画一些间距为2a的平行线束，向此地面随机投长为2l的针，l　　(1.6)
　　其中，x和θ分别在[0，a]和[0，π]两区间上均匀独立地取值.针与线相交概率为
　　(1.7)
　　其中，P( )表示事件 发生的概率.从图1.2也可看出，针与平行线束的相交概率就是曲线下方面积与矩形面积之比.
　　图1.2 曲线下方面积表相交频率
　　设随机投针N次，其中n次针与线束相交，当N足够大，可用频率n/N作为概率P的估计值，从而算得π的值
　　(1.8)
　　历史上，曾有一些学者作了随机投针的实验，得到了π的近似估计值.表1.1列出了部分实验结果.
　　表1.1 圆周率π的实验值
　　我们亦用随机投针实验计算了π值，在图1.3中画出了数值结果与真值误差的绝对值随投针次数的变化规律，进而检验了蒙特卡罗方法误差公式.
　　图1.3 投针求π的误差随的变化关系.其中点为模拟结果，实线是误差公式
　　概率论和随机行走理论的许多进展，可以看做是蒙特卡罗方法的基础.例如，开尔文(LordKelvin)已经用到随机抽样去帮助求气体动能理论中的一些时间积分，他的随机抽样包含了关联系数分布.柯朗(Courant)、弗里德里希斯(Friedrichs)和卢伊(Lewy)认为一定的随机行走与某个偏微分方程的解相等价.20世纪30年代，费米(EnricoFermi)作了一些数值实验，后来被称为蒙特卡罗计算.他在研究当时最新发现的中子的行为时，设计了一个中子可能与凝聚态物质相互作用的抽样实验，这导致了后来更多的中子扩散和输运理论.
　　1.3 蒙特卡罗求积分思想
　　第二次世界大战期间，冯 诺伊曼、费米、乌拉姆(S.Ulam)、米特罗波利斯(Metropolis)以及现代数字式计算机的开拓者，均对蒙特卡罗方法的发展做出了突出的贡献.在20世纪40年代末和50年代早期，蒙特卡罗方法有了一个突飞猛进的发展，大量的文章描写了这种新方法，以及如何用它去解决统计力学、辐射输运、经济模型等领域中的问题.遗憾的是，当时受计算机的功能所限，用计算机模拟的方法，并没有获得比解析结果更为深入的成果.由于计算机的计算速度的快速增长，才使艰巨和复杂的计算得以完成，人们从失败中学到了很多东西.
　　同时，降低误差方法的理论进展，尤其是“重要抽样方法”和1953年米特罗波利斯算法P6]的提出，使蒙特卡罗方法在当时达到了一个很高的水平，也一直影响到现在.
　　在我国，蒙特卡罗方法是伴随着20世纪60年代原子能科学技术的需要而发展起来的W.目前，一些高校和研究所开设了“计算物理学”必修课，仅在其中用一两章来介绍蒙特卡罗方法及应用①.
　　本节用蒙特卡罗方法计算一个定积分的值，来阐述它的基本思想.
　　假设要计算一个光滑函数g(x)在区间[a，b]的定积分
　　(1.9)
　　除了高斯求积公式外，标准的数值方法是计算被积函数在等间隔点Xi上的值，然后求和
　　(1.10)
　　其中，权重Wi并不依赖于g，它们由算法的精度来确定，仅在最简单的矩形求积公式中的每个才等于1.用蒙特卡罗方法求积分，也将用到(1.10)式，不过所有的权重系数叫恒等于1，而且巧是随机选的.
　　若在区间[a，b]上有界，对于[a，b]上的定积分，只需做线性变换，即可化成[0，1]区间上的积分，且被积函数值也能在[0，1]之间.令x=a+(b-a)u，则有
　　(1.11)
　　其中.因，所以.故变为求定义域在[0，1]区间，值域也在[0，1]范围的定积分
　　(1.12)
　　如果积分区间的上下限有一个或两个为无穷大，则可通过一个非线性变换，例如，在[0，oo)区间积分，可考虑拉盖尔(Laguerre)(exp(-x))变换;对于(-oo，oo)区间积分，则实施厄米(Hemiite)(exp(-x2))变换，设法将积分区间转变成为有界，这样一来，才可以产生简单子样.
　　1.随机投点法
　　如图1.4所示，计算定积分值也就是求曲线以下阴影的面积.建立一个随机投点法模型，向正方形，内随机投点，其中和是两个独立的均勻分布随机数.
　　图1.4 在矩形中随机投点计算曲线阴影面积
　　若第i个点落入曲线下方的区域，即满足条件，则称第i次实验成功.随机投点实验成功的概率p为
　　(1.13)
　　重复进行随机投点，记录成功次数M和总实验次数N，显然两者之比等于曲边梯形面积与正方形面积的比，那么定积分(1.12)的MC估计是
　　(1.13)