!============================================================================ ! Name : IsingModelF.f90 ! Author : Georg Bergner ! Version : ! Copyright : ! Description : Simple Fortran program to run Ising model simulations. ! Warning: several features (like random number generator) ! not optimal for precise investigations. !============================================================================ ! Functions to find the neighbouring points in up and down direction ! using periodic boundary conditions. integer function shiftup(num,NS) implicit none integer num integer NS shiftup=num+1 if(num==NS) then shiftup=1 endif end function shiftup integer function shiftdown(num,NS) implicit none integer num integer NS shiftdown=num-1 if(num==1) then shiftdown=NS endif end function shiftdown program IsingModelF implicit none ! The grid size and the parameter for the scan are fixed at compile time. integer, parameter :: XS=40 integer, parameter :: YS=40 ! Minimum and maximum value for the scan. double precision, parameter :: Jmin=0.2 double precision, parameter :: Jmax=0.8 integer, parameter :: nsteps=20 ! The thermalization has to be checked / adapted. integer, parameter :: nthermal=200 integer, parameter :: numconfigs=1000 ! Fulloutput provides the information about the thermalization. integer, parameter :: fulloutput=0 integer configuration(1:XS,1:YS) integer x,y,n,step double precision J,tmpsum, average, average2,averageabs,r integer shiftup, shiftdown double precision correlator(1:XS) double precision correlator2(1:XS) double precision correlatoroutput(2*nsteps+1,1:XS) print *,"#Parameter Nconfigs=",numconfigs," nthermal=",nthermal," XS=",XS," YS=",YS,"\n" ! Scan of parameter range of J. do step=1,nsteps J=Jmin+(step-1)*(Jmax-Jmin)/dble(nsteps-1) ! Initial configuration. do x=1,XS do y=1,YS configuration(x,y)=1 end do end do average=0 average2=0 averageabs=0 ! start the updater. do n=1,numconfigs ! measure magnetization tmpsum=0.0 do x=1,XS do y=1,YS tmpsum=tmpsum+configuration(x,y) end do end do tmpsum=tmpsum/dble(XS*YS) if(fulloutput.EQ.1) then print*, n,tmpsum end if if(n>=nthermal) then average=average+tmpsum average2=average2+tmpsum*tmpsum averageabs=averageabs+dabs(tmpsum) end if ! metropolis update ! This is simple typewriter ordering ! which is not the best for all applications. do x=1,XS do y=1,YS !change of local action tmpsum=configuration(shiftup(x,XS),y)+configuration(shiftdown(x,XS),y)+& configuration(x,shiftup(y,YS))+configuration(x,shiftdown(y,YS)) tmpsum=2.0*J*tmpsum*configuration(x,y) ! Simple random number generator. ! Not very precise results can be expected. call RANDOM_NUMBER(r) if (r < dexp(-tmpsum)) then configuration(x,y) = -1.0*configuration(x,y) end if end do end do end do ! Print out of average value for each J ! Limited number of observables, more can be added. ! Note that the error estimate is drastically underestimated ! since the autocorrelation time is not included. ! The absolute value of the order parameter provides a ! clearer identification of the phase transition. average=average/dble(numconfigs-nthermal+1) average2=average2/dble(numconfigs-nthermal+1) averageabs=averageabs/dble(numconfigs-nthermal+1) print *,J,average,dsqrt((average2-average*average)/dble(numconfigs-nthermal)),averageabs,& (average2-average*average)/dble(XS*YS) end do end program IsingModelF