二维Ising模型的程序设计.doc

二维Ising模型的程序设计一、课题名称：二维Ising模型的程序设计二、班级和姓名： *三、主要内容：1. 研究的内容和算法：Ising模型最初由Lenz提出和用来作为铁磁性的一个模型。后来成为他的研究生Ising的博士论文的题目。1925年，Ising给出了一维情况下的解，该解显示，在一维情况下，Ising模型没有相变解。1944年，Onsager得到了二维Ising模型的准确解，二维时就有了相变。对于三维，至今还没有严格解，需依靠数值计算得到。物质在外磁场H中的磁场强度M为（1）抗磁体，0，数值很小且随温度反比或与温度无关；（3）铁磁体，在一定相变温度Tc（Curie温度）之下，M不随H作线性变化，具有磁滞回线是磁体物质的在磁场中行为的基本特性，磁化率与外磁场有关。在Tc之上时，铁磁性消失，转变为顺磁性。（4）反铁磁体，温度在Tc之上时是顺磁体，之下时随温度下降而降低。对于二维Ising模型，令：G=Ld为一个d维、共有N个格点的体系，在每个格点i上有一个自旋，可以朝上或朝下的方向。用自旋变量i表示，。在外磁场H中，体系的哈密顿量为：，其中J 为交换关联系数，B表示单个自旋的磁矩，表示只对格点i周围最邻近的给点j求和。 J为正时为铁磁体的模型，各个自旋倾向于同向排列； J为负时为反磁体的模型，各个自旋倾向于反向排列。2.模拟二维Ising模型的步骤：为了方便，令为1。(1) 选择任意一个初始位形Xx1,x2, xN；(2) 按1/N等概率的选取一个格点i，将其自旋反向，得到一个 新的位形Xx1,x2, xN；(3) 利用公式 ，计算能量差E=E(X)-E(X) ，若E0，则再产生一个0,1之间的随机数，如果 ；则位形改变有效，否则位形不变；(5) 返回步骤(2)，进行下一次迭代。（）固定外加磁场，磁化强度随温度的变化源程序#include #include #include #include #include #include #include #include using namespace std;const double pi=3.1415926;const double em=1.0e-12;double randx,randnum;unsigned long randxi=1;const unsigned long randa=16807;unsigned long randm=2147483637;/ random number creation *void random_number() if ( randxi = randm ) randxi=randm-1; randxi=(randa*randxi)%randm; randx=randxi; randnum=randx/randm; return;/end random_number */ main *int main(void) FILE *fp8; fp8=fopen(固定外加磁场，随温度.txt,w+); const int lenx=100,leny=lenx; int spinlenxleny,imlenx,iplenx; double mag0=0.0,mag=0.0; double enj=0.5,magmax=0.0; double temp=0.0,den=0.0,trapro=0.0; int l=0,mcsi=0,mcsmax=0; int i=0,j=0,ix=0,iy=0,iz=0,kx=0,ky=0,kz=0,t=0; int dspin=0,imx=0,ipx=0,imy=0,ipy=0; double H=0.0; for (i=0;ilenx;i+) imi=i-1; ipi=i+1; im0=lenx; iplenx-1=0; /设置边界条件 mag0=lenx*leny; magmax=abs(mag0);/磁矩最大绝对值 coutsetw(15)H= H; mcsmax=pow(10,5); for (i=1;i=100;i=i+1) temp=temp+0.1; mag=mag0; for (ix=0;ixlenx;ix+) for (iy=0;iyleny;iy+) spinixiy=1;/初始状态，设置所有初始磁矩为1 for (mcsi=0;mcsimcsmax;mcsi+) random_number(); kx=int(lenx*randnum); random_number(); ky=int(leny*randnum); imx=imkx; ipx=ipkx; imy=imky; ipy=ipky; dspin=spinimxky+spinipxky+spinkximy+spinkxipy;/ 对格点（ix，iy）周围最近邻的磁矩求和 den=2.0*enj*spinkxky*dspin+2*H*spinkxky; /能量变化，由于磁矩从-1到1或1到-1，所以乘系数2 trapro=exp(-den/temp); random_number(); if (den = 0 | randnum = trapro) mag=mag-2.0*spinkxky; /磁矩的变化 spinkxky=-spinkxky; /自旋反转 else goto label1; label1: continue; fprintf(fp8,%15.6f %15.6f n, temp,mag/magmax); coutsetw(10)OKn; goto labelout; labelout: return 0; / end main *计算结果及具体分析讨论取关联系数为0.5，外加磁场分别取0.5，1.0，1.5，得到磁化强度随温度的变化关系如下图：外加磁场由图可得，在不同的外加磁场下，随着温度的增加，磁化强度降低，这应该是温度增加导致的电子的热运动增加，电子运动变得更加无序，所以磁矩变小。（）固定温度，磁化强度随外加磁场的变化源程序#include #include #include #include #include #include #include #include using namespace std;const double pi=3.1415926;const double em=1.0e-12;double randx,randnum;unsigned long randxi=1;const unsigned long randa=16807;unsigned long randm=2147483637;/ random number creation *void random_number() if ( randxi = randm ) randxi=randm-1; randxi=(randa*randxi)%randm; randx=randxi; randnum=randx/randm; return;/end random_number */ main *int main(void) FILE *fp8; fp8=fopen(固定温度，随外加磁场.txt,w+); const int lenx=100,leny=lenx; int spinlenxleny,imlenx,iplenx; double mag0=0.0,mag=0.0; double enj=0.9,magmax=0.0; double temp=0.0,den=0.0,trapro=0.0; int l=0,mcsi=0,mcsmax=0; int i=0,j=0,ix=0,iy=0,iz=0,kx=0,ky=0,kz=0,t=0; int dspin=0,imx=0,ipx=0,imy=0,ipy=0; double h=-10; for (ix=0;ixlenx;ix+) for (iy=0;iyleny;iy+) spinixiy=1;/初始状态，设置所有初始磁矩为1 for (i=0;ilenx;i+) imi=i-1; ipi=i+1; im0=lenx; iplenx-1=0; mag0=lenx*leny; magmax=abs(mag0); mag=mag0; coutsetw(15)Temp= temp; mcsmax=pow(10,5); for (i=1;i=2000;i=i+1) t=t+1; if(t=1000) h=h+0.02;/磁场从-10增加到10 else h=h-0.02; /磁场从10减少到-10 for (mcsi=0;mcsimcsmax;mcsi+) random_number(); kx=int(lenx*randnum); random_number(); ky=int(leny*randnum); imx=imkx; ipx=ipkx; imy=imky; ipy=ipky; dspin=spinimxky+spinipxky+spinkximy+spinkxipy;/ 对格点（ix，iy）周围最近邻的磁矩求和 den=2.0*enj*spinkxky*dspin+2*h*spinkxky;/能量变化，由于磁矩从-1到1或1到-1，所以乘系数2 trapro=exp(-den/temp); random_number(); if (den = 0 | randnum = trapro) mag=mag-2.0*spinkxky;/磁矩的变化 spinkxky=-spinkxky;/自旋反转 else goto label1; label1: contin