4
$\begingroup$

I am working on 2D Ising model, when I plot "specific heat vs temperature", I can't see any discontinuity at critical temperature around Tc ~ 2.7K. I am enclosing results of all other thermodynamic quantities, which are correct except the specific heat.

Edit 1:

As asked by @Daniel Shapero, I ran the simulation for 64 x 64 lattice, results are posted in comments. As you can see the absence discontinuity in both specific heat and susceptibility plots.

Output:

mag eng sus sp

code:

#!/usr/bin/env python
# coding: utf-8

#--------------------------------------------------
## required packages ##
#--------------------------------------------------

import numpy as np
from numpy.random import random
import matplotlib.pyplot as plt

#--------------------------------------------------
## definitions ##
#--------------------------------------------------

def spin_field(rows, cols):
    ''' generates a configuration with spins -1 and +1'''
    return np.random.choice([-1, 1], size = (rows, cols))

def neighbours(x, y, lattice, dim):
    ''' finds the neigbours of a particular lattice point with periodic boundery conditons '''
    top   = (x - 1, y)
    bottom  = ((x + 1) % dim, y)
    left    = (x, y - 1)
    right = (x, (y + 1) % dim)

    return [lattice[left[0], left[1]],
            lattice[right[0], right[1]],
            lattice[top[0], top[1]],
            lattice[bottom[0], bottom[1]]
           ]

def energy_calc(x, y, lattice, dim, J, B):
    ''' calulates the energy of the whole configuration '''
    dE = - J * lattice[x, y] * (np.sum(neighbours(x, y, lattice, dim)) + B)
    return dE

def total_energy(lattice, dim, J, B):
    ''' Energy of the configuration  '''
    TE = 0
    for x in range(dim):
        for y in range(dim):
            TE += - J * lattice[x, y] * (np.sum(neighbours(x, y, lattice, dim)) + B)
    return TE

#--------------------------------------------------
## relavant details about the lattice ##
#--------------------------------------------------

# number of monte carlo simulations
mcs = 100

# sparse averaging
relax_sweeps = 50

# square lattice dimensions
dim = 20

# coupling constant
J = 1.0

# external magnetic field
B = 0.0

# initialisation of lattice with random spins
lattice = spin_field(dim, dim)

n1 = 1.0 / (relax_sweeps * dim * dim)
n2 = 1.0 / (relax_sweeps * relax_sweeps * dim * dim)

#--------------------------------------------------
## main program ##
#--------------------------------------------------

# number of points between temperature range
m = 50

#dummy indices
index = 0

# this array contains thermodynamics quantities
quantities = np.zeros([m,5])

for temp in np.linspace(0.1,5.0,m):
    M1 = M2 = E0 = E1 = E2 = 0
    for sweep in range(0, mcs + relax_sweeps):
        for x in range(0,dim):
            for y in range(0,dim):
                dE = -2 * energy_calc(x, y, lattice, dim, J, B)
                if (dE <= 0):
                    lattice[x,y] *= -1
                elif (np.exp(-1 * dE / temp) >= random()):
                    lattice[x,y] *= -1
                else:
                    continue
        if (sweep > mcs):
            E0 = total_energy(lattice, dim, J, B)
            E1 += E0
            E2 += E0 * E0

    quantities[index,0] = temp
    quantities[index,1] = (n1*E2 - n2*E1*E1) / (temp * temp)
    index += 1

#--------------------------------------------------
## plotting of relevant thermodynamic quantities ##
#--------------------------------------------------

plt.figure(3)
plt.plot(quantities[:,0],quantities[:,1],'k.',label='{0} x {0}'.format(str(dim)))
plt.legend()
plt.xlabel('Temperature')
plt.ylabel('Specific Heat')
plt.title('Ising Model 2D')
plt.show()

```
$\endgroup$
2
  • 2
    $\begingroup$ Can you run the simulation again at larger sizes? 40 x 40, 80 x 80, etc. This could just be a finite size effect and the change will grow sharper as you go to larger system sizes. See also this answer from the physics SE. $\endgroup$ Mar 19 '20 at 16:47
  • $\begingroup$ @DanielShapero I simulated it for 64 x 64 lattice. But I still discontinuity is absent in both specific heat and susceptibility, please look at the results, link $\endgroup$
    – 147875
    Mar 19 '20 at 17:23
5
$\begingroup$

I sincerely thank @Daniel Shapero for directing me towards this answer.

Discontinuity in the specific heat or susceptibility curves to be visible significantly, you should take much more finer measurements for large number of sweeps, say, I ran the simulation for

  1. 1024 steps for system to reach equilibrium
  2. 1024 steps for sparse averaging/to measure specific heat finely
  3. Temperature ranging from 1.6 - 3.2, with 100 points in-between
  4. Lattice I considered is 20 x 20

Output:

$\endgroup$
2
  • $\begingroup$ If your question is answered you can accept your own answer. $\endgroup$
    – cfdlab
    Mar 20 '20 at 3:23
  • $\begingroup$ @cfdlab but I can only accept it after 2 days $\endgroup$
    – 147875
    Mar 20 '20 at 3:25

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.