\documentclass[10pt]{article}
\usepackage{a4wide}
\usepackage{listings}
\usepackage{graphicx}
\usepackage{times}
\usepackage{amsfonts}
\usepackage{amssymb}
\newcommand{\iint}[1]{\int\!\!\!\!\int\limits_{\!\!\!#1}}
\newcommand{\iiint}[1]{\int\!\!\!\!\int\!\!\!\!\int\limits_{\!\!\!\!\!\!\!\!\!#1}}
\newcommand{\text}[1]{\mathrm{#1}}

\begin{document}

\section*{Format for delivery of report and programs}
%
The format of the project is that of a printed file or hand-written report.
The programs should also be included with the report. Write only your candidate number on the first
page of the report and state clearly that this is your report for project 5 of FYS3150/FYS4150, fall 2012. 
There will be a box marked 'FYS3150/FYS4150' at the reception of the Department of Physics (room FV128). 

\section*{Project 5, First and second order phase transitions, deadline  December 10, 12pm}




For this project you can build upon the program you developed for project 4 for the Ising model and 
the Potts model described in the lecture notes, see chapter 13. 

The Potts model has been, in addition to the Ising model, widely used in studies
of phase transitions in statistical physics. The so-called two-dimensional 
$q$-state Potts
model has an energy given by
\[
  E=-J\sum_{<kl>}^{N}\delta_{s_l,s_k},
\]
where the spin $s_k$ at lattice position $k$ can take the values
$1,2,\dots,q$. The Kronecker delta function $\delta_{s_l,s_k}$ equals 
unity if the spins are equal and is zero otherwise.
$N$ is the total number of spins.
For $q=2$ the Potts model corresponds to the Ising model. To see that 
we can rewrite the last equation as
\[
  E=-\frac{J}{2}\sum_{<kl>}^{N}2(\delta_{s_l,s_k}-\frac{1}{2})-\sum_{<kl>}^{N}\frac{J}{2}.
\]
Now, $2(\delta_{s_l,s_k}-\frac{1}{2})$ is +1 when $s_l=s_k$ and $-1$ 
when they are different. This model is thus equivalent to the Ising model  
except a trivial difference in the energy minimum given by a an additional
constant and a factor $J\rightarrow J/2$.
One of the many applications of the Potts model is to helium absorbed on 
the surface of graphite.

The Potts model exhibits a second order phase transition for low values
of $q$ and a first order transition for larger values of $q$.
Using Eherenfest's definition of a phase transition, 
a second order phase transition has second derivatives
of the free energy that diverge or are discontinuous (the heat capacity and susceptibility
in our case) while a first order transition has first derivatives
like the mean energy that are discontinuous.

\begin{enumerate}
\item[a)] Write a program which simulates the $q=2$ Potts model 
for a two-dimensional lattice with $40\times 40$
spins and compute the average energy and specific
heat. Compute also the absolute value of the mean magnetization, the susceptibility and 
the correlation function
\begin{equation}
  G_{ij} = \langle s_is_j \rangle-\langle s_i\rangle\langle s_j \rangle.
\end{equation}
Establish an appropriate temperature range for where you see a sudden change
in the heat capacity and susceptibility. 
To get appropriate statistics you should
allow for at least $10^6$ Monte Carlo cycles. 
In setting up this code you can easily extend your code from project 4.

In the first part of this project we limit ourselves to the Metropolis algorithm.

Compare these results with those obtained with the two-dimensional Ising model.
The exact critical temperature for the Ising model is $T_C=2.269$.
For the Potts model with $q=2$ it is half this value.

\item [b)]  Here we study the time autocorrelation function, discussed in chapter 13.7.2 of the lecture notes.
We use again the square lattice with $L=40$ spins in the $x$ and $y$ directions.
Choose first a temperature below the assumed critical  temperature (in units of $kT/J$) and study the 
mean energy and magnetisation (absolute value) as functions of the number of Monte Carlo cycles.
Use both an ordered (all spins pointing in one direction or have the same value) and a random
spin orientation as starting configuration. 
How many Monte Carlo cycles do you need before you reach an equilibrium situation?
Can you extract an equilibration measure?
Repeat this analysis for a temperature above the critical . Try also different seeds for the random
number generators.  The evaluation of the autocorrelation function should be included
in all your further calculations.
\newline\newline

Near $T_C$ we can characterize the behavior of many physical quantities
by a power law behavior.
As an example the mean magnetization is given by
\begin{equation}
  \langle {\cal M}(T) \rangle \sim \left(T-T_C\right)^{\beta},
\end{equation}
where $\beta=1/8$ is a so-called critical exponent. A similar relation
applies to the heat capacity 
\begin{equation}
  C_V(T) \sim \left|T_C-T\right|^{\alpha},
\end{equation}
and the susceptibility
\begin{equation}
  \chi(T) \sim \left|T_C-T\right|^{\gamma},
\end{equation}
with $\alpha = 0$ and $\gamma = 7/4$.
Another important quantity is the correlation length, which is expected
to be of the order of the lattice spacing for $T>> T_C$. Because the spins
become more and more correlated as $T$ approaches $T_C$, the correlation
length increases as we get closer to the critical temperature. The divergent
behavior of $\xi$ near $T_C$ 
is
\begin{equation}
  \xi(T) \sim \left|T_C-T\right|^{-\nu}.
  \label{eq:xi}
\end{equation}
A second-order phase transition is characterized by a
correlation length which spans the whole system.
Since we are always limited to a finite lattice, $\xi$ will
be proportional with the size of the lattice. 
Through finite size scaling relations, see chapter 8.3 of Newman and Barkema,
it is possible to relate the behavior at finite lattices with the 
results for an infinitely large lattice.
The critical temperature scales then as
\begin{equation}
 T_C(L)-T_C(L=\infty) \sim aL^{-1/\nu},
 \label{eq:tc}
\end{equation}
with  $a$ a constant and  $\nu$ defined in Eq.~(\ref{eq:xi}).
The correlation length is given by
\begin{equation}
  \xi(T) \sim L\sim \left|T_C-T\right|^{-\nu}.
  \label{eq:xi2}
\end{equation}
and if we set $T=T_C$ one obtains
\begin{equation}
  \langle {\cal M}(T) \rangle \sim \left(T-T_C\right)^{\beta}
  \rightarrow L^{-\beta/\nu},
  \label{eq:scale1}
\end{equation}
a heat capacity
\begin{equation}
  C_V(T) \sim \left|T_C-T\right|^{-\gamma} \rightarrow L^{\alpha/\nu},
  \label{eq:scale2}
\end{equation}
and susceptibility
\begin{equation}
  \chi(T) \sim \left|T_C-T\right|^{-\alpha} \rightarrow L^{\gamma/\nu}.
  \label{eq:scale3}
\end{equation}

For all productions runs, you should parallelize the code at this stage. 
When doing this, be sure that you obtain the same
results as with the serial code.


\item [c)]  We will now study the behavior of the $q=2$ Potts  model in two dimensions close to the 
critical temperature as a function of the lattice size $L$.
Calculate the expectation values 
for $\langle E\rangle$ and $|\langle {\cal M}\rangle|$,
the specific heat $C_V$, the susceptibility $\chi$ and the correlation function
$G(r)$
as functions  of $T$ for $L=20$, $L=40$, $L=60$, $L=80$ and $L=100$ for $T$ near the critical temperature
$T_C\approx 1.1345$ for small temperature steps $\Delta T=0.01-0.05$. 
Plot  $\langle E\rangle$, $\langle {\cal M}\rangle$, $C_V$ and $G(r)$
as functions of $T$. Can you see an indication of a phase transition?
Try to extract the correlation length $\xi$ from the formal definition of the correlation function
\begin{equation}
  \xi^{-1} = -\lim_{r \rightarrow \infty} \frac{\partial}{\partial r} ln G(r),
\end{equation}
with $r=|i-j|$.

\item[d)]  Use Eq.~(\ref{eq:tc}) and the exact result
$\nu=1$ in order to estimate $T_C$ in the thermodynamic limit $L\rightarrow \infty$
using your simulations with $L=20$, $L=40$, $L=60$, $L=80$ and $L=100$.

\item[e)]  We use the exact result 
$kT_C/J=1/ln(1+\sqrt{2})\approx  1.1345$ and $\nu=1$.
Determine the numerical values of
$C_V$, $\chi$ and ${\cal M}$ at the exact value $T=T_C$
for  $L=20$, $L=40$, $L=60$, $L=80$ and $L=100$.
Plot $log_{10} ({\cal M})$
and $log_{10}(\chi$  $log_{10})$ as functions of  $log_{10}(L)$
and use the scaling relations of 
Eqs.~(\ref{eq:scale1}) and (\ref{eq:scale3}) in order to determine the constants
$\beta$ and $\gamma$. 
Are your log-log plots close to straight lines?
The exact values are $\beta=1/8$ and $\gamma=7/4$.

\item[f)]  Use the exact value of $T_C$ and compute the correlation time $\tau \sim L^z$ 
for the above lattices and make a log-log plot.
Comment your results for the Metropolis algorithm.



\item[g)] This exercise is optional and gives you an additional 30\% to the final score. 
Extend the calculations to the Potts model with 
$q=3$ and $q=6$. 
Establish first the location of the peak in the 
specific heat and study the behavior of the mean energy and magnetization 
as functions 
of $q$. Do you see a noteworthy change in behavior from the $q=2$ case?
Here you need only to find a region where you see that these variables
change dramatically. 
We study closer the behavior of the $q=3,6$ Potts  models in two dimensions close to the 
critical temperature as a function of the lattice size $L\times L$.
Calculate the expectation values 
for $\langle E\rangle$ and $|\langle {\cal M}\rangle|$,
the specific heat $C_V$ and the correlation function
$G(r)$
as functions  of $T$ for $L=20$, $L=40$, $L=60$, $L=80$ and $L=100$ for $T$ near the established 
critical temperature from the previous step.
Plot  $\langle E\rangle$, $\langle {\cal M}\rangle$, $C_V$ and $G(r)$
as functions of $T$. Can you see an indication of a phase transition? What is its character?
Comment your results.
Compare your critical temperature  with the values from the literature, see
for example the references listed below.

\end{enumerate}
For a reference on the critical temperatures of the Potts Models, see 
J.~L.~Monroe, Physical Review E {\bf 66}, 066129 (2002) and
\url{http://prola.aps.org/abstract/PRE/v66/i6/e066129}.
See also Challa and Landau at \url{http://link.aps.org/abstract/PRB/v34/p1841}.






\end{document}




\end{document}








