关于状态分立系统的改进的Metropolis方法

维普资讯

第１９卷第２期　２Ｏ０２年３月　计算物理　Ｖｅ　１．１９，Ｎｏ．２‘．Ｉ．．　ＣＨ［ＮＥＳＥ　ＪＯＵＲＮＡＬ　ＯＦ　ＣＯＭＰＵＴＡＴＩＯＮＡＬ　ＰＨＹＳＩＣＳ　ｍｍ　（文●壤号】ｌｏｏ１．２４６Ｋ（２ｏｏ２）ｏ２．０１０３．０５　关于状态分立系统的改进的Ｍｅｔｒｏｐｏｌｉｓ方法　吴金闪　，　包景东　，　杨展如　（Ｉ．北京师范大学物理系，北豪Ｉ∞８７５；２．北豪师范大学理论物理研究所．北京ｌｆｉＯｇＴＳ）　［摘薹】对状态分立系统的平衡态统计问题给出了一个改进的Ｍｅｔｒｏｐｏｌｉｓ方法，并在伊辛摸壹上与传统方法进行　了对比，得到　一致的结果，效率提高了约２５％．针对伊辛模型的抽样过程给出了其马尔科夫过程所对应的ＭＩ＿　方程．将所提出的方法应用于伊辛模型，相当于对萁Ｍａｓｔｅｒ方程进行了模拟．作为推广，一大类临界动力学横型的　Ｍｕｔｅｒ方程都可以用本方法进行模拟，以用于其暂态过程的研究．　［关ｔ词】改进Ｍｅｔｒｏｐｏｌｉｓ方法；伊辛模型；Ｍｏｎｔｅ　Ｃａｒｌｏ方法　【中圈分类号１　０５４０　［文献标识码】Ａ　曼分布，其中　为广义相空间坐标，可观测量的平　Ｏ　引言　由处理经典统计系综平衡态问题发展而来的　Ｍｅｔｒｏｐｏｌｉｓ方法，以构造马尔科夫链的方式实现了抽　样．同时避开了计算待定的归一化因子．但是在考虑　诸如伊辛模型等状态分立（“自旋”向上或向下）的系　均值可表为　１　０（　）ｅ－　ｄｘ　＜Ｑ（ｚ））＝Ｌ　———一，　Ｊ　ｅ。　ｄｘ　（１）　式中ｐ＝击，　为玻尔兹曼常数，　为绝对温度．　传统的Ｍｅｔｒｏｐｏｌｉｓ方法的精神是利用细致平衡原理　在相空间中对（　）按玻尔兹曼分布抽样，然后计算　统时，我们发现分母是容易计算的，因此简略了Ｍｅ．　ｔｒｏｐｏｌｉｓ方法通过两两比较其转移率来判断的技巧，　保留了其构造马尔科夫链的精神，直接对所有的态　求分布和抽样．成为一个改进的Ｍｅｔｒｏｐｏｌｉｓ方法．　对于真实系统，其非平衡演化是通过散射、漂移　平均值（Ｄ＞．采用的抽样过程是　①对在　（　）态的系统，随机取一新态ｊ（　）；　②若△　（　—Ｊ）＜０，则跃迁　—Ｊ成功；　等过程实现的，而不是由细致平衡原理决定的热力　学肤迁．但是在近平衡系统．其平衡态的性质是相同　的，而且认为实现平衡的过程和热力学跃迁过程是　等价的．这正是近平衡临界动力学的研究方法之一．　另一方面，除了传统的Ｍｅｔｒｏｐｏｌｉｓ方法之外，完全可　能存在另外一条途径也能实现准演化，而且在这一　方法中，向低能态的跃迁可能被拒绝，向高能态的跃　迁可能被接受．而我们所改进的方法正是如此一个　更加对称的准演化方式．在此基础上，本文利用和发　展了一套描述这种马尔科夫过程的主方程，井用改　进的Ｍｅｔｒｏｐｏｌｉｓ方法进行了模拟．　③若△　（１一』）≥０，则取随机散｝∈（Ｏ，１），若｛　≤ｅ。　，则跃迁卜－　成功，否则系统保持原态；　④重复①．　可以证明，上述抽样的渐近态满足玻尔兹曼分布“　．　现在，我们考虑一个有Ｌ个状态参量｛ｓ　ｌ（１＝　１，２，…，Ｌ）的系统，若每一ｓ　有ｊＩｆ　个分立的状态　Ｉ　：（　＝１，２，…，ｊＩｆ。），则在每一元过程只有　一个　参量发生变化（可能会有ｓ　，ｓ　两个参量相耦台的　情况，此时需对方程略作改动）的近似前提下，其分　布函数向平衡态的演化服从如下Ｍａｓｔｅｒ方程，记Ｓ　的取值为　．　（　，…，　－，　）＝　１　分立状态特征与改进的Ｍｅｔｒｏｐｏｌｉｓ　方法　哈密顿量为　：　（　）的经典系统服从玻尔兹　【收蕞日期】２咖一ｏ７—２０；［修回日期］２００１—０１—２０　Ｊ．ｉ　Ｊ　∑【　Ｐ（‘“，　，…；Ｉ）一　；Ｐ（…，　：，…；‘）】，　（２）　［基量琦置　国家自然科学基盘（１９７７５ＣＣ８）及教育部博士基盘（１９９９００２７０７）资助项目　［作毒柚介】吴金闪（１９７７一），男、浙江兰溪，硬士生，占Ｌ事理论物理及统计物理方面的研究　

维普资讯

ｌ０４　计　算　物　理　第ｌ９卷　其中　是ｓ由　：状态跃迁至　的转移速率　下面对跃迁几率　：的具体形式给一个简捷的　推导：系统哈密顿量为　＝　（Ｓ）．系统处于状态｛　＝　｝的几率由玻尔兹曼分布决定，并满足细致平衡　（　｜，…，一，…，　‘；ｔ）＝　∑［　（一　一０．ｔ）Ｐ（…．一　．…）一　（　一一　）Ｐ（…，　，…）］．　以一维为例，由上述定理得　：（７）　原理．设在某一元过程中Ｓ由　＝　态变为Ｓ　＝　：　惫．而　：可取｝　Ｉ的所有值．为了求出跃迁Ｉ　：一　｝末态的概率分布可先求转移速率的两两比值　，南，　（８）　然后对ｓ的所有可能归一化，而得到　（　：一　），即　山　其中＾＝ｅ２　一　”　．　一＝南　（９）　（＊　一ｘ：）　ｅ　：　．ｅ　‘　）　：～：’　－　ｅ　　）而用传统的Ｍｅｔｒｏｐｏ￣ｓ方法可得　ｅ一　可’　船量．可进一步把Ｅ（　：）写成　Ｅ（　）＝△　（　：）＋Ｅ。，　这里风是总能量中与Ｓ．无关的能量．因此　蚰　（　ｆ）　，¨　¨　”　；ｉ　　：　Ｔ　：　：；＝　：　　：隔　。：　．　．；　其中Ｅ（　：）为系统在　＝　：，其它状态量不变时的　（４）　同理　★　１，　，　～，　＋ｔ＝１．　尽管∑　埘（　一　）≠１，不满足归一化，但埘（　一　≈）已经具有跃迁几率的意义，能够用来模拟演化过　程了．如果要作归一化的话，那么折合跃迁几事为　（　）＝　一，　－　（５）　一　～　ｌ　’　¨　即　’　，　“　丽¨　丽ｌ　’，　（１０）　至此问题转化为对～个可列状态的马尔科夫过　程的模拟．改进Ｍｅｔｒｏｐｏｌｉｓ方法的抽样过程如下：　①对｛　Ｉ：　：｝的系统，对某一状态参量Ｓ．计算　各个状态　：的相关能量△　（　）和状态占有率　（　：），得到各个可取值的分布；　②取一均匀随机数　∈（０．１），若∑　（　）＜　＜　＋　ｉ＿　■ｉ，　・＝　＿　与（８）式和（９）式相同．可见在伊辛模型上两种方法　等价．由（８）式和（９）式可制定新的演化过程Ｍｏｎｔｅ　Ｃａｒｌｏ步，如下：　∑　（　：），则取新状态为　，即对分布实现离散抽　样；　③对其它状态参量抽样重复①　①取格点ｓ．，计算近邻格点的白旋之和，并计　算　得到（８）式和（９）式；　②取一随机数　∈（０，１），若｛＜１０１（　）（郾（８）　式），则ｓ　＝１，否则Ｓ＝一１；　③重复①．　对比前面所述的传统Ｍｅｔｒｏｐｏｌｉｓ方法的抽样过　程可见，在改进后的Ｍｅｔｒｏｐｏｌｉｓ方法中，每一单Ｍｏｎｔｅ　Ｃａｒｌｏ步步调用了一次用来产生预选态的随机数，也　少进行了一次判断，可望提高效率．　作为算例，本文对３７　ｘ３７×３７的简立方格点进　行了计算，其临界温度　满足　Ｔ　一～０．２２１．‘　＝２伊辛模型及计算结果　对于铁磁伊辛模型…，其哈密顿量为　Ⅳ＝一，　ｄ．・　，　（６）　＝（１１）…　在图１（ａ）和（ｂ）中我们计算的是平均磁化强度（Ｊ】ｌｆ）　｛　∑　）－图中　＝　．特别在图ｌ（ｂ）中，同　式中（ｉ，Ｊ）表示最近邻，有Ｌ个状态参量　，每一参　量有两个可能取值　：±１．其Ｇｌａｕｂｅｒ机制　准演化动力学Ｍａｓｔｅｒ方程为　。　的　时对时间与系综进行了平均，这是因为在临界点附　近，涨落较大．我们对临界点附近散了仔缩的搜索，　发现对步长为Ｏ．Ｏ∞５的数据，　在０．００１的精度　

维普资讯

第２期　吴金闪等：关于状态分立系统的改进的Ｍｅｔｒｏｐｏ￣方法　ｌＯ５　内是正确的．当然继续加细计算可提高临界点精度．　另外从图２也可见远离临界点及其附近的两侧的确　存在着定性区　４．如果加细计算也可从另一方面得　到临界点．　一、—＾　～、　－　一　Ｉ　．　＼、．　扭）　０１　０　２　０．３　Ｚ　０　４　０　５　Ｚ　（ａ）ｐ１１ａｓｅ　ｄｉａｇｒａｍ　ｏｆ日ｓｉｎｇｌｅ　ｃｕｌ￣ｉｃ　ｌａｔｉｔｃｅ［ｓｉｎｇ　ｍｏｄｅｌｆｏｒ　ａ　ｓｙｓｔｅｍ　ｏｆ　ｓｉｚｅ　３７　ｘ　３７　ｘ　３７　（ｂ）ｄｅ诅ｉｌｅｄ　ｐｉｃｔｕｒｅ　０ｆ（ａ）ｉｎ　ｔｈｅ　ｒｅ　＿ｍ　ｎｅａｒｔｈｅ　ｃｒｉｉｔｃａｌ　ｐｏｉＩＩｔ　图１　３７　ｘ　３７　ｘ　３７简立方格点伊辛模型的相图及赫界点附近的放大图　Ｆｉｇ．１　Ｐｈａｓｅ　ｄｉ　ｇｒ∞＿ｏｆ　ａ　ｓｉｎｇｌｅ　ｃｕｂｉｃｌ　ｏｅＩｓｌｎｇｍｏｄｅｌ　ａｎｄｉｔｓ　ｄｅｕｆｉｌｅｄ　Ｐ　¨ｅ　ｎｃｇｔｒｔｈｅ　ｃｒｉｉｔｃａｌ　ｐｏｉｎｔ　｝　【　图２不同温度下平均磁化强度的暂态图像　Ｆｉｇ．２　ＴＩ．ｅ　Ｉ￣ｍｓｉｅｎｔ－ｐ１１ｅ力０Ｉｎｅ螂Ｉ　ｄｉａｇｒａｍ　０ｆ　ｔｈｅ　啪暑ｅ呻　珥ｔ　ｕｎｄｅｒ　ｄｉｆｅｍｔ蛔叩ｒＩｌＩ驿・　

维普资讯

计　算　物　理　第ｌ９卷　图２是　值分别为０．２００，０．２２０，０．２２２，０．２２５　的平均磁化强度的暂态过程图像，对应温度由高到　低．在图２（ａ）中ｚ在临界点左劁略远，即较高温时　曲线在零点附近涨落，完全可以视为噪声；在图２　（ｄ）中ｚ的临界点右侧略远，即较低温时涨落较小，　有非零的稳定值；而在图２（ｂ）和图２（ｃ）中，在临界点　附近，涨落较大，同时存在小平台，且驰豫时间非常　大．这说明本方法可以成功地用于暂态过程的研究．　与传统算法对比，我们发现对于多态系统，必须　增加调用随机数的次数，这样不仅增加了机时．而且　对随机散提出了更高的要求．而在我们的方法中每　一加，而本方法的调用随机数与参量效无关。可以预　期，状态参量越多，每一状态参量的可取值麓多，本　方法的优越性就会越显著．而且改进的Ｍｅ，￣ｉｓ方　法对２维正方格点给出了与理论值ｋ　１　＝　２．２６９Ｊ（即Ｚ　＝　—　＝ｏ．４４０７）符合的临界点ｚ。　０．４４０，见图３（ａ）；对３维简立方格点的３态模型　也给出了与物理图景一致的相图：由于有Ｓ＝０的　态的参与，临界温度降低，见图３（ｂ）．　裹１　两种方法在３个模型上计算时闫的对比　Ｔａｂｌｅ１　Ｃ．￣￣ｐａｒｉｓｏｎ　ｏｆｔｈｅｄｍｎｌ＿岫ｉｔｍｅ　Ｏｆ曲ｅ　ｉｍｐｃｏｖｅｄ　ａｎｄ　ｔｈｅ　ｌｒａｄｉ￣ｉｏｎａｌ　Ｍｅｔｒｏｐａ￣ｍｅｆｌｍｄｓ　蚰ｔｈ眦ｄ胁ｒｅｎｔ　ｍｏｄｅｌｓ　单Ｍｏｎｔｅ　Ｃａｒｌｏ步只需要调用一次随机数．在对３　个模型的模拟中，我们发现整体效率有了很大的提　高．表１是在３个伊辛模型：２维２态（２ｄ２ｓ），３维２　态（３ｄ２ｓ）与３维３态（３ｄ３ｓ，一个自旋可以取零的“伊　辛”模型）上两种方法的效率对比．一般说来，由于传　统方法调用随机数的次数随参量数目的增加而增　Ｚ　Ｚ　（ａ）ｐｈａｓｅ　ｄｉａｇｒａｍ　ｏｆ　ａ　ｓｑｕａｒｅ　ｌａｔｔｉｃｅ　Ｉｓｉｎｇ　ｍｏｄｅｌ　ｆｏｒ　ａ　ｓｙｓｔｅｍ　ｏｆ　ｓｉｚｅ　１００　ｘ　１００　（ｂ）ｐｈａｓｅ　ｄｉａｇｒａｍ　ｏｆ　ａ　ｔｈｒｅｅ．８Ｌ４　ｔｈｒｅｅ．￣ｍｅｎｓｉｏｎａｌ　ｓｙｓｔｅｍ　图３　２维２态与３维３态“伊辛”模型的相图　Ｆｉｇ．３　Ｐｈａ￥￣ｄｉａｇｒａｍ　ｏｆ　ｌｓｉｎｇ　ｍｏｄｅｌ　此外，再看方程（５），当系统的状态可取连续值　时，求和化为积分，此时分母不能简单求出，重新回　到了原来的Ｍｅｔｒｏｐｏｌｉｓ方法．可见原来的Ｍｅｔｒｏｐｏｌｉｓ　方法只取相邻两态之间比较的方法正是避过了求分　母的积分．众所周知，Ｍｅｔｏｒｏｌｐｉｓ方法的精髓有两点：　其一是构造一个马尔科夫链来实现抽样；其二是用　。　

维普资讯

第２期　吴金闲等：关于状态分立系统的改进的№ｔｌ　方珐　ｌ０７　其特殊的两两比较的方法避免了分母下无穷多个态　的求和．但是对于状态分立系统，我们却可以直接把　分母算出来而不必局限于原来的Ｍｅｔｒｏｐｏｌｉｓ方法．而　用元过程来构造准演化的思想却正是源于Ｍｅｔｏｐｏ－ｒ　【参考文●ｊ　Ｂｉｎｄｅｒ　Ｋ．Ｈｅｅｎｕａｎｎ　Ｄ　Ｗ．Ｍｏｎｔｅ　Ｃａｄｅ　ｓｉ－柚＿ｌｉ呻ｉｎ　Ｓ　・　ｉｆｓｔｌｃｅｌ　ＰＩＩｐ．啷［Ｍ］．Ｓｐｄ．ｅｒｇ－Ｖｅｒｌｑ￣；１９ａＢ．２７．　扬晨如．分形物理学［Ｍ］．上海ｔ上海科学教育出版　社，１９９６．９５．　Ｇ］ａｖｂｅｒ　Ｒ　Ｊ．Ｔｉｍｅ－ｄｅｐｅｎｄｅｎｔ　ｌ噍　∞０ｆ　ｔｈｅ　ｈｂｓｓ　ｍｏｄｅｌ　Ⅱｓ方法的精神．　３结语　本文以伊辛模型为例，用本文改进的Ｍｅｔｒｏｐｏｌｉｓ　方法处理了一类描述准平衡演化的Ｍａｓｔｅｒ方程，得　到了合理的暂态过程和临界现象．模型的研究对象　可以扩展为一大类近平衡临界动力学模型．对于连　【Ｊ］．Ｊ　Ｍｎｔｈ　Ｐｈｙ，，１９６３．４：２９４．　Ｋａｗａｓａｋｉ　Ｋ．Ｄｉｆｕｓｉｏｎ　ｅｏｎｓｔａｎｔ　ｎｅａｒｔｈｅ口　ｐｏｌｎｔ　ｔｉｍｅ・ｄｅｐｅｎｄｅｎｔ［ｓｉｎ８　ｍｏｄｅｌ［Ｊ］．Ｈ　Ｐ，ｅｖ，１９６６，１４５：２２４．　Ｚｈｕ　Ｊｉａｎ－ｙａｎｇ，Ｙａｎ８　Ｚ　Ｒ．Ｅｘｉｔｄ　０ｆｔｈｅ　ｋｉｎｅｔｉｃ　Ｇａｕｓｓｉａｎ　ｍｏｄｅｌ［Ｊ］．Ｐ　Ｂ　Ｒｅｖ，１９９９，髓：１５５１．　ＺＩＩＩｌ　Ｊｉａｎ・ｙａｎｇ．Ｙａｎｇ　Ｚ　Ｒ．Ｓｏｌａｂｌｅ　ｌ・　ｅ　Ｃ．￣ｍｄａｎ　ｍｏｄｅｌ　续取值模型，本文所提出的方法可以退化到传统的　Ｍｅｔｒｏｐｏｌｉｓ方法；对离散的多态模型，本方法能大量　节约随机数，进而提高了计算的效率．可以预期本方　法对多参量模型将显示出更大的优越性．　；｝［［　ｉｎ　ｅｘｔｅｒｎａｌｆｉｅｌｄ［Ｊ］．Ｐｈｙｓ　ｖ．硼学．２ＯＯＯ．　，聪：２１０．　包景束蘩特卡罗方｛击教程［Ｍ］．北京：北京舜范大　王克纲．Ｇｉｂｂｓ　Ｅｎｓｅｍｂｌｅ　Ｍｏｎｔｅ－Ｃａｒｌｏ方珐及应用［Ｊ］．　８　ｌ　２　３　４　５　６　７　计算物理，２０Ｃ０（４）．　ＡＮⅡ　０ＶＥＤ　Ｍ　ＴＲｏｐ０ＬＩＳ　Ａ　Ｐ】Ｐ】　ＣＨ　Ｔｏ　Ｔ丑置　ＤＩＳＣＲＥＴＥ－ＳＴＡＴＥ　ＳＹ＿ｓＴＥＭ　ＷＥ　Ｊｉｎ－ｓｈａｎ‘　，ＢＡＯ　Ｊｉｎｓｏｄｏ　、ＹＡＮＧ驰ｔｌ卜ｒｕ　（１　忡－眦ｏｆＰｈｙｓｉｃｓ，Ｂｅｒｇ　脯　椭时，＆曲　时，　撕１００８７５，ＰＲ　ｃｈ抽｛　１００８７５．ＰＲ　Ｃｈ／ｎａ）　２　Ｉｎｓ￣／ｔｔａｅ　ｏｆＴｈｅｏｒｅｔｉｃａｌＰｈｙｓｉｃｓ，＆　Ｎｏｒｍａｌ　［Ａｌｍｔｒａｅｔ］Ａｎ　ｉｍｐｒｏｖｅｄ　Ｍｅｔｒｏｐｏｌｌｓ　ａｐｐｒｏａｃｈ　ｉｓ　ｐｒｏｐｏｓｅｄ　ｅｓｐｅｃｉａｌｌｙ　ｏｔ　ｓｌａｖｅ　ｔｈｅ　ｅｑｕｉｌｉｂｒｉｕｍ　ｐ￣ｏｂ［ｓｍ　０ｆ　Ｈ　ｄｉｓｃｔｘ￣ｅ－ｓｔｔｎｅ　ｓｙｓｔｅｍ．ｕｓｅｄ　ｔｏ　ｄｅａｌ　ｗｌｌｈ　Ｉｈｅ［ｓｉｎｇ　ｍｏ＆．１　ａｎｄ　ｐｒｏｖｅｎ　ｔｏ　ｂｅ　ｗｅｌｌ　ｃｏｎｓｉｓｔｅｎｔ　ｗｌｌｈ　ｔｈｅ　ｔｒａｄｉｔｉｏｎａｌ　Ｍｅｔｒｏｐｏｌｉｓ　ｍｅｌｈｅｄ．Ｂｕｔ　ｉｔ’ｓ　ｍｏｌｅ　ｅｆｉｃｉｆｅｎｔ∞ｉｔ　ｒｅｄｕｃｅｓ　２５　ｐｅｒｃｃａｌｌ　ｔｈｅ　ｓｌｍ，ｄ￣ｔｉｍｅ．Ｅ￣ｐｅｃｉａｌｌｙｆｏｒｔｈｅ［ｓｉｎｇｍｏｄｅｌ，ｔｈｅＭａｓｔｅｒ　ｅｑｕａｔｉｏｎｉｓ　ｇｅｎｅｒａｌｉｚｅｄｉｎｔｈｅ　ｓｉｍｕｌａｔｉｏｎ．Ｔｈｉｓ　ｋｉｎｄ　ｏｆＭａｓｔｅｒ　如　ｃｒｉｔｉ・　ｃａ］ｋｉｎｅｌｉｃｓ　ｎｄ　ａｉｔｓ’￣ａｎｓｉｅｎｔ　ｐｈｅｎｏｍｅｎａ　ｃ丑１１　ｅｆｆｅｃｔｉｖ由ｂｅ　ｓｔｕｄｉｅｄ　ｂｙ　ｕｓｉｎｇ　ｔｈｅ　ｍｏｄｉｉｆｄ　ｅＭｅｔｒｏｐｏｌｉ￣ｍｅ￣ｏｄ　ｐ￣ｗｅｎｌｄ　ｏｈｅｌｅ．　［Ｋｅｙ　ｗｅｌｄｓ］ｈｎｐｒｏｖｄ　ｅＭｅｔｌｏｐｏｌｌｓ　ａｐｐｒｏｅｃｈ；Ｉｓｉｎｇ　ｍｏｄｅ］；Ｍｏｎｔｅ　Ｃａｒｌｏ　［ｉｈ￣ｄｗｄ　ｄａｔｅ］２０００—０７一如ｉ［Ｒｅｖｌ￣ｄａｔｅ］２００１一Ｏｌ一如