\documentclass{article}
\usepackage{amsmath,amssymb,pictex,fullpage,cite}

\allowdisplaybreaks

\makeatletter
\let\over\@@over
\let\atop\@@atop
\let\above\@@above
\let\overwithdelims\@@overwithdelims
\let\atopwithdelims\@@atopwithdelims
\let\abovewithdelims\@@abovewithdelims
\let\frac\undefined
\let\genfrac\undefined
\makeatother

\def\smallskip{\vskip\smallskipamount}
\def\medskip{\vskip\medskipamount}
\def\bigskip{\vskip\bigskipamount}

\let\{=\lbrace
\let\}=\rbrace
\let\[=\lbracket
\let\]=\rbracket

\long\def\items#1{\begin{itemize}\parsep=0pt\itemsep=0pt#1\end{itemize}}
\long\def\enum#1{\begin{enumerate}\parsep=0pt\itemsep=0pt#1\end{enumerate}}

\def\eq#1{\begin{equation}#1\end{equation}}
\def\Align#1{\begin{align}#1\end{align}}
\def\Aligned#1{\begin{aligned}#1\end{aligned}}
\def\Gather#1{\begin{gather}#1\end{gather}}
\def\Gathered#1{\begin{gathered}#1\end{gathered}}
\def\Multline#1{\begin{multline}#1\end{multline}}
\def\Subeq#1{\begin{subequations}#1\end{subequations}}
\def\Matrix#1{\begin{matrix}#1\end{matrix}}
\def\pMatrix#1{\begin{pmatrix}#1\end{pmatrix}}
\def\Cases#1{\begin{cases}#1\end{cases}}

\def\d{\partial}
\def\bd{{\bar\d}}
\def\bpsi{{\bar\psi}}
\def\bvarphi{{\bar\varphi}}
\def\Re{\mathop{\rm Re}\nolimits}
\def\Im{\mathop{\rm Im}\nolimits}
\def\sh{\mathop{\rm sh}\nolimits}
\def\ch{\mathop{\rm ch}\nolimits}
\def\tg{\mathop{\rm tg}\nolimits}
\def\ctg{\mathop{\rm ctg}\nolimits}
\def\tr{\mathop{\rm tr}\nolimits}
\def\Tr{\mathop{\rm Tr}\nolimits}
\def\id{\mathop{\rm id}\nolimits}
\def\arctg{\mathop{\rm arctg}\nolimits}
\def\bn{{\boldsymbol n}}
\def\bp{{\boldsymbol p}}
\def\br{{\boldsymbol r}}
\def\bR{{\boldsymbol R}}
\def\bS{{\boldsymbol S}}
\def\bJ{{\boldsymbol J}}
\def\bsigma{{\boldsymbol\sigma}}
\def\Z{{\mathbb Z}}
\def\R{{\mathbb R}}
\def\C{{\mathbb C}}
\def\cH{{\cal H}}
\def\bz{{\bar z}}
\def\bw{{\bar w}}
\def\bw{{\bar w}}
\def\tphi{{\tilde\phi}}
\def\trho{{\tilde\rho}}
\def\ta{{\tilde a}}
\def\ve{\varepsilon}
\def\vq{{\vec q}}
\def\vx{{\vec x}}
\def\vy{{\vec y}}
\def\vJ{{\vec J}}
\def\va{{\vec a}}
\def\vb{{\vec b}}
\def\eff{{\rm eff}}
\def\const{{\rm \,const\,}}
\def\eF{{\epsilon_F}}
\def\Ech{{E_{\text{ch}}}}
\def\Esp{{E_{\text{sp}}}}
\def\Espel{{E_{\text{sp}}^{\text{el}}}}
\def\Eim{{E_{\text{im}}}}
\def\el{{\text{el}}}
\def\im{{\text{im}}}

\def\lcolon{\mathopen:}
\def\rcolon{\mathclose:}

\def\aux{{\tilde1}}

%%%%%%% For figures: %%%%%%%
\def\ar#1 #2, #3 #4 {\arrow <0.3truecm> [0.1,0.3] from #1 #2 to #3 #4 }
\catcode`!=11
\def\ln#1 #2, #3 #4 {\!start(#1,#2)\!ljoin(#3,#4)}
\catcode`!=12
\def\rl#1 #2, #3 #4 {\putrule from #1 #2 to #3 #4 }
\def\trl#1 #2, #3 #4 {\linethickness .8pt
	\putrule from #1 #2 to #3 #4 \linethickness .4pt}
\def\point#1 #2 {\put{\hbox{\kern -1pt .}} [Bl] at #1 #2 }
\def\bfpoint#1 #2 {%
\put{\hbox{\kern -2pt \raise -2.3pt%
\hbox{$\textstyle\bullet$}}} [Bl] at #1 #2 }
\def\figcap#1#2\par{\medskip{\narrower\eightpoint\noindent
	Fig.~#1. \ignorespaces#2\par}\medskip}
%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\makeatletter
\def\section{\@startsection{section}{1}{\z@}%
                                   {-3.5ex \@plus -1ex \@minus -.2ex}%
                                   {2.3ex \@plus.2ex}%
                                   {\normalfont\normalsize\bfseries}}
\def\subsection{\@startsection{subsection}{2}{\z@}%
                                     {-3.25ex\@plus -1ex \@minus -.2ex}%
                                     {1.5ex \@plus .2ex}%
                                     {\normalfont\normalsize\bfseries}}
\makeatother


\begin{document}

\centerline{\bf Lecture~1}
\centerline{\bf Eight-vertex model and commuting transfer matrices}
\bigskip

Since we shall speak about some models of physical importance, let us
formulate first the general physical framework. Consider a classical
system of interacting particles. Let $C$ be a state of this system,
which can run some (generally infinite) set of admissible states. Let
$E(C)$ be the energy of the state~$C$, defined by usual Hamiltonian
mechanics. Suppose that the system weakly interacts with a thermal bath
of the temperature~$T$. The most fundamental postulate of statistical
mechanics, the {\it Gibbs law}, says that the probability of the state
$C$ is given by
\eq{
w(C)={1\over Z}\,e^{-E(C)/T}.
\label{GibbsLaw}
}
It is easy to understand, what is the proportionality coefficient $1/Z$.
As the total probability for the system to be in any state is unity, we
have for the {\it partition functions\/} of the system
\eq{
Z=\sum_{C}e^{-E(C)/T}.
\label{PartitionFunction}
}
Since the space of configurations can be continuous, the sum may turn
out to be an integration, but we shall not consider this general case
in these lectures.

The partition function encodes the most fundamental observable
thermodynamic function of the model, the {\it free energy\/}:
\eq{
F=E-TS=-T\log Z.
\label{FreeEnergy}
}
Here $E$ is the {\it total energy\/} of the system,
$$
E=\sum E(C)w(C)={T^2\over Z(T)}{d\over dT}Z(T)
=-T^2{d\over dT}{F\over T}
=F-T{dF\over dT},
$$
while
$$
S=-\sum w(C)\log w(C)=-{dF\over dT}
$$
is the {\it entropy}. The first equality in~(\ref{FreeEnergy}) is
the thermodynamic definition of the free energy, which, in principle,
was established before Gibbs, while the second equality provides its
statistical interpretation.

Let $f_i(C)$ be some functions of the state of the system. The {\it
correlation functions\/} are generic expectation values
$$
\langle f_1\ldots f_n\rangle
=\sum_Cw(C)f_1(C)\ldots f_n(C).
$$
Note, that the correlation functions can be expressed in terms of
derivatives of a generalized partition function of the system with
external fields $F_i$ coupled to the variables $f_i$:
$$
Z(F_1,\ldots,F_N)=\sum_C e^{-E(C;F_1,\ldots,F_N)/T},
\qquad
E(C;F_1,\ldots,F_N)=E(C)-\sum_i F_if_i.
$$

Usually, from the physical point of view the most interesting objects
are local correlation functions, i.e.\ the correlation functions of the
amounts $f_i(C)$ that can be measured at some space points. We shall
discuss examples of such functions in detail later.

Naturally, evaluation of the partition function and correlation
functions is a difficult task, except some trivial examples solved by
classics of the science. Most these solved examples are systems of
independent particles, each of which possesses a finite set or, at
least, a finite-dimensional space of states and a simple function for
$E(C)$ (like a quadratic function in the Boltzmann gas). But how to do
in the case of interacting particles? The usual approach in physics is
developing some approximate methods, based either on the perturbation
theory or on some experimentally supported assumptions. It is a very
effective way, but sometimes physicists need some additional support to
their assumptions.

How mathematicians can help physicists? First, they can propose some
rigorous estimates, which can prove the effects predicted by
physicists. Many of such important estimates were proposed (see, e.~g.\
\cite{Ruelle,Sinai}). Second, they may propose some sophisticated
methods to solve exactly some particular nontrivial examples. Though
these methods can be not completely rigorous, they better convince
physicists, because they are more in the way of physical thinking, and
can be used by physicists themselves.

Let us now slightly specify the problem. Forget about motion of
particles. We can do it for some problems, for example, if positions of
particles are fixed by the crystalline lattice. The set of states of
the system ({\it configurations\/}) is reduced to the direct product of
sets of internal states of particles (`spins'). We shall assume these
`spins' to be discrete variables. Physically the discrete `spins' may
originate in quantum states of atoms (e.g.\ physical spin states) with
diagonal interaction or in other sources, e.g. in types of atoms in the
substitutional solutions or in positions of atoms in different quantum
wells as in the ice-type systems.

In these lectures we shall discuss two-dimensional models of statistical
mechanics. Why two-dimensional? Because in one dimension the lattice
statistical models are trivially reduced to quantum mechanics of small
systems without phase transitions or other interesting features. On the
other hand, in three dimensions the problem is too complicated and the
classes of solvable models are too narrow. We restrict our attention
onto two classes of two-dimensional solvable lattice models, which are
in a sense basic: eight-vertex model and solid-on-solid (SOS)
solvable model.

Let us start with the ice model on the square lattice called also the
six-vertex model~\cite{Lieb}. Consider a square lattice made of oxygen
ions ($\rm O^{2-}$) in the vertices and hydrogen ions ($\rm H^+$) on
the vertical and horizontal bonds of the lattice, which we shall call
edges. (Do not think of all this too seriously, because
experimentalists are unable to produce any two-dimensional ice-like
systems, but this picture of ice can help you both to remember the
formulation of the model and to better understand the physical origin
of such kind of problems.) We know that the hydrogen ions can form with
the oxygen ions two types of bonds: strong and narrow polar bond and
weak and long `hydrogen' bond. It means that a hydrogen ion lying on an
edge of the lattice must be positioned near one of the oxygen ion and
far from another oxygen ion. There are two such positions on each edge.
It forms discrete `spin' state on each edge. We shall think that the
`spin' is equal to `$+$' or `$+1$' if the hydrogen ion is positioned
near the right end of a horizontal edge or near the upper end of the
vertical edge, and is equal to `$-$' or `$-1$' if it is positioned near
the left or lower end of the edge.

Further, from the neutrality condition we conclude, that in the vicinity
of each oxygen ion there must be just two hydrogen ions. It imposes a
restriction onto values of `spins' at edges surrounding each vertex.
Namely, denote these spins by $\ve_1$, $\ve_2$, $\ve'_1$, $\ve'_2$:
$$
\vcenter{
\beginpicture
\setcoordinatesystem units <2cm,2cm> point at 0 -0.5
%
\rl 2.5 2 , 2.5 1
\rl 3 1.5 , 2 1.5
\put{$\ve_1$} [Bl] at 2.6 1.0
\put{$\ve_2$} [Bl] at 2.0 1.6
\put{$\ve'_1$} [Bl] at 2.6 1.9
\put{$\ve'_2$} [Bl] at 3.0 1.6
%
\endpicture
}
$$
Then this restriction, called the {\it ice rule}, can be expressed as
\eq{
\ve_1+\ve_2=\ve'_1+\ve'_2
\label{IceRule}
}
There are six configurations around each vertex that satisfy the ice rule:
$$
\underbrace{
\vcenter{
\beginpicture
\setcoordinatesystem units <1cm,1cm> point at 0 0
\rl -0.5 0.0 , 0.5 0.0
\rl 0.0 -0.5 , 0.0 0.5
\put{$+$} [tc] <0pt,-3pt> at 0.0 -0.5
\put{$+$} [Br] <-3pt,-2pt> at -0.5 0.0
\put{$+$} [Bc] <0pt,3pt> at 0.0 0.5
\put{$+$} [Bl] <3pt,-2pt> at 0.5 0.0
\endpicture}
\qquad
\vcenter{
\beginpicture
\setcoordinatesystem units <1cm,1cm> point at 0 0
\rl -0.5 0.0 , 0.5 0.0
\rl 0.0 -0.5 , 0.0 0.5
\put{$-$} [tc] <0pt,-3pt> at 0.0 -0.5
\put{$-$} [Br] <-3pt,-2pt> at -0.5 0.0
\put{$-$} [Bc] <0pt,3pt> at 0.0 0.5
\put{$-$} [Bl] <3pt,-2pt> at 0.5 0.0
\endpicture}
}_{\atop{\textstyle a}}
\qquad
\underbrace{
\vcenter{
\beginpicture
\setcoordinatesystem units <1cm,1cm> point at 0 0
\rl -0.5 0.0 , 0.5 0.0
\rl 0.0 -0.5 , 0.0 0.5
\put{$+$} [tc] <0pt,-3pt> at 0.0 -0.5
\put{$-$} [Br] <-3pt,-2pt> at -0.5 0.0
\put{$+$} [Bc] <0pt,3pt> at 0.0 0.5
\put{$-$} [Bl] <3pt,-2pt> at 0.5 0.0
\endpicture}
\qquad
\vcenter{
\beginpicture
\setcoordinatesystem units <1cm,1cm> point at 0 0
\rl -0.5 0.0 , 0.5 0.0
\rl 0.0 -0.5 , 0.0 0.5
\put{$-$} [tc] <0pt,-3pt> at 0.0 -0.5
\put{$+$} [Br] <-3pt,-2pt> at -0.5 0.0
\put{$-$} [Bc] <0pt,3pt> at 0.0 0.5
\put{$+$} [Bl] <3pt,-2pt> at 0.5 0.0
\endpicture}
}_{\atop{\textstyle b}}
\qquad
\underbrace{
\vcenter{
\beginpicture
\setcoordinatesystem units <1cm,1cm> point at 0 0
\rl -0.5 0.0 , 0.5 0.0
\rl 0.0 -0.5 , 0.0 0.5
\put{$+$} [tc] <0pt,-3pt> at 0.0 -0.5
\put{$-$} [Br] <-3pt,-2pt> at -0.5 0.0
\put{$-$} [Bc] <0pt,3pt> at 0.0 0.5
\put{$+$} [Bl] <3pt,-2pt> at 0.5 0.0
\endpicture}
\qquad
\vcenter{
\beginpicture
\setcoordinatesystem units <1cm,1cm> point at 0 0
\rl -0.5 0.0 , 0.5 0.0
\rl 0.0 -0.5 , 0.0 0.5
\put{$-$} [tc] <0pt,-3pt> at 0.0 -0.5
\put{$+$} [Br] <-3pt,-2pt> at -0.5 0.0
\put{$+$} [Bc] <0pt,3pt> at 0.0 0.5
\put{$-$} [Bl] <3pt,-2pt> at 0.5 0.0
\endpicture}
}_{\atop{\textstyle c}}
$$
Let us think that each of this vertex configuration may have its own
energy $E_1,\ldots,E_6$. The energy of the system is equal to the sum of
energies of vertices. It defines the model nearly completely. The only
thing to specify is the boundary conditions. We shall consider three
possibilities:

1.~The model of $M$ columns and $N$ rows with cyclic (toroidal) boundary
conditions.

2.~The model with fixed spins at the boundary.

3.~The model on an infinite lattice, considered as a limit of any of
these models as $M,N\to\infty$.

The six-vertex model is known to be solvable if
$$
E_1=E_2,
\qquad
E_3=E_4,
\quad
E_5=E_6.
$$
It means that the configurations braced together on the picture possess
the same energies and the same Boltzmann weights:
\eq{
a=e^{-E_1/T}=e^{-E_2/T},
\qquad
b=e^{-E_3/T}=e^{-E_4/T},
\qquad
c=e^{-E_5/T}=e^{-E_5/T}.
\label{abcdef}
}
This model is called the homogeneous six-vertex model without external
field.

The disadvantage of the six-vertex model is that it possesses some
pathologic physical properties, related to severeness of the ice
condition. Namely, consider its phase diagram:
$$
\Aligned{
{\rm AF}:
&\quad
c>a+b,
\\
{\rm F_1}:
&\quad
a>b+c,
\\
{\rm F_2}:
&\quad
b>a+c,
\\
{\rm D}:
&\quad
{\textstyle{1\over2}}(a+b+c)\ge a,b,c.
}
$$
The antiferroelectric (AF) region is the region of antiferroelectric
order considered as excitations above the following ground states:
\eq{
\vcenter{
\beginpicture
\setcoordinatesystem units <1cm,1cm> point at -4 0
\rl 1.0 3.3 , 1.0 0.7
\rl 2.0 3.3 , 2.0 0.7
\rl 3.0 3.3 , 3.0 0.7
\rl 0.7 3.0 , 3.3 3.0
\rl 0.7 2.0 , 3.3 2.0
\rl 0.7 1.0 , 3.3 1.0
\put{$+$} [lb] at 1.1 3.4
\put{$+$} [lb] at 0.4 3.1
\put{$+$} [lb] at 3.1 3.4
\put{$+$} [lb] at 2.4 3.1
\put{$+$} [lb] at 2.1 2.4
\put{$+$} [lb] at 1.4 2.1
\put{$+$} [lb] at 1.1 1.4
\put{$+$} [lb] at 0.4 1.1
\put{$+$} [lb] at 3.4 2.1
\put{$+$} [lb] at 3.1 1.4
\put{$+$} [lb] at 2.4 1.1
\put{$+$} [lb] at 2.1 0.4
\put{$-$} [lb] at 2.1 3.4
\put{$-$} [lb] at 1.4 3.1
\put{$-$} [lb] at 1.1 2.4
\put{$-$} [lb] at 0.4 2.1
\put{$-$} [lb] at 3.4 3.1
\put{$-$} [lb] at 3.1 2.4
\put{$-$} [lb] at 2.4 2.1
\put{$-$} [lb] at 2.1 1.4
\put{$-$} [lb] at 1.4 1.1
\put{$-$} [lb] at 1.1 0.4
\put{$-$} [lb] at 3.4 1.1
\put{$-$} [lb] at 3.1 0.4
%
\setcoordinatesystem units <1cm,1cm> point at -8 0
\rl 1.0 3.3 , 1.0 0.7
\rl 2.0 3.3 , 2.0 0.7
\rl 3.0 3.3 , 3.0 0.7
\rl 0.7 3.0 , 3.3 3.0
\rl 0.7 2.0 , 3.3 2.0
\rl 0.7 1.0 , 3.3 1.0
\put{$-$} [lb] at 1.1 3.4
\put{$-$} [lb] at 0.4 3.1
\put{$-$} [lb] at 3.1 3.4
\put{$-$} [lb] at 2.4 3.1
\put{$-$} [lb] at 2.1 2.4
\put{$-$} [lb] at 1.4 2.1
\put{$-$} [lb] at 1.1 1.4
\put{$-$} [lb] at 0.4 1.1
\put{$-$} [lb] at 3.4 2.1
\put{$-$} [lb] at 3.1 1.4
\put{$-$} [lb] at 2.4 1.1
\put{$-$} [lb] at 2.1 0.4
\put{$+$} [lb] at 2.1 3.4
\put{$+$} [lb] at 1.4 3.1
\put{$+$} [lb] at 1.1 2.4
\put{$+$} [lb] at 0.4 2.1
\put{$+$} [lb] at 3.4 3.1
\put{$+$} [lb] at 3.1 2.4
\put{$+$} [lb] at 2.4 2.1
\put{$+$} [lb] at 2.1 1.4
\put{$+$} [lb] at 1.4 1.1
\put{$+$} [lb] at 1.1 0.4
\put{$+$} [lb] at 3.4 1.1
\put{$+$} [lb] at 3.1 0.4
\endpicture
}
\label{GroundStates}
}
The excitations can be considered as some loops of changed spins and
can have finite energy even on the infinite lattice. The correlation
functions of spins decrease as $e^{-r/\xi}$ as the distance
$r\to\infty$. The constant $\xi$ is called {\it correlation length}. It
is a normal behavior of correlation functions out of special point
called critical points.

In the ferroelectric regions ($\rm F_1$ and $\rm F_2$) the situation is
strange: all excitations consist of the lines of flipped spins in the
SW--NE direction, which are infinite in the infinite volume limit! It
means that these excitations possess large energy, and their
contribution tends to zero as $M,N\to\infty$. We have the situation of
a {\it frozen order}. The free energy of such system per site tends to
zero.

Another pathological feature is related to the disordered (D) region.
The whole region turns out to be {\it critical}. It means that the
correlations between local variables like spins decrease with distance
$r$ like $r^{-2d}$ with an appropriate {\it scaling dimension\/}~$d$.
From physics we know that critical points always lie on surfaces that
separate phases in the system. But here we have a situation where the
critical points form a region on the phase diagram. One can expect,
that if we add to the model new configurations that break the ice
conditions, this critical region will become a surface that separates
two phases.

Unfortunately, physically reasonable solvable generalizations of the
six-vertex model are unknown. If we relate some energy to disbalance of
charge at a vertex, we shall lose solvability. But it is possible to
generalize the model in a `mathematical' way without lost of
solvability as follows~\cite{Sutherland,FanWu}. Suppose the
configurations around a vertex to be admissible if
\eq{
\ve_1+\ve_2=\ve'_1+\ve'_2\bmod4.
\label{EightVertexRestriction}
}
It means that we admit two more vertex configurations:
$$
\underbrace{
\vcenter{
\beginpicture
\setcoordinatesystem units <1cm,1cm> point at 0 0
\rl -0.5 0.0 , 0.5 0.0
\rl 0.0 -0.5 , 0.0 0.5
\put{$+$} [tc] <0pt,-3pt> at 0.0 -0.5
\put{$+$} [Br] <-3pt,-2pt> at -0.5 0.0
\put{$-$} [Bc] <0pt,3pt> at 0.0 0.5
\put{$-$} [Bl] <3pt,-2pt> at 0.5 0.0
\endpicture}
\qquad
\vcenter{
\beginpicture
\setcoordinatesystem units <1cm,1cm> point at 0 0
\rl -0.5 0.0 , 0.5 0.0
\rl 0.0 -0.5 , 0.0 0.5
\put{$-$} [tc] <0pt,-3pt> at 0.0 -0.5
\put{$-$} [Br] <-3pt,-2pt> at -0.5 0.0
\put{$+$} [Bc] <0pt,3pt> at 0.0 0.5
\put{$+$} [Bl] <3pt,-2pt> at 0.5 0.0
\endpicture}
}_{\atop{\textstyle d}}
$$
Such model is called {\it eight-vertex}. If, in addition to
(\ref{abcdef}), the corresponding energies $E_7$ and $E_8$ are equal,
$$
d=e^{-E_7/T}=e^{-E_8/T},
$$
the model turns out to be solvable~\cite{Baxter72}. The properties
described below can be found in~\cite{Baxter}.

The phase diagram of the eight-vertex model looks like:
\eq{
\Aligned{
{\rm AF_1}:
&\quad
c>a+b+d,
\\
{\rm AF_2}:
&\quad
d>a+b+c,
\\
{\rm F_1}:
&\quad
a>b+c+d,
\\
{\rm F_2}:
&\quad
b>a+c+d,
\\
{\rm D}:
&\quad
{\textstyle{1\over2}}(a+b+c+d)>a,b,c,d.
}\label{8vPhaseDiagram}
}
In this case the disordered region D is not critical. The correlation
length is finite and correlation functions decrease exponentially. The
critical points lie on the boundaries of the regions
$$
a=b+c+d;
\quad
b=a+c+d;
\quad
c=a+b+d;
\quad
d=a+b+c
$$
and at the special surfaces
$$
\Aligned{
&a=0,\quad {\textstyle{1\over2}}(b+c+d)\ge b,c,d;
\\
&b=0,\quad {\textstyle{1\over2}}(a+c+d)\ge a,c,d;
\\
&c=0,\quad {\textstyle{1\over2}}(a+b+d)\ge a,b,d;
\\
&d=0,\quad {\textstyle{1\over2}}(a+b+c)\ge a,b,c;
}
$$
There are special maps (`dualities') that can map each of the regions
(\ref{8vPhaseDiagram}) on to another~\cite{FanWu}. So it is enough to
study only one region, e.g.\ $\rm AF_1$. In this region the ground
states again look like~(\ref{GroundStates}).

I repeated many times the words `solvability', `solvable'. What does
they mean? Though there is no rigorous definition of this term, it must
mean approximately the following: there are some quantities of physical
importance (partition function, correlation functions) that could be
found exactly in these models. What are the conditions for solvability
of this or that model? Let us make several steps to see such (presumably
sufficient) condition.

First introduce the weight matrix (called also $R$ matrix) of the
eight-vertex model:
\eq{
R^{\ve'_1\ve'_2}_{\ve_1\ve_2}=
\vcenter{\beginpicture
\setcoordinatesystem units <1cm,1cm> point at 0 0
\rl -0.5 0.0 , 0.5 0.0
\rl 0.0 -0.5 , 0.0 0.5
\ar -0.49 0.0 , -0.5 0.0
\ar 0.0 -0.49 , 0.0 -0.5
\put{$\ve_1$} [tc] <0pt,-3pt> at 0.0 -0.5
\put{$\ve_2$} [Br] <-3pt,-2pt> at -0.5 0.0
\put{$\ve'_1$} [Bc] <0pt,3pt> at 0.0 0.5
\put{$\ve'_2$} [Bl] <3pt,-2pt> at 0.5 0.0
\endpicture}
\label{RmatrixDef}
}
We added the arrows to the lines to define the orientation on the lattice,
so that we could rotate or distort it. The $R$ matrix can be written as 
\eq{
R=\bordermatrix{&++&+-&-+&--\cr
++&a&&&d\cr+-&&b&c&\cr-+&&c&b&\cr--&d&&&a}.
\label{RmatrixForm}
}

Now introduce the $L$ matrix
\eq{
(L^{\mu'}_\mu)^{\ve'_1\ldots\ve'_N}_{\ve_1\ldots\ve_N}
=\vcenter{
\beginpicture
\setcoordinatesystem units <1cm,1cm> point at 0 0
\rl 0.0 -2.5 , 0.0 -1.3
\setdashes
\rl 0.0 -1.3 , 0.0 -0.2
\setsolid
\rl 0.0 -0.2 , 0.0 0.5 \ar 0.0 -2.49 , 0.0 -2.5
\rl -0.5 0.0 , 0.5 0.0 \ar -0.49 0.0 , -0.5 0.0
\rl -0.5 -1.5 , 0.5 -1.5 \ar -0.49 -1.5 , -0.5 -1.5
\rl -0.5 -2.0 , 0.5 -2.0 \ar -0.49 -2.0 , -0.5 -2.0
\put{$\mu$} [tc] <0pt,-3pt> at 0.0 -2.5
\put{$\mu'$} [Bc] <0pt,3pt> at 0.0 0.5
\put{$\ve_N$} [Br] <-3pt,-2pt> at -0.5 0.0
\put{$\ve'_N$} [Bl] <3pt,-2pt> at 0.5 0.0
\put{$\ve_2$} [Br] <-3pt,-2pt> at -0.5 -1.5
\put{$\ve'_2$} [Bl] <3pt,-2pt> at 0.5 -1.5
\put{$\ve_1$} [Br] <-3pt,-2pt> at -0.5 -2.0
\put{$\ve'_1$} [Bl] <3pt,-2pt> at 0.5 -2.0
\endpicture}
=R_{\mu_{N-1}}^{\mu'}{\vphantom{R}}_{\ve_N}^{\ve'_N}\ldots
R_{\mu_1}^{\mu_2}{\vphantom{R}}_{\ve_2}^{\ve'_2}
R_\mu^{\mu_1}{\vphantom{R}}_{\ve_1}^{\ve'_1}.
\label{LDef}
}
We shall consider this object as a matrix in indices $\mu$, $\mu'$ and
an operator in the product $\C^2\otimes\C^2\otimes\cdots\otimes\C^2$
($N$ times) spanned on the vectors $v_{\ve_1}\otimes
v_{\ve_2}\otimes\cdots\otimes v_{\ve_N}$, where $v_+$ and $v_-$ form
the natural basis in $\C^2$. The product associated to
$\ve_1,\ldots,\ve_N$ is called {\it quantum\/} space, while the space
$\C^2$ associated to $\mu$ is called {\it auxiliary\/} space. We shall
always omit the `quantum' indices and sometimes omit the auxiliary
indices, substituting them by a numeric subscript labeling the space,
e.g.~$L_1$.

Now we are ready to define the {\it transfer matrix}
\eq{
T=\tr_1L_1\equiv\sum_\mu L^\mu_\mu.
\label{TDef}
}
Consider the eight-vertex model of $M$ columns and $N$ rows with cyclic
boundary conditions. It is easy to see that the partition function is
equal to
$$
Z=\Tr T^M,
$$
where the trace $\Tr$ is taken over the quantum space. Let $\Lambda_1\ge
\Lambda_2\ge\cdots\ge\Lambda_{2^N}$ be the eigenvalues of the transfer
matrix. Then
$$
Z=\sum_i\Lambda_i^M.
$$
In the large $M$ limit the leading contribution is given by
$\Lambda_1^M$ and we have
$$
\kappa\equiv Z^{1/MN}\to\Lambda_1^{1/N},
\qquad
f\equiv{F\over MN}\to-{T\over N}\log\Lambda_1
\qquad\text{as $M\to\infty$},
$$
where the partition function per site $\kappa$ and the partial free
energy $f$ are introduced. Surely, it is necessary to study the behavior
of subleading contributions at large $M$ and $N$ to substantiate these
formulas, and it can be done, but in these lectures we shall assume
these formulas to be correct without a proof.

Anyhow the problem is reduced to that of the quantum mechanics: we
have an operator $T$ of evolution by one step in the `time' direction
and we have to diagonalize it. When this can be done exactly? From
classical mechanics we know that the system is solvable (more
precisely, integrable) if we have sufficiently many integrals of motion in
involution (the Liouville theory). Though there is no `quantum
Liouville theorem' it can be expected that in quantum mechanics the
situation is similar. We have to find some other commutative integrals of motion,
i.e.\ operators that commute with the transfer matrix and with each
other. Let us look for them in the same form as the transfer matrix.
Namely, let $R'$ be the matrix of the form (\ref{RmatrixForm}) with some
new `weights' $a'$, $b'$, $c'$, $d'$. Define the $L$ operator $L'$ and
the transfer matrix $T'$ in terms of $R'$ according to (\ref{LDef}) and
(\ref{TDef}). Besides, for any product $V_1\otimes\cdots\otimes
V_n=\C^2\otimes\cdots\C^2$ we shall denote by $R_{ij}$ the matrix $R$
acting on the $i$th and $j$th components of the product. Then there is
a

{\itshape{\bfseries Theorem~\cite{Baxter71}.} If there exists an
invertible matrix $R''$ of the form (\ref{RmatrixForm}), such that the
{\bfseries Yang--Baxter equation} holds
\eq{
R''_{12}R'_{13}R^{\vphantom{\prime}}_{23}
=R^{\vphantom{\prime}}_{23}R'_{13}R''_{12}
\label{YBE}
}
or, graphically,
$$
\vcenter{
\beginpicture
\setcoordinatesystem units <1cm,0.87cm> point at 0 0
\rl 1.0 -0.5 , -1.0 -0.5 \ar -0.99 -0.5 , -1.0 -0.5
\ar 0.25 1.0 , -0.75 -1.0
\ar -0.25 1.0 , 0.75 -1.0
\put{$1$} [tc] <0pt,-3pt> at 0.75 -1.0
\put{$2$} [tc] <0pt,-3pt> at -0.77 -1.0
\put{$3$} [Br] <-3pt,-2pt> at -1.0 -0.5
\put{$R'$} [tr] <-1pt,-3pt> at 0.5 -0.5
\put{$R^{\vphantom{\prime}}$} [tl] <1pt,-3pt> at -0.5 -0.5
\put{$R''$} [Bl] <3pt,-2pt> at 0.0 0.5
\endpicture}
\quad=\quad
\vcenter{
\beginpicture
\setcoordinatesystem units <1cm,0.87cm> point at 0 0
\rl 1.0 0.5 , -1.0 0.5 \ar -0.99 0.5 , -1.0 0.5
\ar 0.75 1.0 , -0.25 -1.0
\ar -0.75 1.0 , 0.25 -1.0
\put{$1$} [tc] <0pt,-3pt> at 0.25 -1.0
\put{$2$} [tc] <0pt,-3pt> at -0.25 -1.0
\put{$3$} [Br] <-3pt,-2pt> at -1.0 0.5
\put{$R$} [Br] <-1pt,3pt> at 0.5 0.5
\put{$R'$} [Bl] <1pt,3pt> at -0.5 0.5
\put{$R''$} [Bl] <3pt,-2pt> at 0.0 -0.5
\endpicture}
$$
then the transfer matrices $T$ and $T'$ commute.
\eq{
TT'=T'T.
\label{TT'commut}
}
}

The proof is elementary in graphical form:
$$
T'T=\
\vcenter{
\beginpicture
\setcoordinatesystem units <0.5cm,0.5cm> point at 0 0
\rl 0 4.5 , 0 3.5 \rl 0 0.5 , 0 -0.78 \ar 0 -0.78 , 0 -0.8
\rl 1 4.5 , 1 3.5 \rl 1 0.5 , 1 -0.78 \ar 1 -0.78 , 1 -0.8
\setdashes
\rl 0 3 , 0 1
\rl 1 3 , 1 1
\setsolid
\rl 1.5 0 , -0.78 0 \ar -0.78 0 , -0.8 0
\rl 1.5 4 , -0.78 4 \ar -0.78 4 , -0.8 4
\put{$L^{\vphantom{\prime}}$} [tc] <0pt,-3pt> at 0.0 -0.8
\put{$L'$} [tc] <0pt,-3pt> at 1.0 -0.8
\endpicture}
\ =\
\vcenter{
\beginpicture
\setcoordinatesystem units <0.5cm,0.5cm> point at 0 0
\ln 0 4.5 , 1 5.5 \ln 1 5.5 , 0 6.5
\ln 1 4.5 , 0 5.5 \ln 0 5.5 , 1 6.5
\put{$R''$} [Bl] <6pt,-2pt> at 0.5 5
\put{\rlap{$R^{\prime\prime-1}$}} [Bl] <6pt,-2pt> at 0.5 6
\rl 0 4.5 , 0 3.5 \rl 0 0.5 , 0 -0.78 \ar 0 -0.78 , 0 -0.8
\rl 1 4.5 , 1 3.5 \rl 1 0.5 , 1 -0.78 \ar 1 -0.78 , 1 -0.8
\setdashes
\rl 0 3 , 0 1
\rl 1 3 , 1 1
\setsolid
\rl 1.5 0 , -0.78 0 \ar -0.78 0 , -0.8 0
\rl 1.5 4 , -0.78 4 \ar -0.78 4 , -0.8 4
\put{$L^{\vphantom{\prime}}$} [tc] <0pt,-3pt> at 0.0 -0.8
\put{$L'$} [tc] <0pt,-3pt> at 1.0 -0.8
\endpicture}
\ =\
\vcenter{
\beginpicture
\setcoordinatesystem units <0.5cm,0.5cm> point at 0 0
\ln 1 4.5 , 0 5.5 \ln 0 4.5 , 1 5.5
\ln 0 -0.5 , 1.18 -1.68 \ar 1.18 -1.68 , 1.2 -1.7
\ln 1 -0.5 , -0.18 -1.68 \ar -0.18 -1.68 , -0.2 -1.7
\put{$R''$} [Bl] <6pt,-2pt> at 0.5 -1
\put{\rlap{$R^{\prime\prime-1}$}} [Bl] <6pt,-2pt> at 0.5 5
\rl 0 4.5 , 0 3.5 \rl 0 0.5 , 0 -0.5
\rl 1 4.5 , 1 3.5 \rl 1 0.5 , 1 -0.5
\setdashes
\rl 0 3 , 0 1
\rl 1 3 , 1 1
\setsolid
\rl 1.5 0 , -0.78 0 \ar -0.78 0 , -0.8 0
\rl 1.5 4 , -0.78 4 \ar -0.78 4 , -0.8 4
\put{$L'$} [Br] <-3pt,3pt> at 0.0 -0.0
\put{$L^{\vphantom{\prime}}$} [Bl] <3pt,3pt> at 1.0 -0.0
\endpicture}
\ =\
\vcenter{
\beginpicture
\setcoordinatesystem units <0.5cm,0.5cm> point at 0 0
\ln 0 4.5 , 1 5.5 \ln 1 5.5 , 0 6.5
\ln 1 4.5 , 0 5.5 \ln 0 5.5 , 1 6.5
\put{\rlap{$R^{\prime\prime-1}$}} [Bl] <6pt,-2pt> at 0.5 5
\put{$R''$} [Bl] <6pt,-2pt> at 0.5 6
\rl 0 4.5 , 0 3.5 \rl 0 0.5 , 0 -0.78 \ar 0 -0.78 , 0 -0.8
\rl 1 4.5 , 1 3.5 \rl 1 0.5 , 1 -0.78 \ar 1 -0.78 , 1 -0.8
\setdashes
\rl 0 3 , 0 1
\rl 1 3 , 1 1
\setsolid
\rl 1.5 0 , -0.78 0 \ar -0.78 0 , -0.8 0
\rl 1.5 4 , -0.78 4 \ar -0.78 4 , -0.8 4
\put{$L'$} [tc] <0pt,-3pt> at 0.0 -0.8
\put{$L^{\vphantom{\prime}}$} [tc] <0pt,-3pt> at 1.0 -0.8
\endpicture}
\ =\
\vcenter{
\beginpicture
\setcoordinatesystem units <0.5cm,0.5cm> point at 0 0
\rl 0 4.5 , 0 3.5 \rl 0 0.5 , 0 -0.78 \ar 0 -0.78 , 0 -0.8
\rl 1 4.5 , 1 3.5 \rl 1 0.5 , 1 -0.78 \ar 1 -0.78 , 1 -0.8
\setdashes
\rl 0 3 , 0 1
\rl 1 3 , 1 1
\setsolid
\rl 1.5 0 , -0.78 0 \ar -0.78 0 , -0.8 0
\rl 1.5 4 , -0.78 4 \ar -0.78 4 , -0.8 4
\put{$L'$} [tc] <0pt,-3pt> at 0.0 -0.8
\put{$L^{\vphantom{\prime}}$} [tc] <0pt,-3pt> at 1.0 -0.8
\endpicture}
\ =TT'.
$$
Note, that the third equality means that
\eq{
R''_{12}L'_1L^{\vphantom{\prime}}_2=L^{\vphantom{\prime}}_2L'_1R_{12}.
\label{RLL}
}
which generalizes the Yang--Baxter equation to the situation where the
space $3$ is substituted by the whole quantum space.

What are the solutions to the Yang--Baxter equation? It turns out that they
can be written in the form
$$
R=R(u_2-u_3),
\qquad
R'=R(u_1-u_3),
\qquad
R''=R(u_1-u_2)
$$
with some analytic function $R(u)$. The {\it spectral parameters\/}
$u_i$ can be associated with the spaces $V_i$. The Yang--Baxter equation
takes the form
\eq{
R_{12}(u_1-u_2)R_{13}(u_1-u_3)R_{23}(u_2-u_3)
=R_{23}(u_2-u_3)R_{13}(u_1-u_3)R_{12}(u_1-u_2).
\label{YBEu}
}
The solution matrix elements $a(u),\ldots,d(u)$ of the matrix $R(u)$ are
written in terms of the Jacobi theta functions $\theta_i(u;\tau)$
($i=1,\ldots,4$) with quasiperiods $1$ and $\tau$ ($\Im\tau>0$). Namely,
in the region $\rm AF_1$ we have
\eq{
\Aligned{
a(u)
&=\rho(u;\epsilon,r)s(1-u;\epsilon,r),
\cr
b(u)
&=\rho(u;\epsilon,r)s(u;\epsilon,r),
\cr
c(u)
&=\rho(u;\epsilon,r)s(1;\epsilon,r),
\cr
d(u)
&=\rho(u;\epsilon,r)s(1-u;\epsilon,r)s(u;\epsilon,r)s(1;\epsilon,r),
}\label{abcd-ell}
}
where
$$
s(u;\epsilon,r)={\theta_1\!\left({u\over2r};{i\pi\over2\epsilon r}\right)
 \over\theta_2\!\left({u\over2r};{i\pi\over2\epsilon r}\right)}.
$$
The parameters $\epsilon>0$, $r>1$ are fixed numbers. The scalar
function $\rho(u;\epsilon,r)$ is arbitrary.

In the limit $r\to\infty$ we obtain the six-vertex model:
$$
\Aligned{
a(u)
&\sim\sh\epsilon(1-u),
\\
b(u)
&\sim\sh\epsilon u,
\\
c(u)
&\sim\sh\epsilon,
\\
d(u)
&=0.
}
$$

The presence of a continuous set of solutions $R(u)$ means that there is
an infinite family of commuting transfer matrices:
\eq{
T(u_1)T(u_2)=T(u_2)T(u_1).
\label{TTu}
}
The matrix $R(0)$ is proportional to the transposition matrix
$$
P=\pMatrix{1&&&\\&1&1&\\&1&1&\\&&&1}.
$$
It means that the transfer matrix $T(0)$ is proportional to the shift
operator. Define now a set of Hamiltonians $H_1,H_2,\ldots$ as follows
\eq{
T^{-1}(0)T(u)=1+H_1u+H_2u^2+\ldots.
\label{HnDef}
}
They all commute with each other and with the shift operator
$$
[H_m,H_n]=0,
\qquad
[H_m,T(0)]=0.
$$
Not all of them are independent. In a finite system only a finite number
of them are independent. It turns out that the model is indeed solvable.
It is not easy to prove this, and it demands some other ideas to find
even the partition function. You can find the solution in Baxter's
book~\cite{Baxter}.

What are the Hamiltonians $H_n$? The general $H_n$ is a complicated
operator, but the simplest one is given by
$$
H_1=\const
-{1\over2}\sum^N_{k=1}(J_x\sigma^x_k\sigma^x_{k+1}
+J_y\sigma^y_k\sigma^y_{k+1}+J_z\sigma^z_k\sigma^z_{k+1})
$$
with $\sigma^a_k$ are the Pauli sigma-matrices,
$$
\sigma^x=\pMatrix{&1\\1&},
\qquad
\sigma^y=\pMatrix{&i\\-i&},
\qquad
\sigma^z=\pMatrix{1&\\&1},
$$
acting on the $k$th component of the tensor product. The coefficients
$J_a$ are functions of $\epsilon$, $r$. The common factor is not so
interesting, but the ratios of these coefficients are important $u$
independent combinations of weights:
$$
\Delta={2J_z\over J_x+J_y}={a^2+b^2-c^2-d^2\over2ab},
\qquad
\Gamma={J_x-J_y\over J_x+J_y}={cd\over ab}.
$$
In the case of the six-vertex model $\Gamma=0$ and we have
the XXZ model ($J_x=J_y$).

It can be said that the XYZ chain is much more physical model than the
eight-vertex model itself. The largest eigenvalue of $T$ corresponds to
the lowest eigenvalue of $H_1$. The next-to-large eigenvalues of $T$
correspond to the first excitations above the ground state of the XYZ
model. In the infinite-volume limit the lowest two states are
degenerate, while the gap between these lowest states and the other
excitations remains finite and only vanishes at the critical points.
This gap is the inverse correlation length in the time dimension, while
the spectrum above the gap defines the correlation length in the spatial
dimension.

Generally, the two-dimensional lattice models of classical statistical
mechanics are related to the one-dimensional models of quantum
mechanics.

If the function $\rho$ satisfies the conditions
\eq{
\rho(u)\rho(-u)=(s(u;\epsilon,r)s(-u;\epsilon,r)+s^2(1;\epsilon,r))^{-1},
\qquad
\rho(u)=\rho(1-u),
\label{rhoConditions}
}
the $R$ matrix satisfies two additional conditions:
\Align{
\span
R_{12}(u)R_{21}(-u)
=\id
\qquad
&&\text{\bf(Unitarity)},
\label{Unitarity}
\\
\span
R(u)_{\ve_1\ve_2}^{\ve'_1\ve'_2}
=R(1-u)_{\ve'_2,-\ve_1}^{\ve_2,-\ve'_1}
\qquad
&&\text{\bf(Crossing symmetry)}.
\label{CrossingSymmetry}
}

It is possible to find the solution to the equations
(\ref{rhoConditions}), such that $R(u)$ has the minimal number of poles
on the strip $0<\Re u<1$. It reads
\Align{
\rho(u;\epsilon,r)
&=x^{1-r/2}
{(x^{2r+2};x^{4r})_\infty(x^{2r-2};x^{4r})_\infty
\over(x^{2r};x^{4r})^2_\infty}
{(x^{2r}z;x^{4r})_\infty(x^{2r}z^{-1};x^{4r})_\infty
\over(x^{4r-2}z;x^{4r})_\infty(x^2z^{-1};x^{4r})_\infty}
{g(z^{-1})\over g(z)},
\notag
\\
g(z)
&\equiv g(z;\epsilon,r)
={(x^2z;x^4,x^{2r})_\infty(x^{2r+2}z;x^4,x^{2r})_\infty
\over
(x^4z;x^4,x^{2r})_\infty(x^{2r}z;x^4,x^{2r})_\infty},
\label{rhoFactor}
}
with the brace function
$$
(z;p_1,\ldots,p_N)
=\prod_{n_1,\ldots,n_N=0}^\infty(1-zp_1^{n_1}\ldots p_N^{n_N})
$$
and the `multiplicative' parameters $x$, $p$, $z$ defined as
$$
x=e^{-\epsilon},
\qquad
p=x^{2r},
\qquad
z=x^{2u}.
$$
It turns out that this solution gives just the $R$ matrix for which the
partition function per site is equal to~1 according to Baxter's
solution. It means that solutions of such simple `reflection equations'
(\ref{Unitarity}) and (\ref{CrossingSymmetry}) make it possible to
easily reproduce the result of tedious and involved calculations based
on the Bethe equations!


\begin{thebibliography}{9}
\itemsep=0pt\parskip=0pt
\bibitem{Ruelle}
D.~Ruelle, Statistical mechanics. Rigorous results.
World Scientific, 1999
\bibitem{Sinai}
Ya.~G.~Sinai, Theory of phase transitions, Pergamon Press, 1983
\bibitem{Lieb}
E.~H.~Lieb, Phys.\ Rev.\ {\bf162} (1967) 162
\bibitem{Sutherland}
B.~Sutherland, J.~Math.\ Phys.\ {\bf11} (1970) 3183
\bibitem{FanWu}
C.~Fan, F.~Y.~Wu, Phys.\ Rev.\ {\bf B2} (1970) 723
\bibitem{Baxter72}
R.~J.~Baxter, Annals Phys.\ {\bf70} (1972) 193
\bibitem{Baxter}
R.~J.~Baxter, Exactly solved models in statistical mechanics,
Academic Press, 1982
\bibitem{Baxter71}
R.~J.~Baxter, Stud.\ Appl.\ Math.\ (Mass.\ Inst.\ of Technology) {\bf50}
(1971) 51
\end{thebibliography}












\end{document}

