
\section{Statistical Mechanics}



\subsection{ Statistical Physics} 

Some physical processes, especially quantum processes such as nuclear 
reactions, are inherently probabilistic, and can only be described 
statistically.

Many large classical systems (e.g. a volume of gas) have so many
variables, or {\em degrees of freedom} 
(positions and velocities of every molecule), that
an exact treatment is completely intractable, and in any case, not
really useful.  A statistical approach works very well, since the 
number of degrees of freedom is extremely large.

{\em Statistical mechanics} can be used to model both these types of systems.






\subsection{ Stochastic (probabilistic) Processes}

{\bf Example: Flux of neutrons in a nuclear reactor}

Reactor is simulated by tracking a neutron with random initial
conditions (position, momentum) through the system.  
There are different {\em probabilities} 
for absorbtion, scattering, triggering nuclear reactions to produce more
neutrons, etc., for the different materials in the reactor. 
A statistical ensemble of such neutrons gives the required flux.

% Figure ??






\subsection{ Systems with Many Degrees of Freedom}

{\bf Example: Thermodynamic properties of gases}

It is totally impractical and not useful to know the exact
microstate (positions and momenta of all molecules) of a large volume
of gas.

Useful properties are {\bf statistical}:
average energy of particles (temperature), average momentum change
from collisions with walls of container (pressure), etc.

The law of large numbers implies the error in the averages decreases 
as the number of particles increases.  
Macroscopic volume (state) of gas has O($10^{23}$) molecules. 
Thus a statistical approach works very well!

A given state or {\em configuration} of the system occurs with a 
particular {\em probability}, subject to a {\em probability distribution}.






\subsection{ Applications of Statistical Mechanics}


Statistical Mechanics is applicable to:

\begin{itemize}

\item Radioactive processes, nuclear reactions

\item Quantum processes in semiconductors, superconductors, etc.

\item High energy physics, reactions in particle accelerator collisions

\item Thermal properties of matter (solids, liquids, gases), 
      phase transitions, heat conduction

\item Magnetic and electronic properties of matter

\item Metallic alloys

\item Diffusion and percolation

\end{itemize}






\subsection{ Applications of Statistical Mechanics (cont.)}

Statistical mechanics techniques are also applicable to:

\begin{itemize}

\item Quantum field theories of electromagnetism (electrons, photons)
		and strong nuclear force (quarks, gluons)

\item String theories and quantum gravity

\item Quantum chemistry (structure and interactions of molecules)

\end{itemize}







\subsection{ Calculations in Statistical Mechanics}

{\bf Analytic methods:} almost always approximations.  Very few exact
solutions.  Good for understanding basic physics, but often break down
in regions of interest such as phase transitions.  Little knowledge or
control of errors, not easily systematically improvable.

\bigskip

{\bf Monte Carlo methods:} calculations usually involve sums or integrals
over very large number of dimensions or configurations.  Monte Carlo
is the obvious choice.  Can compute quantities for any parameter values
(including where analytic methods don't work).  Errors can be
estimated, and are systematically improvable by using more sample
configurations and larger systems (i.e. more computer power!).




\section{Spin Models}



\subsection{ Example: Spin Models of Magnetism}

This example will be used as the basis for most of the discussion of
simulation in statistic physics, for a number of reasons:

\begin{itemize}

\item Simple models of real physical systems

\item Relevant to a variety of other systems linked by the
		theory of critical phenomena and phase 
		transitions, e.g. thermodynamics and quantum
		field theory

\item Analytic solutions exist for some simple models, which is
      useful for checking numerical techniques

\item Classic example of the use of Monte Carlo methods

\item Different algorithms exist

\item Interesting application for parallel computation -- 
      a lot of work has been done on parallel computers 
      (and has helped drive their development) 

\end{itemize}






\subsection{ Magnetic Phase Transitions}

Many magnetic crystals have a phase transition from a magnetized
(ordered) state to a demagnetized (disordered) state.  This is not
just a simple change in magnetization -- it is characterized by a
divergence of certain quantities (for an infinite system) at the Curie
Temperature $T_c$.

% Figures ??






\subsection{ Simple Model of Magnetism}

Unpaired electron spins couple and align. 
The sum of their magnetic fields gives macroscopic magnetism.

\begin{tabbing}
indenting \= high temperature \=   \kill
\> low temperature   \>  $\Rightarrow \,$  order	\\
\>           \>  $\rightarrow \,$  alignment \\
\> 	     \>  $\rightarrow \,$  large magnetization  \\
\>           \>  \\
\> high temperature  \>  $\Rightarrow \,$  disorder (thermal fluctuations) \\
\>           \>  $\rightarrow \,$  random alignment \\
\> 	     \>  $\rightarrow \,$  spins and fields cancel \\
\> 	     \>  $\rightarrow \,$  no magnetization 
\end{tabbing}






\subsection{ The Ising Model}

Spins have only two states, $+1 \, (\uparrow)$ or $-1 \, (\downarrow)$.
Energy is given by
%
\begin{displaymath}
E = -J \sum_{\langle i,j \rangle} S_i S_j
\end{displaymath}
%
\hspace*{2cm} $S_i$ =  spin at site $i$ of crystalline lattice \\
\hspace*{2cm} $\langle i,j \rangle$ =  nearest neighbors in the lattice \\
\hspace*{2cm} $J$ =  interaction strength 

\begin{tabbing}
indentation \= bigspace \= spins not aligned- \=   \kill

\> $\uparrow \uparrow$, $\downarrow \downarrow$  \>  spins aligned 
                                                 \> $\Rightarrow \, E = -J$ \\
\> $\uparrow \downarrow$, $\downarrow \uparrow$  \>  spins not aligned 
                                                 \> $\Rightarrow \, E = +J$
\end{tabbing}

If $J > 0$, the state of lowest energy is when all spins are aligned.
This is called a ferromagnet.  

$J < 0$ is an antiferromagnet, which also occurs in real systems  
(will come back to this later).






\subsection{ Ising Ferromagnet}

For this model, statistical mechanics tells us

\begin{tabbing}
indent \= high temperature plus space \=   \kill

\>  low temperature ($T \rightarrow 0$)  \>  $\rightarrow \,$  low $E$	\\
\>                             \>  $\rightarrow \,$  spins aligned  \\
\>                             \>  $\rightarrow \,$  large magnetization  \\
\>                             \>   \\
\>  high temperature ($kT \gg J$)  \>  $\rightarrow \,$ large thermal energy \\
\>                      \>  $\rightarrow \,$  spins randomly aligned  \\
\>                      \>  $\rightarrow \,$  no magnetization
\end{tabbing}

Here magnetization $M = {\displaystyle \sum_i} \, S_i$.

So low $T$ and high $T$ behavior for the Ising model is the same 
as a real ferromagnet.






\subsection{ Solution of the Ising Model}

This model was proposed by Lenz in 1923 as a problem for his graduate
student, Ising, who showed that the 1-d model had no phase transition,
i.e. the magnetization varied continuously and slowly from $T=0$ to
$T=\infty$, with the susceptibility $\chi$ staying finite.

Onsager solved the 2-d problem exactly in 1944, showing that it had a
phase transition of the type seen in real ferromagnets.  This led to
great advances in the theory of phase transitions and critical
phenomena.

The 3-d problem also has a phase transition, but is not exactly soluble.
Results come from analytic approximations and Monte Carlo simulation,
which is the usual case for more complicated spin models.






\subsection{ Statistical Mechanics of Spin Models}

How does one go about solving the Ising model?

Want to work out thermodynamic averages of physically measurable
quantities such as the magnetization $M$ and the energy $E$ as a function
of the temperature $T$.  The thermal average of $M$ is
%
\begin{displaymath}
\langle M \rangle = \sum_{C} \, p(C) \, M(C)
\end{displaymath}
%
where $M(C)$ is the magnetization of the configuration $C$, $p(C)$ is the
probability distribution for the configuration as a function of $T$, and
the sum is over all configurations.

A fundamental result in statistical mechanics, which we will not
derive here, but accept as a physical law, is that $p$ is the Boltzmann
distribution.






\subsection{ The Boltzmann Distribution}

The Boltzmann probability distribution is defined by
\begin{eqnarray*}
p(C) & = & {1 \over Z} e^{-E(C)/kT} \\
   Z & = & \sum_C e^{-E(C)/kT}
\end{eqnarray*}

Here $E$ is the energy (or, more generally, the Hamiltonian) of the
configuration, $T$ is the temperature, and $k$ is Boltzmann's constant.
For convenience we often choose units so that $k = 1$, in which case $T$
is measured in energy units.

$Z$ is called the partition function, and is the fundamental quantity in
statistical mechanics.  All quantities of interest can be extracted
from $Z$, so an analytic formula for $Z$ implies an exact solution of the
model.  The sum will actually be an integral if the space of possible
configurations is continuous rather than discrete.





\section{Monte Carlo Methods}




\subsection{ Review of Monte Carlo Integration}

% Figure ?? 

$
\begin{array}{rl}
\hspace*{4in} I & \, = \, {\displaystyle \int_{a}^{b} } \, f(x) \, dx \\
\hspace*{4in}   & \, \approx \, {(b-a) \over N} 
{\displaystyle \sum_{i=1}^{N} } \, f(x_i) 
\end{array}
$

\bigskip

Instead of choosing $x_i$ at regular intervals, just choose them at
random.  Error is purely statistical, $\propto 1 / \sqrt{N}$ 
{\em independent of the dimension of the integral}.
Other integration methods have error $\propto 1 / N^{1/d}$.  

Example is the numerical estimate of pi, using simple 2-d Monte Carlo 
(dartboard method). 

$
\begin{array}{rl}
\hspace*{2cm} \pi & = {\displaystyle \int_{-1}^{1} \int_{-1}^{1} } 
\, p(x,y) \, dx \, dy \\
\hspace*{2cm}     & \approx {4 \over N} {\displaystyle \sum_{i=1}^{N} }
\, p(x_i,y_i) \\
\hspace*{2cm} p(x,y) & = 1 \qquad x^2 + y^2 \leq 1 \\
\hspace*{2cm}        & = 0 \qquad \mbox{otherwise}
\end{array}
$

% Figure ?? 






\subsection{ Monte Carlo Simulation}

To calculate sums (or integrals, if the spins are continuous-valued)
over such a large number of degrees of freedom, we resort to Monte
Carlo techniques.

As with the Monte Carlo integration mentioned earlier, we could just
generate configurations at random, and approximate the real thermal
averages, such as $\langle M \rangle$, by Monte Carlo averages:
%
\begin{displaymath}
\overline{M} = 
 {1 \over N} \, { { {\displaystyle \sum_{i=1}^N} \, M(C_i) \, e^{-E(C_i)/kT} } 
\over { {\displaystyle \sum_{i=1}^N} e^{-E(C_i)/kT} } }
\end{displaymath}

The problem here is that because of the rapidly varying exponential
function in the Boltzmann distribution, most randomly chosen
configurations will make a negligible contribution to the sum,
since $E$ will be relatively large. 






\subsection{ Calculating the Partition Function}

Suppose we want to calculate the exact partition function $Z$ 
numerically.  Need to do this for all $T$, but let's start with just
one temperature.  For a real system, want O($10^{23}$) spins, but let's
start small and try to solve for a $32 \times 32$ lattice, 
i.e. O($10^3$) spins.

Number of configurations in the sum = $2^{32 \times 32} \sim 10^{300}$.

Suppose we had a gigantic parallel supercomputer, with 10 million
processors.  Each processor could generate a configuration $C$, calculate
$E(C)$ and the Boltzmann factor $e^{-E(C)/kT}$
and add it to the sum over configurations 
in one nanosecond (less than the time for a single instruction for the
fastest modern computer).  
We run this calculation for the age of the universe. This gives

\hspace*{1cm} $10^7$ procs $\times \, 10^9$ configs/proc/sec 
              $\times \, 10^{14}$ sec/yr $\times \, 10^{10}$ yrs \\
\hspace*{1cm} $\sim 10^{40}$ configs \qquad (not real close to what is needed)






\subsection{ Importance Sampling}

In order to get sensible, accurate results when simulating statistical
systems with a rapidly varying Boltzmann distribution, it is vital to
use the idea of importance sampling in Monte Carlo integration.

Clearly the ideal situation would be to sample configurations with a
probability given by their Boltzmann weight $p(C)$.  Then the Monte
Carlo average for $M$ would just be:
%
\begin{displaymath}
\overline{M} = {1 \over N} \sum_{i=1}^N \, M(C_i) 
\end{displaymath}

This is great!  Except that the sampling probability 
$p(C_i) = e^{-E(C_i)/kT} / Z$ 
depends on the partition function $Z$, which is basically what we are trying 
to calculate in the first place!  

If we don't know what $p(C)$ is, how can we do this? 






\subsection{ Markov Processes}

Let us set up a so-called {\em Markov chain} of configurations $C_t$ by the
introduction of a fictitious dynamics.  The ``time'' $t$ is computer time
(marking the number of iterations of the procedure), NOT real time --
our statistical system is considered to be in equilibrium, and thus
time invariant.

Let $P(A,t)$ be the probability of being in configuration $A$ at time $t$.

Let $W(A \! \rightarrow \! B)$ be the probability per unit time, or 
{\em transition probability}, of going from $A$ to $B$. Then:
%
\begin{displaymath}
P(A,t \! + \! 1) = P(A,t) + \sum_B \, 
[ \, W(B \! \rightarrow \! A) \, P(B,t) - W(A \! \rightarrow \! B) \, P(A,t) \, ] 
\end{displaymath}

At large $t$, once the arbitrary initial configuration is ``forgotten,''
want $P(A,t) \rightarrow p(A)$.






\subsection{ Detailed Balance}

Clearly a sufficient (but not necessary) condition for an equilibrium
(time independent) probability distribution is the so-called detailed
balance condition
%
\begin{displaymath}
W(A \! \rightarrow \! B) \, P(A,t) = W(B \! \rightarrow \! A) \, P(B,t) 
\end{displaymath}

This method can be used for any probability distribution, but if we
choose the Boltzmann distribution
%
\begin{eqnarray*}
{W(A \! \rightarrow \! B) \over W(B \! \rightarrow \! A) } 
          & = & {p(B) \over p(A)} = { e^{-E(B)/kT} \over e^{-E(A)/kT} } \\
          & = & e^{- \Delta E / kT} \\
\Delta E  & = & E(B) - E(A)
\end{eqnarray*}

N.B.  $Z$ does not appear in this expression! It only involves quantities 
that we know $(k, T)$ or can easily calculate $(E)$.






\subsection{ The Metropolis Algorithm}

This dynamic method of generating an arbitrary probability
distribution was invented by Metropolis, Teller$^2$, and Rosenbluth$^2$ 
in 1953 (supposedly at a Los Alamos dinner party).

There are many possible choices of the $W$'s which will satisfy detailed
balance. They chose a very simple one:
  
$
\begin{array}{clc}
\hspace*{4.1cm} W(A \! \rightarrow \! B) 
& = \, e^{- \Delta E / kT} & \mbox{if} \; \Delta E > 0 \\
\hspace*{4.1cm}
& = \, 1                   & \mbox{if} \; \Delta E \leq 0 
\end{array}
$

So, if $E(B) > E(A)$
%
\begin{displaymath}
{W(A \! \rightarrow \! B) \over W(B \! \rightarrow \! A) } 
          \, = \, { e^{-[E(B)-E(A)]/kT} \over 1} 
          \, = \, e^{- \Delta E / kT} 
\end{displaymath}

and if $E(B) \leq E(A)$
%
\begin{displaymath}
{W(A \! \rightarrow \! B) \over W(B \! \rightarrow \! A) } 
          \, = \, {1 \over e^{-[E(A)-E(B)]/kT} } 
          \, = \, e^{- \Delta E / kT} 
\end{displaymath}






\subsection{ Monte Carlo Algorithms}

So we have a valid Monte Carlo algorithm if:

\begin{itemize}

\item  We have a means of generating a new configuration $B$ 
		from a previous configuration $A$ such that 
		the transition probability $W(A \! \rightarrow \! B)$ 
		satisfies detailed balance

\item  The generation procedure is {\it ergodic}, i.e. every configuration
		can be reached from every other configuration
		in a finite number of iterations

\end{itemize}

The Metropolis algorithm satisfies the first criterion for all
statistical systems.  The second criterion is model dependent, and not
always true (e.g. at $T=0$).






\subsection{ Monte Carlo References}



Here are some references to books with introductions to
Monte Carlo methods, especially for statistical physics:
%
\begin{itemize}

\item
S. Koonin, 
{\em Computational Physics}

\item
H. Gould and J. Tobochnik,
{\em An Introduction to Computer Simulation Methods, Vol. 2}

\item
O. Mouritsen,
{\em Computer Studies of Phase Transitions and Critical Phenomena}

\item
M.H. Kalos and P.A. Whitlock,
{\em Monte Carlo Methods, Vol. I. Basics}

\end{itemize}


These are a bit more advanced:
%
\begin{itemize}

\item
K. Binder ed., 
{\em Monte Carlo Methods in Statistical Physics}

\item
K. Binder ed., 
{\em Applications of the Monte Carlo Method in Statistical Physics} 

\item
D.P. Landau et. al, eds.
{\em Computer Simulation Studies in Condensed Matter Physics: 
Recent Developments}

\item
K. Kremer and K. Binder, {\it Comput. Phys. Reports} {\bf 7}, 259 (1988).

\end{itemize}





% \end{document}


