晶粒长大fortran程序
![晶粒长大fortran程序](https://img.360docs.net/img6b/037s7wyy20enniy4twhe-b1.webp)
- 1、下载文档前请自行甄别文档内容的完整性,平台不提供额外的编辑、内容补充、找答案等附加服务。
- 2、"仅部分预览"的文档,不可在线预览部分如存在完整性等问题,可反馈申请退款(可完整预览的文档不适用该条件!)。
- 3、如文档侵犯您的权益,请联系客服反馈,我们会尽快为您处理(人工客服工作时间:9:00-18:30)。
!原子的三维跳跃
integer x,y,z,xyz(1:1000,1:1000),xn(1:6),yn(1:6),zn(1:6),rn real xyzm(1:1000)
write(*,*)"实验天数Jmax,实验次数Imax"
read(*,*) jmax,imax
xn=(/0,0,0,0,1,-1/)
yn=(/0,0,-1,1,0,0/)
zn=(/1,-1,0,0,0,0/)
iseed=rtc()
do j=1,jmax
x=0
y=0
z=0
do i=1,imax
rn=6*ran(iseed)+1
x=x+xn(rn)
y=y+yn(rn)
z=z+zn(rn)
xyz(j,i)=x*x+y*y+z*z
end do
end do
open(1,file="F:\089024352yi.dat")
do i=1,imax
xyzm=0
xyzm(i)=1.0*sum(xyz(1:jmax,i))/jmax
write(1,*) i,xyzm(i)
end do
close(1)
end
!MC单晶长大
use msflib
parameter ir=400,jr=400
integer is(0:ir+1,0:jr+1),tmax,isn(1:8),nstate,t,nr,ix,iy write(*,*)"please input the time step"
read(*,*)tmax
iseed=rtc()
irc=ir/2
jrc=jr/2
r=min(irc,jrc)-10
is=10
is(irc,jrc)=2
open(1,file="f:\089024352er.dat")
do t=1,tmax
do x=1,ir
do y=1,jr
ix=ir*ran(iseed)+1
jy=jr*ran(iseed)+1
isn=(/is(ix-1,jy-1),is(ix-1,jy),is(ix-1,jy+1),is(ix,jy-1), !is(ix,jy+1),is(ix+1,jy-1),is(ix+1,jy),is(ix+1,jy+1)/) e0=count(isn.ne.is(ix,jy))
if (e0.eq.0) cycle
nr=8*ran(iseed)+1
nstate=isn(nr)
e=count(isn.ne.nstate)
rd=ran(iseed)
de=e-e0+7*(nstate-is(ix,jy))+2.5*rd-1.25
if (de.lt.0.0) is(ix,jy)=nstate
isre=setcolor(is(ix,jy))
isre=setpixel(ix,jy)
end do
end do
write(1,*) t, sqrt(1.0*count(is.eq.2))
end do
close(1)
end
!MC多晶长大
use msflib
parameter ir=400,jr=400,nmax=100
integer is(0:ir+1,0:jr+1),tmax,isn(1:8),nstate,
!t,nr,ix0,iy0,ix,iy
integer igv(0:nmax)
write(*,*)"please input the time step"
read(*,*)tmax
iseed=rtc()
! is=0
do i=1,nmax
ix0=ir*ran(iseed)+1
jy0=jr*ran(iseed)+1
is(ix0,jy0)=i
end do
!is(0,1:jmax)=is(imax,1:jmax)
!is(imax+1,1:jmax)=is(1,1:jmax)
!is(0:imax+1,0)=is(0:imax+1,jmax)
!is(0:imax+1,jmax+1)=is(0:imax+1,1)
open(1,file="F:\089024352san.dat")
igv=1
igv(0)=10
do t=1,tmax
! iarea=0
is(0,0:jr+1)=is(ir,0:jr+1)
is(ir+1,0:jr+1)=is(1,0:jr+1)
is(0:ir+1,0)=is(0:ir+1,jr)
is(0:ir+1,jr+1)=is(0:ir+1,1)
do x=1,ir
do y=1,jr
ix=ir*ran(iseed)+1
jy=jr*ran(iseed)+1
isn=(/is(ix-1,jy-1),is(ix-1,jy),is(ix-1,jy+1),is(ix,jy-1),
! is(ix,jy+1),is(ix+1,jy-1),is(ix+1,jy),is(ix+1,jy+1)/)
e0=count(isn.ne.is(ix,jy))
if(e0.eq.0)cycle
nr=8*ran(iseed)+1
nstate=isn(nr)
rd=ran(iseed)
e=count(isn.ne.nstate)
ig=igv(isn(nr))-igv(is(ix,jy))
de=ig+e-e0+2.5*rd-1.25
if(de.lt.0.0) is(ix,jy)=nstate
iii=mod(is(ix,jy),15)
if(iii==0)iii=iii+1
isre=setcolor(iii)
isre=setpixel(ix,jy)
! if(is(ix,jy).ne.10)iarea=iarea+1
end do
end do
! write(1,*)t,1.0*iarea
write(1,*)t,sqrt(1.0*count(is.ne.0))
end do
close(1)
end
!CA单晶长大
use msflib
parameter ir=400,jr=400
integer is(0:ir+1,0:jr+1),tmax,isn(1:8),nstate,t,nr,is0,iy0,ix,jy integer is1(0:ir+1,0:jr+1),isn1(1:8)
write(*,*)"please input the time step"
read(*,*)tmax
iseed=rtc()
irc=ir/2