二维随机格点伊辛模型中相变的蒙特卡罗模拟.pdf

2023 年 6 月伊犁师范大学学报（自然科学版）Jun.2023第 17 卷 第 2 期Journal of Yili Normal University（Natural Science Edition）Vol.17 No.2二维随机格点伊辛模型中相变的蒙特卡罗模拟刘洋，赵薇，周恒为*，黄以能*（伊犁师范大学 物理科学与技术学院，新疆凝聚态相变与微结构实验室，新疆 伊宁835000）摘要：针对二维随机格点伊辛模型（2D-RSIM）至今仍未获得体系热平衡态的精确解，以及已有的蒙特卡罗（MC）模拟结果存在自旋空位的浓度点不够多、系统性不够足等问题，采用MC方法，对模型体系的比热C随温度T变化，特别是所表现的相变行为随的演化，进行了细致的MC模拟.结果表明，对所有的值，C随着T变化存在一个峰，并且随着变大，峰形由尖锐变得圆滑，说明2D-RSIM中发生了弥散相变，并且相变温度Tc随的增大总体上单调降低，细致分析首次发现存在3个Tc随线性变化的区域（I区、II区、III区），其中I区下降斜率最大，III区下降斜率最小.关键词：随机格点伊辛模型；相变；蒙特卡罗模拟中图分类号：O414文献标识码：A文章编号：2097-0552（2023）02-0025-070引言引言为了描述顺磁-铁磁相变，Lenz和Ising提出了现在一般称之为伊辛模型（Ising Model，IM）1-3的著名理论模型4.IM的两个主要假设为：1）伊辛自旋（spin）假设，即自旋只能处于向上或向下两个状态，这里的自旋表示是材料晶体点阵上的原子、分子、离子中电子的自旋磁矩和轨道磁矩的总磁矩，准确地说是赝自旋（Pseudo-spin），一般简称为自旋5；2）最近邻相互作用假设，即晶体点阵上的自旋，只有在最近邻的情况下，才存在相互作用.后来的研究表明，IM模型也是描述顺电-铁电相变最为成功的模型6.Klein-Brout提出的长程相互作用随机格点（random-site）伊辛自旋模型7，是为了描述磁性金属固溶体材料（CuMn1-、AuFe1-等）而提出的.Klein-Brout 模型的 3 个主要假设为：1）自旋和自旋空位（spinvacant）假设，即晶体点阵的任一格点上由自旋（如AuFe1-中的Fe原子）或自旋空位占据，自旋空位定义为无自旋或自旋为0（如AuFe1-中的Au原子）；2）点阵上自旋随机分布假设，即在保证每个点阵格点上的平均浓度为1-条件下，自旋在点阵上随机分布；3）自旋之间长程相互作用假设，即模型体系中，任意一对自旋之间的相互作用为RKKY势8-10.Binder等5用IM的最近邻自旋相互作用假设，对Klein-Brout模型的长程相互作用进行简化，即得到现在一般所谓的随机格点伊辛模型（random-site-IM，RSIM）.Zhang-Huang6发现，在RSIM中再引入随机场收稿日期：2023-04-03基金项目：新疆维吾尔自治区重点实验室开放课题（2021D04015）；伊犁师范大学研究生科研创新项目（YS2022G016）.作者简介：刘洋（1997），男，河南南阳人，伊犁师范大学在读研究生，研究方向：凝聚态物理.*通信作者：周恒为（1968），女，博士，教授，研究方向：凝聚态物理；黄以能（1965），男，南京大学博士生导师，天池特聘教授，研究方向：凝聚态物理.伊犁师范大学学报（自然科学版）2023年（random-field，RF），即RFRSIM，可以对弛豫铁电性（relaxor-ferroelectricity）进行较好的描述.目前，虽然已经获得一维RSIM（1D-RSIM）体系热平衡态的精确解，但是2D-RSIM、3D-RSIM仍未获得相应的精确解5，6.已有的蒙特卡罗（Monte Carlo，MC）模拟结果11-14，则存在所模拟的自旋空位浓度点既不够多、系统性也不够足等问题.针对该问题，本文用MC方法15，对2D-RSIM体系的比热随温度变化，特别是所表现的相变行为，随的演化，进行了细致的计算机模拟.1二维随机格点伊辛模型与蒙特卡罗方法介绍二维随机格点伊辛模型与蒙特卡罗方法介绍2D-RSIM模型的哈密顿量H为H=-J limn i=1n-1j=1ni,ji+1,jri,jri+1,j+i=1nj=1n-1i,ji,j+1ri,jri,j+1.（1）其中，i,j表示二维点阵中第i行第j列格点上的自旋，i,j=1或-1，表示自旋的状态；J为最近邻自旋之间的相互作用能常数；为自旋空位浓度(0 1)；r是0到1之间的随机数；ri,j为随机函数，当r 时，ri,j=0；r 时，ri,j=1；n为沿着晶轴方向的自旋数目；n 表示热力学极限（thermodynamic limit）；0为真空磁导率；为自旋磁矩.对伊辛自旋体系，H就是体系中任意一个自旋取1或-1的自旋构型（spin configuration）的能量.另外，上述自旋体系处于可以交换能量的温度为T的热浴（heat bath）中，即自旋体系和热浴共同构成一个正则系综（canonical ensemble）.当模型体系处于热平衡时，基于Boltzmann-Gibbs统计得，体系中单个自旋的平均内能U、平均比热C为U 11-limn Hn2,(2A)C UT=1kBT21()1-2limn -H2-H2n4.（2B）其中，H1Z He-H,（3A）-H21Z H2e-H,（3B）Z为与H对应的配分函数：Z exp()-H.（3C）表示对体系中所有自旋状态（i,j=1或-1）求和，()kBT-1.将方程（3A）（3B）（3C）代入方程（2B），原则上可以严格计算2D-RSIM的热平衡态的C，但是如上文所述，迄今为止仍然未获得其精确解.MC模拟是一种以概率和统计理论为基础，用计算机实现统计抽样，以获得问题近似解的方法.下文具体介绍利用MC方法，对2D-RSIM体系中，晶体点阵、点阵上自旋随机分布、体系中自旋翻转、体系的热平衡自旋构型、体系的C等模拟过程.因为热力学实际系统包含巨大的粒子数，所以往往可以看作是一个无限大的系统.但是当进行模拟时，因为计算机内存和运行时间的限制，只能选择有限体系.本文选择模拟2D-RSIM体系中所包含的自旋数目26刘洋等：二维随机格点伊辛模型中相变的蒙特卡罗模拟第2期n2()1-=2002.由于有限体系必然存在边界，本文选用的是自由性边界条件，即边界上自旋与模拟体系之外不存在相互作用.二维正方点阵的模拟：生成晶格常数为a、原胞数目为n n的二维正方点阵，如图1所示.点阵上随机自旋的模拟：对应每个格点位置，生成0至1之间的随机数r，若r ，则在格点上模拟放置自旋；若r r时，自旋发生偏转；当p r时，自旋不发生偏转.热平衡自旋构型的模拟：本文选择自旋全部向上的自旋构型为模拟的初始构型，即i,j=1.体系每经过一个MC步后，一般会得到一个新的自旋构型，对第k步体系的一个确定自旋构型的能量为Hk=-J i=1n-1j=1ni,ji+1,jri,jri+1,j+i=1nj=1n-1i,ji,j+1ri,jri,j+1.（5）因为本文模拟的系综是正则系综，所以选择随MC步数变化的Hk作为热平衡的判据参量，即一定的MC步数内Hk的平均值，基本不随MC步数变化，就认为模拟体系趋于温度为T的热平衡态.如图2所示，当=0.0时，对T=2.50J=1.11Tc（Tc将在下文讨论），k 102步以后，Hk的平均值基本不随k变化，即可以认27伊犁师范大学学报（自然科学版）2023年为体系趋于热平衡态；对T=2.26J=1.00Tc，当k 103步以后，体系趋于热平衡态；对T=2.03J=0.90Tc，T=1.83J=0.81Tc，T=1.42J=0.63Tc，当k 20步以后，体系趋于热平衡态.图 2系列温度（T）下，蒙特卡罗模拟的二维随机格点伊辛模型（2D-RSIM）（自旋空位浓度=0.0）的能量（Hk）随蒙特卡罗步数的变化结果图3的相应模拟结果表明，对不同的，当T=0.81Tc时，k 103以后，体系基本都趋于热平衡态.图 3特定约化温度（T Tc=0.81，Tc为相变温度）下，蒙特卡罗模拟的系列自旋空位浓度的二维随机格点伊辛模型（2D-RSIM）的能量（Hk）随蒙特卡罗步数的变化结果由于本文重点模拟2D-RSIM的相变特征，即研究温区主要在Tc附近，k 103步后认为体系基本趋于热平衡是可行的.但是，考虑到不同以及不同T下，体系到达热平衡所需的MC步数不同，步数越大模拟结果会越好.因此本文选择的到达热平衡的MC步数（kE）为50 000步，总的MC步数（kF）为150 000.体系C的模拟：方程2-3所需计算的H、-H2可以用如下的方程近似得到：H11-1n21kF-kE+1k=kEkFHk,（6A）-H21()1-21n41kF-kE+1k=kEkFHk2.（6B）28刘洋等：二维随机格点伊辛模型中相变的蒙特卡罗模拟第2期2模拟结果与分析讨论模拟结果与分析讨论图4为系列自旋空位浓度、自旋数目为n n的2D-RSIM中，单个自旋的平均比热C随温度T变化的计算机模拟结果.当=0.0时，2D-RSIM简化为2D-伊辛模型（2D-IM）.由于2D-IM的C的精确解已经由昂萨格（Onsager）给出3，17，这里作为验证本文模拟结果的依据.首先，C的精确解表明，2D-IM的相变温度Tc=2.27J，本文的模拟结果Tc=2.26J（图4a），比精确解低0.44%，微小的偏差源于本文模拟的是n n=200 200的有限体系的边界效应；其次，C的精确解还表明，其在Tc附近对数发散，即C -ln|T-Tc，除去极为靠近Tc点C值较小外，本文的模拟结果也给出了基本相同的结果（图4a）.因此，本文对2D-IM的计算机MC模拟结果（方程4-6）与精确解结果（方程2-3）偏差较小，可以预期本文对2D-RSIM的模拟结果也具有较高的可信度.图 4系列自旋空位浓度、自旋数目为n n的2D-RSIM中，单个自旋的平均比热C随温度T变化的计算机模拟所得结果.其中，黑线：=0.0，n 的精确解；红色圆点：=0.0，n n=200 200；蓝色圆点：=0.1，n n=211 211；洋红色圆点：=0.2，n n=224 224；海军蓝色圆点：=0.3，n n=239 239；紫色圆点：=0.4，n n=258 258；橄榄绿色圆点：=0.5，n n=283 283；深青色圆点：=0.6，n n=316 316；皇家蓝圆点：=0.7，n n=365 365；橙色圆点：=0.8，n n=447 447；紫罗兰圆点：=0.9，n n=633 633图4的结果还表明，对所有的值，C随着T变化明显存在一个峰，并且随着变大，峰形由尖锐变得越来越圆滑，并且在不断变宽.虽然按照早期的Ehrenfest相变分类法18，圆滑的比热峰表示体系中不存在相变；但是随着研究的深入，尽管仍然存在争议19，20，现在一般认为体系中发生了弥散相变（diffuse phase transi29伊犁师范大学学报（自然科学版）2023年tion）6或者初始相变（incipient phase transition）19.本文作者认为2D-RSIM中发生了弥散相变，并将C的峰值对应的温度称为相变温度Tc（图4a）.由图5可见，Tc随的增大总体上单调降低，细致分析则可以发现存在3个Tc随线性变化的区域（I区、II区、III区），其中I区下降斜率最大，III区下降斜率最小.上述Tc随变化的3个线性区域是本文首次发现，对应的微观机制仍需深入研究.图 5计算机模拟所得的二维随机格点伊辛模型（2D-RSIM）的相变温度Tc随自旋空位浓度的变化结果参考文献：1 HAGGKVIST R，ROSENGREN A，LUNDOW P H，et al.On the Ising Model for the Simple Cubic Lattice J.Adv.Phys.，2007，56（5）：653-755.2 ISING E.Report on the Theory of Ferromagnetism J.Z.Phys.，1925，31（1）：253-258.3 MCCOY B M，WU T T.The Two-Dimensional Ising Model M.New York：Dover Publications，Inc.，2014.4 TARONIA.90 Years of the Ising Model J.Nat.Phys.，2015，11（12）：997-997.5BINDER K，YOUNG A P.Spin-Glasses：Experimental Facts，Theoretical Concepts，and Open QuestionsJ.Rev.Mod.Phys.，1986，58（4）：801-976.6 ZHANG L L，HUANG Y N.Theory of Relaxor-Ferroelectricity J.Sci.Rep.，2020，10（1）：5060-5078.7 Klein M W，BROUT R.Statistical Mechanics of Dilute Copper Manganese J.Phys.Rev.，1963，132（6）：2412-2426.8 RUDERMAN M A，KITTEL C.Indirect Exchange Coupling of Nuclear Magnetic Moments by Conduction Electrons J.Phys.Rev.，1954，96（1）：99-102.9 KASUYAT.ATheory of Metallic Ferromagnetism and Anti-ferromagnetism on Zeners Model J.Prog.Thero.Phys.，1956，16（1）：45-57.10 YOSIDAK.Magnetic Properties of Cu-MnAlloys J.Phys.Rev.，1957，106（5）：893-898.11 DEALCANTARAAA，DE SOUZAAJ F，MOREIRAF G B.Microcanonical Simulation of the Site-Diluted 3-DimensionalIsing-Model J.Phys.Rev.B，1994，49（13）：9206-9209.12 SELKE W，SHCHUR L N，VASILYEV O A.Specific Heat of Two-Dimensional Diluted Magnets J.Physica A，1998，259（3-4）：388-396.13 FULCO U L，NOBRE F D，DA SILVA L R，et al.Investigation of Critical Properties in the Two-Dimensional Site-DilutedIsing Ferromagnet J.PhysicaA，2001，297（1-2）：131-141.14 KUTLU B，GENC A E.The Estimation of a（kTC（P）/J，P）Phase Diagram for a Two-Dimensional Site-Diluted Ising ModelUsing a MicrocanonicalAlgorithm J.PhysicaA，2013，392（3）：451-457.30刘洋等：二维随机格点伊辛模型中相变的蒙特卡罗模拟第2期15 NEWMAN M E J，BARKERMAG T.Monte Carlo Methods in Statistical Physics M.New York：Clarendon Press，2012.16 GLAUBER R J.Time-Dependent Statistics of Ising Model J.J Math Phys.，1963，4（2）：294-307.17 ONSAGER L.Crystal Statistics I：A Two-Dimensional Model with an Order-Disorder Transition J.Phys.Rev.，1944，65（3）：117-149.18 EHRENFEST P.Phase Conversions in a General and Enhanced Sense，Classified According to the Specific Singularities ofthe Thermodynamic Potential J.Proc.Konin.Aka.Weten.Amst.，1933，36（1）：153-157.19 LAVIS DA，KUHN R，FRIGG R.Becoming Large，Becoming Infinite：TheAnatomy ofThermal Physics and PhaseTransitionsin Finite Systems J.Found Phys.，2021，51（1）：90-159.20 BANGU S.Understanding Thermodynamic Singularities：Phase Transitions，Data，and Phenomena J.Phil.Sci.，2009，76（4）：488-505.【责任编辑：张建国】Monte Carlo Simulation of Phase Transitions in a Two-DimensionalRandom Lattice Ising ModelLiu Yang,Zhao Wei,Zhou Hengwei*,Huang Yineng*（Xinjiang Laboratory of Phase Transitions and Microstructures in Condensed Matters,College of Physical Science and Technology,Yili Normal University,Yining,Xinjiang 835000,China）Abstract:In view of the fact that the two-dimensional random lattice Ising model(2D-RSIM)has not yet obtained an accuratesolution to the thermal equilibrium state of the system,and the existing Monte Carlo(MC)simulation results have problems such asthe concentrationpoints of spin vacancy being not sufficient and systematical,this paper uses the MC method to conduct a detailedMC simulation of the specific heat C of the model system as a function of temperature T,especially the evolution of the phasetransition behavior as a function of.The results show that for allvalues,C has a peak as T changes,and asincreases,the peakshapechangesfromsharptosmooth,indicatingthatadiffusephasetransitionoccursin2D-RSIM,andthephasetransitiontemperatureTcdecreases monotonously asincreases.Detailed analysis has found for the first time that there are three regions(regions I,II,andIII)whereTclinearly changes with,with the largest decline slope in region I and the smallest decline in region III.Key words:random lattice Ising model;phase transition;Monte-Carlo simulation31