https://june-yeop.tistory.com/entry/Monte-Carlo-simulation-of-2-dim-Ising-model-using-Metropolis-algorithm  
 
  If you had studied quantum statistical physics, then you must be familiar to the spontatneous magnetisation of which most representative system is the Ising model, the most typical example of second order phase transition.

  Even though the way to solve the 2-dim Ising model exactly is well-known (Majorana fermionic operator + Wigner-Jordan transform, or C.N. Yang's paper), it demands solid backgrounds in quantum field theory or operational analysis to understand. Hence, we've taught the spontaneous magnetisation mostly via Landau's Mean Field Approach in undergraduate level statistical physics course.

  However, we can derive the tendancy of this phase transition in approximate level from a combination of numerical methods of which names are Monte-Carlo simulation and Metropolis algorithm.

  In this post, I'm gonna adumbrate basic concept of Monte-Carlo simulation and explain the codes that I wrote using Fortran90 and Python line by line.

0. Theoretical Set-Up

 
ARRRRRRRRR

1. Evolutions of spin-configuration and magnetisation w.r.t. iteration

Firstly, I wrote a Fortran code that shows how spin configuration and magnetisation of the system vary as the Monte-Carlo time-step. I used Fortran 90 for the main program and Python for visualisations.

1.1 Main Program (Fortran 90)

 
Let us begin with "PROGRAM ~" as usual.


    PROGRAM ISING_MODEL
  


Let us declare variables. For the place marked by question marks, just put any values as you want.

    !---DECLARE VARIABLES---------------------------

    IMPLICIT NONE


    INTEGER(4) :: I_ROW, I_COL 
    INTEGER(4) :: SPIN_LOW, SPIN_COL
    INTEGER(4) :: ITR
    iNTEGER(4) :: PRE_MAGNETISATION

    INTEGER, PARAMETER :: DIM_ROW = ???, DIM_COL = ???
    INTEGER, PARAMETER :: NUMBER_TOTAL_SPIN = DIM_ROW * DIM_COL
    INTEGER, PARAMETER :: ITERATION = ???

    INTEGER, DIMENSION(DIM_ROW , DIM_COL) :: SPIN_CONFIG


    DOUBLE PRECISION, DIMENSION(ITERATION + 1) :: NORMALISED_MAGNETISATION

    ! SPIN CORRELATION
    DOUBLE PRECISION, PARAMETER :: J = ???                  
    
    ! NUMBER OF NEAREST LATTICE SITES (assume square lattice)
    DOUBLE PRECISION, PARAMETER :: GAMMA = 4.               

    ! CRITICAL TEMPERATURE OF MEAN FIELD THEORY
    DOUBLE PRECISION, PARAMETER :: BETA_C = 1/(J*GAMMA)     

    ! BETA_RATIO = BETA / BETA_C ... THEN BETA = BETA_C * BETA_RATIO
    DOUBLE PRECISION, PARAMETER :: BETA_RATIO = ???         

    DOUBLE PRECISION :: DEL_ENERGY, E_I, E_F
    DOUBLE PRECISION :: RND, RND_ROW, RND_COL


    !--- DIRECTORY RELATED VARIABLES
    
    CHARACTER(100) FOLDER_NAME

    INTEGER :: FU
    CHARACTER(100) FILE_NAME
    CHARACTER(100) FILE_ID
    
    INTEGER :: IERR, UNIT_NUM, IOSTAT


    !--- VARIABLES FOR THE PROGRAM EXECUATION TIME
    REAL :: START_TIME, END_TIME, TOT_TIME
  


Let's call the cpu time to estimate the execuation time, cuz im curious kkk.

    !---LOAD CPU_TIME TO EVALUATE THE EXECUATION TIME---------------------------------
    
    CALL CPU_TIME(START_TIME)

    !---------------------------------------------------
  


Let us set up the initial spin configuration.

We are gonna consider only two cases for the initial configurations,
  i) high temperature limit (completely random configuration)
  ii) 0 temperature limit (completely ordered configuration)

     !---HIGH TEMPERATURE LIMIT ---------------------
    
      DO I_COL = 1, DIM_COL, +1

          DO I_ROW = 1, DIM_ROW, +1

              !---CHOOSE A REAL VALUE IN (0,1) RANGE
              CALL RANDOM_NUMBER(RND)

              !---COMPARE RND W/ 0.5
              !---IF RND < 0.5, GIVE UP-SPIN. ELSE, GIVE DOWN=SPIN
              IF ( RND < 0.50000000000000 ) THEN
                  SPIN_CONFIG(I_ROW, I_COL) = -1
              ELSE   
                  SPIN_CONFIG(I_ROW, I_COL) = +1
              END IF

          END DO
    
      END DO


      !---LOW TEMPERATURE LIMIT

      !DO I_COL = 1, DIM_COL, +1
      !    DO I_ROW = 1, DIM_ROW, +1
      !        SPIN_CONFIG(I_ROW, I_COL) = +1
      !    END DO
      !END DO
  


Let us calculate the normalised magnetisation of the initial spin configuration.

    !---CALCULATE THE MAGNETISATION OF INITIAL STATE-------------------------------
    
    ITR = 0

    PRE_MAGNETISATION = 0
    
    DO I_COL = 1, DIM_COL, +1
        
        DO I_ROW = 1, DIM_ROW, +1
        
            PRE_MAGNETISATION = PRE_MAGNETISATION + SPIN_CONFIG(I_ROW, I_COL)
        
        END DO
    
    END DO
    
    NORMALISED_MAGNETISATION(ITR) = 1. * PRE_MAGNETISATION / NUMBER_TOTAL_SPIN 
    !--- I PUT (1.*) TERM TO PREVENT FROME THE MAGNETISATION BECOMING INTEGER...!
  


Actually, the following procedure is unneccesary if we programmed using other advanced languages such as Python or Julia because these langueges support to make a list consists of a bunch of images. However, there's no package for Fortran that disdplay images as far as i know, we need the follwoing procedure.

I'm gonna save the spin configuration for each Monte-Carlo time step as a txt file while running a DO LOOP in a new folder in a folder that present code exists. The code below conducts these stuffs.

     !--- WRITE THE SPIN CONFIGURATION AS A TXT FILE, and MAGNETISATION DATA IN A NEW FOLDER

      ! Specify the folder name you want to create

        ! Write the integer into a string:
        FU = ITR + 20 !왜냐면 unit은 10부터 출발하자나...!!
    
        FOLDER_NAME = "spin_config"

      ! Use a system-specific command to create the folder
      ! Replace the following line with the appropriate command for your system (e.g., mkdir on Unix-like systems)
      ! On Windows, you can use "mkdir" as well.
        CALL SYSTEM("mkdir " // TRIM(FOLDER_NAME), IERR)

        IF (IERR == 0) THEN

        ! Specify the name of the text file to create inside the folder
            WRITE(FILE_ID, '(I0)') ITR
            FILE_NAME = TRIM(FOLDER_NAME) // "/config" // TRIM(ADJUSTL(FILE_ID)) // ".txt"
        
        ! Open the text file for writing
            OPEN( UNIT = FU , FILE = TRIM(FILE_NAME) , STATUS = "replace", ACTION = "write", IOSTAT = IOSTAT )

            IF (IOSTAT == 0) THEN
          ! Write data to the file
                WRITE(FU, *) SPIN_CONFIG
          ! Close the file
                CLOSE(FU)
            END IF

        END IF


        OPEN( 10, FILE = "magnetisation.txt", STATUS = "new" )    
        WRITE(10 ,*) ITR, NORMALISED_MAGNETISATION(ITR)
  

The last two lines are for making a new text file and openning it to jot down the magnetisation date.
The file is gonna be saved in the same directory of the present Fortran code.

1.2 Visualisation (Python)

 
ARRRRRRRRR

2. Magnetisation w.r.t. Temperature

1.1 Main Program (Fortran 90)

 
ARRRRRRRRR

1.2 Visualisation (Python)

 
ARRRRRRRRR

Leave a comment