当前位置:文档之家› 计算材料学Fortran程序汇总

计算材料学Fortran程序汇总

INTEGER M(1:10000), NUMBER1(0:364), NUMBER2REAL X,YISEED=RTC()DO J=1, 10000NUMBER1=0X=RAN(ISEED)NUMBER1(0)=INT(365*X+1)JJJ=1DO I=1,364Y=RAN(ISEED)NUMBER2=INT(365*Y+1)ETR=COUNT(NUMBER1.EQ.NUMBER2)IF (ETR= =1) THENEXITELSEJJJ=JJJ+1M(J)=JJJNUMBER1(I)=NUMBER2END IFEND DOEND DODO I=1,10000IF(M(I).LE.23) SUM=SUM+1END DOPRINT *,SUM/10000END=================================================Monte Carlo Simulation of One Dimensional DiffusionINTEGER X,XX(1:1000,1:1000)REAL XXM(1:1000)! X:INSTANTANEOUS POSITION OF ATOM! XX(J,I):X*X ,J:第几天实验,I:第几步跳跃! XXM(I): THE MEAN OF XXWRITE(*,*) "实验天数JMAX,实验次数IMAX"READ(*,*) JMAX,IMAXISEED=RTC()DO J=1,JMAX !第几天实验X=0 !!!DO I=1,IMAX !第几步跳跃RN=RAN(ISEED)IF(RN<0.5)THENX=X+1ELSEX=X-1END IFXX(J,I)=X*XEND DOEND DOOPEN(1,FILE="C:\DIF1.DAT")DO I=1,IMAXXXM=0.0XXM(I)=1.0*SUM(XX(1:JMAX,I))/JMAX !!WRITE(1,*) I, XXM(I)END DOCLOSE(1)END============================================================ ! Monte Carlo Simulation of Two Dimensional DiffusionINTEGER X,Y,XY(1:1000,1:1000)REAL XYM(1:1000)! X:INSTANTANEOUS POSITION OF ATOM! XY(J,I):X*Y ,J:第几天实验,I:第几步跳跃! XYM(I): THE MEAN OF XYWRITE(*,*) "实验天数JMAX,实验次数IMAX"READ(*,*) JMAX,IMAXISEED=RTC()DO J=1,JMAX !第几天实验X=0 !!!Y=0 !!!DO I=1,IMAX !第几步跳跃RN=RAN(ISEED)IF(RN.LT.0.25)THENx=xy=y-1END IFIF(RN.LT.0.5.AND.RN.GE.0.25)THENx=xy=y+1END IFIF(RN.LT.0.75.AND.RN.GE.0.5)THENx=x-1y=yEND IFIF(RN.GE.0.75)THENx=x+1y=yEND IFXY(J,I)=X*X+Y*YEND DOEND DOOPEN(1,FILE="C:\DIF2.DAT")DO I=1,IMAXXYM=0.0XYM(I)=1.0*SUM(XY(1:JMAX,I))/JMAX !!WRITE(1,*) I, XYM(I)END DOCLOSE(1)END===============================================! Monte Carlo Simulation of One Dimensional Diffusion INTEGER X,XY(1:1000,1:1000),y,XN(1:4),YN(1:4),RNREAL XYM(1:1000)! X:INSTANTANEOUS POSITION OF ATOM! XY(J,I):X*Y ,J:第几天实验,I:第几步跳跃! XYM(I): THE MEAN OF XYWRITE(*,*) "实验天数JMAX,实验次数IMAX"READ(*,*) JMAX,IMAXXN=(/0,0,-1,1/)YN=(/-1,1,0,0/)ISEED=RTC()DO J=1,JMAX !第几天实验X=0 !!!Y=0 !!!DO I=1,IMAX !第几步跳跃RN=4*RAN(ISEED)+1X=X+YN(RN)Y=Y+YN(RN)XY(J,I)=X*X+Y*YEND DOEND DOOPEN(1,FILE="C:\DIF2.DAT")DO I=1,IMAXXYM=0.0XYM(I)=1.0*SUM(XY(1:JMAX,I))/JMAX !!WRITE(1,*) I, XYM(I)END DOCLOSE(1)END做三维空间随机行走??留作业! Monte Carlo Simulation of One Dimensional Diffusion INTEGER X,XY(1:1000,1:1000),y,XN(1:6),YN(1:6),ZN(1:6),RN REAL XYM(1:1000)! X:INSTANTANEOUS POSITION OF ATOM! XY(J,I):X*Y ,J:第几天实验,I:第几步跳跃! XYM(I): THE MEAN OF XYWRITE(*,*) "实验天数JMAX,实验次数IMAX"READ(*,*) JMAX,IMAXXN=(/0,0,-1,1,0,0/)YN=(/-1,1,0,0,0,0/)ZN=(/0,0,0,0,1,-1/)ISEED=RTC()DO J=1,JMAX !第几天实验X=0 !!!Y=0 !!!Z=0DO I=1,IMAX !第几步跳跃RN=6*RAN(ISEED)+1X=X+XN(RN)Y=Y+YN(RN)Z=Z+ZN(RN)XY(J,I)=X*X+Y*Y+Z*ZEND DOEND DOOPEN(1,FILE="C:\DIF2.DAT")DO I=1,IMAXXYM=0.0XYM(I)=1.0*SUM(XY(1:JMAX,I))/JMAX !!WRITE(1,*) I, XYM(I)END DOCLOSE(1)==============================================通过该程序了解fortran语言如何画图(通过像素画图)USE MSFLIBINTEGER XR,YR !在的区域中画一个圆PARAMETER XR=400,YR=400INTEGER R, S(1:XR,1:YR)X0=XR/2 ! 圆心位置X0,YOY0=YR/2R=MIN(X0-10,Y0-10) !圆半径S=0 !像素的初始状态(颜色)DO I=1,XRDO J=1,YRIF((I-X0)**2+(J-Y0)**2<=R**2)S(I,J)=10IER=SETCOLOR(S(I,J))IER=SETPIXEL(I,J)END DOEND DOEND=================================================! 画一个圆(1、如何选出晶界区域;2、进一步加深对画图的理解)USE MSFLIBINTEGER XR,YR !在的区域中画一个圆PARAMETER XR=400,YR=400INTEGER R, S(0:XR+1,0:YR+1), XN(1:4), YN(1:4), SNSXN=(/0,0,-1,1/)YN=(/-1,1,0,0/)X0=XR/2 ! 圆心位置X0,Y0Y0=YR/2R=MIN(X0-10,Y0-10) !圆半径S=0 !像素的初始状态(颜色)DO I=1,XRDO J=1,YRIF((I-X0)**2+(J-Y0)**2<=R**2)S(I,J)=10IER=SETCOLOR(S(I,J))IER=SETPIXEL(I,J)END DOEND DODO I=1,XR !画晶界DO J=1,YRNDS=0DO K=1,4IF(S(I,J).NE.S(I+XN(K),J+YN(K)))NDS=NDS+1END DOIF(NDS>0)THENIER=SETCOLOR(9)ELSEIER=SETCOLOR(8)END IFIER=SETPIXEL(I,J)END DOEND DOEND如何画有一定宽度的晶界?====================================================!MC模拟一个晶粒的缩小(1、如何定义基体和晶粒;(2、如何寻找边界;(3、如何计算能量(构造或描述问题的概率过程)步骤:1、在基体上画一个原始晶粒,并赋予基体和原始晶粒不同的状态(取向号或体积能)2、寻找晶界3、MC的一个时间步(晶粒长大过程中一种随机性)4、计算晶界上网格点的相互作用,每个格点的相互作用导致晶粒长大①随机选取一个初始格点;②若此点属于晶界,那么可以随机转变为它邻居的状态,若不是晶界,则跳出循环,不发生转变;③计算转变前后的能量变化,⊿E=⊿E v +⊿E s +⊿E q(自由能=体积能+表面能+能量起伏)④若⊿E小于0,则新晶相被接受,网格状态发生转变。

USE MSFLIBPARAMETER IR=400,JR=400INTEGER IS(0:IR+1,0:JR+1),TMAX,ISN(1:8),NSTA TE,T,NR,IX,IYWRITE(*,*)"PLEASE INPUT THE TIME STEP "READ(*,*)TMAXISEED=RTC()! 定义圆心和半径JRC=JR/2R=MIN(IRC,JRC)-10! 定义基体和圆晶粒分别为状态1、状态2IS=1DO I=1,IRDO J=1,JRDISTANCE=SQRT(1.0*(I-IRC)**2+1.0*(J-JRC)**2)IF(DISTANCE.LT.R)IS(I,J)=2ISE=SETCOLOR(IS(I,J))ISE=SETPIXEL(I,J)END DOEND DOOPEN(1,FILE="E:\LUKE.DAT")! 寻找晶粒边界,计算能量,改变状态。

相关主题