Ising Model
Shaun Lahert 08668680
January 15, 2015
Figure 1: Flowchart of the Program - Task 1
Task 2
The storage for the array is in the main program Ising.py and is simply a
numpy array of ints of 1 and 1. I investigated using booleans as storage which
I found to be much more efficient using the .nbytes command but I could not
figure out a way to preform multiplication correctly to simulate taking the inner
product.
Task 3
Initially I used an or statement for the periodic boundary conditions but I got
strange behaviour in the magnetisation so switched to using the mod operator
which thankfully worked
Task 4
This function is contained in metropolis.py
Task 5 & 6
No real comments on these tasks
Task 7
Under the timeskip condition I built an array of magnetization from each time
step array and used numpy.sum to average it. I also averaged over the size of
the array to get a value between 1 and 1.
Task 8
Using the sys module with the command sys.argv[x] Ising.py was changed to
accept commands from the shell.
Task 9
This part of the program is contained in array save pbm.py. Using np.where()
which I discovered by research into vectorisation the 1s were changed into zeros then the whole array was converted in boolean to be saved as a pbm
Task 10
No real comment on this task.
Task 11
The shell script is contained in Control loop.sh. It passes the looped temperature T , the interaction constant J, the magnetic field h and the change
in temperature DdeltaT as its system arguments. After the loop is completed
it calls another python function MagTempPlot.py with a system arguement
Magnetisation Temp.csv, this is the .csv file of magnetisation and temperature
created by Ising.py.
Task 12
Ising.py was used to create plots of magnetisation vs time for each temperature to show an equilibrium is reached. MagTempPlot.py was called to plot
Magnetisation vs temperature this showed a phase transition at around T = 2.4
for J = 2. It was also found that increasing J would raise the critical temperature which is apparent from how the metropolis algorithm works i.e depending
on the ratio of TJ .
Images of the array state around this temperature show a transition from
a random distribution i.e zero magnetization to domain formation to nearly
maximum magnetization. All the images can be found in the images and graphs
folder.
Temperature: 3.0
Temperature: 2.7
Temperature: 2.4
Temperature: 2.9
Temperature: 2.6
Temperature: 2.3
Temperature: 2.1
4
Temperature: 2.8
Temperature: 2.85
Temperature: 2.2
Extra Investigation
Continous Cooling
I worked on two extra features of the program, the first was completely pointless
physically for algorithm but it helped with my bash scripting skills. I was
having issues with no domains forming. After playing with some premade Ising
simulators I convinced myself that the system needed to be continuously cooled
or heated i.e no reinitialisation of the array. To achieve this I used my shell
script to pass the final array state of the previous temperature into python to
be used as the new initial array state.
Given the fact the we skip the first 1000 iterations I realise now this is completely
pointless. But if we didnt and we had a small temperature increment I think
it would be a good model for continuously heating or cooling a system.
Optimization
My second extra feature was figuring out how to get the model to run faster as
in its basic nested conditional loops formulation its terribly slow. I read up on
vectorisation of code and using numpy I got rid of all the loops except for the
timestep loop.
The result of this is in the Vectorisation folder. This was pretty much an
unmitigated disaster. There was two main sections to be optimized, the first
was the calculation of energy difference. This was achieved in a sense by using
numpy.roll() to shift the matrix up, down, left and right. These were added up
and then multiplied by J and the original matrix elementwise to give an array
of energy change doing away with the for loops that iterated over the matrix
and admittedly it was about 20 times faster.
This array of energy change was then passed to a vectorised metropolis algorithm
using np.logical or() which returned an array of True or False values based on
the condition. Using np.where() the matrix elements would flip based on this
True/False array.
Neither of these optimizations worked. I tried out both separately with the
standard code and they both gave the behaviour of magnetization going to 0
exponentially.
