Home Page > > Details

Ising model Assignment, phase transitions Assignment,Help With Python Course Assignment,Python Programming AssignmentDebug Haskell Programming|Debug

The Ising model and phase transitions
Remarks on completing the module
This assignment is summatively assessed. It is imperative that you submit the notebook on time.
Ising model
The task for this assignment is to implement the Ising Model introduced in the lecture. The structure in terms of
a code skeleton provided below needs to be followed. Otherwise the automatic tests, which allow you to test
different parts of you implementation, will not work.
We consider an Ising model, in which the interaction energy of spin is calculated from
where the sum is over the 4 nearest neighbours of .
We will restrict the calculation to a 2 dimensional grid throughout, hence 4 nearest neighbours, ,
. Notice that the expression for is different from the form considered in the lecture.
To simplify the calculations, we will use the dimensionless interaction energy
where
in the following. Here, is Boltzmann's constant, and is the temperature.
Given all spin states, we calculate the ensenble-averaged macroscopic magnetization , as
The brackets denote the ensemble average. The parameter has the dimensions of energy, and
consequently is dimensionless.
Follow the numbered steps in the following cells.
The cells below describe how to proceed, step-by-step. Begin by reading through all steps, comparing the
instructions to the Ising model discussed in the lecture. Complete the assignment using the cells below.
Several cells allow you to test your implementation step by step.
1. Set up the regular grid
Set up a 2D grid in the form of a python class, called Grid . The class should contain
the spin array
the value of J
2019/11/7 phase_transition
https://notebooks.dmaitre.phyip3.dur.ac.uk/miscada-sc/user/czjs88/notebooks/Phase_Transition_19/phase_transition.ipynb 2/10
The spin array should be a 2D array of dimension , with . We will address a particular spin with its 2D
Cartesian coordinates , where and are the indices in the 2D array. So, for
example, spin refers to the spin located at vertex of the grid.
Initialize the spins on the grid randomly, meaning each spin can be either up, , or down, , with
equal probability.
When performing calculations on the grid, we will assume periodic boundary conditions
no marks
2. Calculate the energy
Write a method, energy , as part of the class Grid , which calculates the interaction energy of a given spin,
, by summing over its four nearest neighbours. The function should take the grid array, , and the cell
indices and as parameters. It should return a python tuple containing two dimensionless energies
corresponding to the energy of the current spin state of cell , and the energy of the flipped
spin state .
This means that for a cell with spin state , the method should return and
vice versa.
Remember to account for periodic boudnary conditions on the grid.
You can test the implementation of this method using the test cells provded. What are the interaction energies
of cells (6,6), (15,0) and (31, 17) of the assignment grid given below (please include the answer to this question
in the PDF you hand in).
3. Calculate the probability of flipping a spin
The probability that the spin at vertex is flipped depends on the spin states of its neighbours and the
value of as explained in the lecture.
Write a method prob_flip which calculates the probability that spin is flipped, given the (dimensionless)
interaction energies for the current state and the flipped state .
The probability for a flip is given by
You can test the implementation of this method using the test cells provided. What are the probabilities for cells
(12, 12), (18,0) and (31, 12) of the assignment grid given below (please include the answer to this question in
the PDF you hand in)?
2 marks
4. Calculate the macroscopic magnetisation,
Write a method which calculates the current macroscopic magnetisation of a given grid, and add it to the
Grid class. The function should take the grid-array as a parameter and return the mean, macroscopic
magnetisation,
M
2019/11/7 phase_transition
https://notebooks.dmaitre.phyip3.dur.ac.uk/miscada-sc/user/czjs88/notebooks/Phase_Transition_19/phase_transition.ipynb 3/10
You can test the implementation of this method using the test cells provded. Calculate the magnetisation of the
assignment grid. State the answer to 3 significant digits on the PDF file you hand in.
5. Red-black sweep
Write a method to sweep over all spins in turn, first in (say) and then in in a red-black pattern. Red-black
means, first loop over all the red cells on a red-black chessboard pattern (looping over them first, then ).
Once all the red cells are done, then update all the black cells. Add this method to the Grid class.
For each spin in turn, flip the spin or not following the criterion discussed in the lecture. This means that the
spin in each cell in turn should be flipped with a probability (as discussed in step 3).
You can use the methods implemented in step 2 and 3.
6. Thermalisation and magnetisation measurement
Starting from a random configuration, the system needs to be evolved over a number of full red-black sweeps
in order to reach thermal equilibrium. This thermalization is part of the method you develop in this step.
Write a method that starts by sweeping the grid times to allow for the system to thermalize before you
carry out any measurements.
Next, the method should perform a further sweeps, while in addition computing and recording the
value of after every sweep. Use the method you developed in step 4 and the sweep implementaton of step
5.
and are input parameters of the method. The method should return a numpy array of length
, containing the magentisations measured after each measurement sweep.
Add this method to the Grid class.
7. Thermalisation
Plot the magnetisation over time for 1000 full mesh sweeps for and (include the thermalisation
period in the plot). Include this plot in the PDF you hand in. Save the plot to a file 'Thermalization.pdf'
4 marks
β = 0.1, 0.8 1.6
8. Ensemble-averaged magnetisation, , as a function of
Once the system has thermalized (a good choice is thermalisation sweeps), measure the timeaveraged
magnetisation over 1000 sweeps. From this, estimate the ensemble-averaged magnetisation, . Plot
as a function of for . Include this plot in the PDF you
hand in. Save the plot to a file 'Magnetisation.pdf'.
M¯ β
Ntherm = 400

M¯ 1/β β = 1.6, 1.3, 1.1, 1.0, 0.9, 0.8, 0.7, 0.6, 0.4, 0.1
2019/11/7 phase_transition
https://notebooks.dmaitre.phyip3.dur.ac.uk/miscada-sc/user/czjs88/notebooks/Phase_Transition_19/phase_transition.ipynb 4/10
p g p
4 marks
9. Critical temperature
The critical temperature, , in the mean field approximation, is given by
where h is the number of nearest-neighbours, as discussed in the lecture. Use this to calculate the critical value
which corresponds to the critical temperature and mark it in the plot produced in step 8 (If you could not
produce the plot in step 8 you may also mention the numerical value of on the solution you hand in).
2 marks
10. Mean field approximation
The ensemble-averaged value of the mean magnitization in the mean field approximation, is the solution of
as discussed in the lecture. This equation can not be solved analytically for , for given .
Rewrite the equation in terms of using the relation between and derived before, and solve the resulting
formula numerically for . You may want to use the root finding
methods implemented in previous exercises. Redo the plot of step 8, but on top of the numerical result, over
plot the mean field approximation. Add a legend to the plot which allows to distinguish the solution obtained in
step 8 from the mean field approximation.
Save the plot as 'MeanField.pdf'
4 marks
β = 1.6, 1.3, 1.1, 1.0, 0.9, 0.8, 0.7, 0.6, 0.4, 0.1
Solution
Use the cells below to complete the assignment. Use the test cells provided to make sure you are on the right
track.
In [1]:
%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
import sys
import os
sys.path.append(os.getcwd())
from scipy.interpolate import CubicSpline
import pickle
2019/11/7 phase_transition
https://notebooks.dmaitre.phyip3.dur.ac.uk/miscada-sc/user/czjs88/notebooks/Phase_Transition_19/phase_transition.ipynb 5/10
In [ ]:
You can use the cells below to test your implementation
class Grid:
def __init__(self, size, beta):
'''This function initialises the grid, i.e. it sets the
grid size, the value of beta and initialises the cells of the
grid with randomly chosen 'plus' (1) or 'minus' (-1) states.'''
# set self.size, self.beta, and self.cells
# self.cells is a 2D array, so that self.cells[i,j]=+1 or -1, the spin in gr
# YOUR CODE HERE
raise NotImplementedError()

def energy(self, i, j, beta, grid):
'''This function calculates the energies 'e_plus' and 'e_minus'
corresponding to the two possible states 'plus' and 'minus' for the spin at
of a given grid with a given value of beta.
returns: the two energy states 'e_current' and 'e_flip' as a tuple.
'e_current' is the energy of the spin (i,j) in its current spin state
'e_flip' is the energy of the spin (i,j) if you were to flip its spin
'''
# YOUR CODE HERE
raise NotImplementedError()
return e_current, e_flip

def prob_flip(self, e_current, e_flip):
'''This function calculates the probability of a spin flip
for a given spin, given the energies e_current and e_flip
of the current and the flipped state for the cell.
returns: the probability for the flip'''
# YOUR CODE HERE
raise NotImplementedError()
return probability_for_flip

def sweep(self):
'''This function carries out a single red-black sweep.
returns: nothing. For each spin in turn, it compute the probability for
flipping the spin, using the prob_flip function. Comparing the probablity,
it draws a random number to decide whether or not the flip the spin '''
# YOUR CODE HERE
raise NotImplementedError()

def magnetisation(self, grid):
'''This function calculates the mean magnetisation of all the spin in the gr
returns: the mean magnetisation M'''
# YOUR CODE HERE
raise NotImplementedError()
return M

def do_sweeps(self, n_therm, n_measure):
'''This function carries out n_therm thermalisation sweeps and n_measure mea
At the end of each measurement sweep the average magnetisation is computed a
returns: an array of length 'n_measure' containing the recorded magnetisatio
It uses the sweep function, and the magnitization function'''
# YOUR CODE HERE
raise NotImplementedError()
return magnetisation

2019/11/7 phase_transition
https://notebooks.dmaitre.phyip3.dur.ac.uk/miscada-sc/user/czjs88/notebooks/Phase_Transition_19/phase_transition.ipynb 6/10
Test and assignment grids
The cell below loads the test grids, test_grid_1 and test_grid_2 , and their corresponding values for
, test_beta_1 and test_beta_2 .
In addition, it loads assignement_grid and assignment_beta .
Use the first two grids to test your implementation. Use the third grid to answer the assignment questions.
β
In [3]:
2. Interaction energy calculation
The cell below allows you to test your interaction energy calculation method. If it does not return an error your
implementation might be correct.
In [ ]:
grid 1 loaded, size= 32 , beta= 1.6
grid 2 loaded, size= 32 , beta= 0.8
grid 3 loaded, size= 32 , beta= 0.1
filename = 'test_grid_1.pickle'
f = open(filename, 'rb')
(test_grid_1, test_beta_1) = pickle.load(f)
f.close()
filename = 'test_grid_2.pickle'
f = open(filename, 'rb')
(test_grid_2, test_beta_2) = pickle.load(f)
f.close()
filename = 'assignment_grid.pickle'
f = open(filename, 'rb')
(assignment_grid, assignment_beta) = pickle.load(f)
f.close()
print(" grid 1 loaded, size=", len(test_grid_1)," , beta= ", test_beta_1)
print(" grid 2 loaded, size=", len(test_grid_2)," , beta= ", test_beta_2)
print(" grid 3 loaded, size=", len(assignment_grid)," , beta= ", assignment_beta)
g1 = Grid(len(test_grid_1), test_beta_1)
g2 = Grid(len(test_grid_2), test_beta_2)
cells = [(6,6), (15,0), (31,17)]
energies_1 = [(0.0, 6.4), (0.0, 6.4), (0.0, 6.4)]
energies_2 = [(0.8, 2.4000000000000004), (1.6, 1.6), (2.4, 0.8)]
for c, cell in enumerate(cells):
i = cell[0]
j = cell[1]
e_1 = g1.energy(i, j, test_beta_1, test_grid_1)
e_2 = g2.energy(i, j, test_beta_2, test_grid_2)
assert(np.isclose(e_1, energies_1[c]).all())
assert(np.isclose(e_2, energies_2[c]).all())
print("Your implementation might be correct!")
2019/11/7 phase_transition
https://notebooks.dmaitre.phyip3.dur.ac.uk/miscada-sc/user/czjs88/notebooks/Phase_Transition_19/phase_transition.ipynb 7/10
3. Probability calculation
The cell below allows you to test your probability calculation method. If it does not return an error your
implementation might be correct.
In [ ]:
4. Magnetisation calculation
The cell below allows you to test your magnetisation method. If it does not return an error your implementation
might be correct.
In [ ]:
The following hidden cell uses the assignment grid to test
the calculation of energies 2 marks
the calculation of the probabilities 2 marks
the calculation of the magnetization 2 marks
5 Implement the red-black sweep
6 Implement the thermalization step
7. Thermalisation
Plot the magnetisation over time for 1000 full mesh sweeps for and (include the thermalisation
period in the plot). Include this plot in the PDF you hand in. Save the plot as file 'Thermalisation.pdf'
4 marks
β = 0.1, 0.8 1.6
g1 = Grid(len(test_grid_1), test_beta_1)
energies = [(0.1, 0.3), (0.2, 0.2), (0.3, 0.1), (1.5, 1.6), (0.1, 1.6), (0.8, 2.4)]
probabilities = [0.8187307530779818, 1, 1, 0.9048374180359595, 0.22313016014842982,
this_prob = []
for i, e in enumerate(energies):
this_prob.append(g1.prob_flip(e[0], e[1]))
assert(np.isclose(this_prob, probabilities).all())
print("Your implementation might be correct!")
g1 = Grid(len(test_grid_1), test_beta_1)
assert(np.isclose(g1.magnetisation(test_grid_1),0.193359375))
g2 = Grid(len(test_grid_2), test_beta_2)
assert(np.isclose(g2.magnetisation(test_grid_2),-0.3203125))
print("Your implementation might be correct!")
2019/11/7 phase_transition
https://notebooks.dmaitre.phyip3.dur.ac.uk/miscada-sc/user/czjs88/notebooks/Phase_Transition_19/phase_transition.ipynb 8/10
In [ ]:
8. Ensemble-averaged magnetisation, , as a function of
Once the system has thermalized (a good choice is thermalisation sweeps), measure the timeaveraged
magnetisation over 1000 sweeps. From this, estimate the ensemble-averaged magnetisation, . Plot
as a function of for . Include this plot in the PDF you
hand in. Save the plot to a file 'Magnetisation.pdf'.
Perform the following steps:
for each of the listed values of beta:
create a random grid
sweep the for N=400 initial 'thermalization' steps
then sweep for another 1000 steps, calculating after every sweep
use this to compute the ensemble-averaged magnitization, for that value of
plot as a function of
4 marks
# for each value of beta=[0.1, 0.8, 1.6]
# generate a grid of size=32, with the given value of beta
# perform N=1000 red-black sweeps (each sweep runs over the full 32x32 grid)
# calculate the mean magnetization, M for each sweep
# Plot M as a function of sweep number.
# You may want to use some of the plotting commands below.
# set-up the figure
fig, ax = plt.subplots(1,1, figsize = (8, 5))
file = "Thermalisation.pdf"
# caculate mag, the average magnetization, for N=1000 sweeps
# pass the value of beta into the plot command to generate the label, as in
# ax.plot(np.arange(len(mag[b])), mag[b], label='beta=%.2f'%beta)
# YOUR CODE HERE
raise NotImplementedError()
print("Calculation finished")
# plot the result, annotate the file, and save the file
ax.set_xlabel(r'$N_{steps}$')
ax.set_ylabel(r'$M$')
ax.set_ylim([-1.05, 1.05])
ax.legend()
plt.savefig(file)
fig.show()
2019/11/7 phase_transition
https://notebooks.dmaitre.phyip3.dur.ac.uk/miscada-sc/user/czjs88/notebooks/Phase_Transition_19/phase_transition.ipynb 9/10
In [ ]:
9. Critical temperature
The critical temperature, , in the mean field approximation, is given by
where h is the number of nearest-neighbours, as discussed in the lecture. Use this to calculate the critical value
which corresponds to the critical temperature and mark it in the plot produced in step 8 (If you could not
produce the plot in step 8 you may also mention the numerical value of on the solution you hand in).
Enter your answer in the box below. If you don't see a box, execute the hidden cell below.
2 marks
10. Mean field approximation
The ensemble-averaged value of the mean magnitization in the mean field approximation, is the solution of
as discussed in the lecture. This equation can not be solved analytically for , for given .
Enter your analytically calculated value of beta_c to 3 sig figs
beta_c= 0
# Step 8: Magnetisation
# set-up the figure
file = "Magnetisation.pdf"
fig, ax = plt.subplots(1,1, figsize = (8, 5))
# the range of values of beta
betas = [1.6,1.3,1.1,1.0,0.9,0.8,0.7,0.6,0.4,0.1]
# Loop over values of beta, computing the ensemble-averaged M for each beta
# name the resulting ensemble-averaged M mean_mags
# It is an array with the same dimension as betas
# YOUR CODE HERE
raise NotImplementedError()
# make the plot
ax.set_xlabel(r'$\beta^{-1}$')
ax.set_ylabel(r'$\bar M$')
ax.set_ylim([-1.05, 1.05])
ax.set_xlim(0,4)
plt.plot(1./np.array(betas), mean_mags)
plt.savefig(file)
fig.show()
2019/11/7 phase_transition
https://notebooks.dmaitre.phyip3.dur.ac.uk/miscada-sc/user/czjs88/notebooks/Phase_Transition_19/phase_transition.ipynb 10/10
Rewrite the equation in terms of using the relation between and derived before, and solve the resulting
formula numerically for . You may want to use the root finding
methods implemented in previous exercises. Redo the plot of step 8, but on top of the numerical result, over
plot the mean field approximation. Add a legend to the plot which allows to distinguish the solution obtained in
step 8 from the mean field approximation.
The numerical value of versus shows that the transition from at hight to at low is not
infinitely sharp. To quantify where the transition occurs, it might be useful to compute , the value of where
. Calculate this for both the numerical and mean field approximation, and indicate the point
on the plot.
Save the plot as 'MeanField.pdf'
4 marks
# Implement the mean field calculation here: calculate mean_mag_MF, the mean field a
# for a given value of beta. Also implement the calculation of the critical value, b
# MFA
# YOUR CODE HERE
raise NotImplementedError()
# YOUR CODE HERE
raise NotImplementedError()

Contact Us - Email:99515681@qq.com    WeChat:codinghelp
Programming Assignment Help!