二维Ising模型的程序设计.docx
《二维Ising模型的程序设计.docx》由会员分享，可在线阅读，更多相关《二维Ising模型的程序设计.docx（12页珍藏版）》请在冰点文库上搜索。
二维Ising模型的程序设计
二维Ising模型的程序设计
一、课题名称：
二维Ising模型的程序设计
二、班级和姓名：
***
三、主要容：
1.研究的容和算法：
Ising模型最初由Lenz提出和用来作为铁磁性的一个模型。
后来成为他的研究生Ising的博士论文的题目。
1925年，Ising给出了一维情况下的解，该解显示，在一维情况下，Ising模型没有相变解。
1944年，Onsager得到了二维Ising模型的准确解，二维时就有了相变。
对于三维，至今还没有严格解，需依靠数值计算得到。
物质在外磁场H中的磁场强度M为
（1）抗磁体，χ<0，数值很小且是常数，不随温度变化；
（2）顺磁体，χ>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）选择任意一个初始位形X{x1,x2,…xN}；
（2）按1/N等概率的选取一个格点i，将其自旋反向，得到一个
新的位形X’{x1,x2,…xN}；
（3）利用公式，计算能量差△E=E（X’）-E（X），若△E<0，则改变有效，位形改变X→X’；
（4）如果△E>0，则再产生一个[0,1]之间的随机数ξ，如果；则位形改变有效，否则位形不变；
（5）返回步骤
（2），进行下一次迭代。
（Ⅰ）固定外加磁场，磁化强度随温度的变化
①源程序
#include
#include
#include
#include
#include
#include
#include
#include
usingnamespacestd;
constdoublepi=3.1415926;
constdoubleem=1.0e-12;
doublerandx,randnum;
unsignedlongrandxi=1;
constunsignedlongranda=16807;
unsignedlongrandm=2147483637;
//randomnumbercreation***********************
voidrandom_number（）
{
if（randxi==randm）
{
randxi=randm-1;
}
randxi=（randa*randxi）%randm;
randx=randxi;
randnum=randx/randm;
return;
}
//endrandom_number*****************************
//main*****************************************
intmain（void）
{
FILE*fp8;
fp8=fopen（"固定外加磁场，随温度.txt","w+"）;
constintlenx=100,leny=lenx;
intspin[lenx][leny],im[lenx],ip[lenx];
doublemag0=0.0,mag=0.0;
doubleenj=0.5,magmax=0.0;
doubletemp=0.0,den=0.0,trapro=0.0;
intl=0,mcsi=0,mcsmax=0;
inti=0,j=0,ix=0,iy=0,iz=0,kx=0,ky=0,kz=0,t=0;
intdspin=0,imx=0,ipx=0,imy=0,ipy=0;
doubleH=0.0;
for（i=0;i{
im[i]=i-1;
ip[i]=i+1;
}
im[0]=lenx;
ip[lenx-1]=0;
//设置边界条件
mag0=lenx*leny;
magmax=abs（mag0）;//磁矩最大绝对值
cout<cin>>H;
mcsmax=pow（10,5）;
for（i=1;i<=100;i=i+1）
{
temp=temp+0.1;
mag=mag0;
for（ix=0;ix{
for（iy=0;iy{
spin[ix][iy]=1;//初始状态，设置所有初始磁矩为1
}
}
for（mcsi=0;mcsi{
random_number（）;
kx=int（lenx*randnum）;
random_number（）;
ky=int（leny*randnum）;
imx=im[kx];
ipx=ip[kx];
imy=im[ky];
ipy=ip[ky];
dspin=spin[imx][ky]+spin[ipx][ky]+spin[kx][imy]+spin[kx][ipy];//对格点（ix，iy）周围最近邻的磁矩求和
den=2.0*enj*spin[kx][ky]*dspin+2*H*spin[kx][ky];//能量变化，由于磁矩从-1到1或1到-1，所以乘系数2
trapro=exp（-den/temp）;
random_number（）;
if（den<=0||randnum<=trapro）
{
mag=mag-2.0*spin[kx][ky];//磁矩的变化
spin[kx][ky]=-spin[kx][ky];//自旋反转
}
else
{
gotolabel1;
}
label1:
continue;
}
fprintf（fp8,"%15.6f%15.6f\n",\
temp,mag/magmax）;
}
cout<gotolabelout;
labelout:
return0;
}
//endmain*************************************
②计算结果及具体分析讨论
取关联系数为0.5，外加磁场分别取0.5，1.0，1.5，得到磁化强度随温度的变化关系如下图：
外加磁场
由图可得，在不同的外加磁场下，随着温度的增加，磁化强度降低，这应该是温度增加导致的电子的热运动增加，电子运动变得更加无序，所以磁矩变小。
（Ⅱ）固定温度，磁化强度随外加磁场的变化
①源程序
#include
#include
#include
#include
#include
#include
#include
#include
usingnamespacestd;
constdoublepi=3.1415926;
constdoubleem=1.0e-12;
doublerandx,randnum;
unsignedlongrandxi=1;
constunsignedlongranda=16807;
unsignedlongrandm=2147483637;
//randomnumbercreation***********************
voidrandom_number（）
{
if（randxi==randm）
{
randxi=randm-1;
}
randxi=（randa*randxi）%randm;
randx=randxi;
randnum=randx/randm;
return;
}
//endrandom_number*****************************
//main*****************************************
intmain（void）
{
FILE*fp8;
fp8=fopen（"固定温度，随外加磁场.txt","w+"）;
constintlenx=100,leny=lenx;
intspin[lenx][leny],im[lenx],ip[lenx];
doublemag0=0.0,mag=0.0;
doubleenj=0.9,magmax=0.0;
doubletemp=0.0,den=0.0,trapro=0.0;
intl=0,mcsi=0,mcsmax=0;
inti=0,j=0,ix=0,iy=0,iz=0,kx=0,ky=0,kz=0,t=0;
intdspin=0,imx=0,ipx=0,imy=0,ipy=0;
doubleh=-10;
for（ix=0;ix{
for（iy=0;iy{
spin[ix][iy]=1;//初始状态，设置所有初始磁矩为1
}
}
for（i=0;i{
im[i]=i-1;
ip[i]=i+1;
}
im[0]=lenx;
ip[lenx-1]=0;
mag0=lenx*leny;
magmax=abs（mag0）;
mag=mag0;
cout<cin>>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;mcsi{
random_number（）;
kx=int（lenx*randnum）;
random_number（）;
ky=int（leny*randnum）;
imx=im[kx];
ipx=ip[kx];
imy=im[ky];
ipy=ip[ky];
dspin=spin[imx][ky]+spin[ipx][ky]+spin[kx][imy]+spin[kx][ipy];//对格点（ix，iy）周围最近邻的磁矩求和
den=2.0*enj*spin[kx][ky]*dspin+2*h*spin[kx][ky];//能量变化，由于磁矩从-1到1或1到-1，所以乘系数2
trapro=exp（-den/temp）;
random_number（）;
if（den<=0||randnum<=trapro）
{
mag=mag-2.0*spin[kx][ky];//磁矩的变化
spin[kx][ky]=-spin[kx][ky];//自旋反转
}
else
{
gotolabel1;
}
label1:
continue;
}
fprintf（fp8,"%15.6f%15.6f\n",\
h,mag/magmax）;
}
cout<gotolabelout;
labelout:
return0;
}
//endmain*************************************
②计算结果及具体分析讨论
取关联系数为0.9，画出温度分别为0.25，0.5时的磁化强度随外加磁场的变化关系，如下图：
温度越低，磁滞回线的面积越大，低温时磁矩更大，所以矫顽力更大，导致磁滞回线面积更大。
升级会员 