I couldn't find any simple texts explaining the Kosterlitz-Thouless transition. More specifically can someone explain the role of vortices in the transition.
edit: links explaining the transition in a simple manner are also appreciated. Also the explanation does not have to be mathematically rigorous, a qualitative explanation is just fine.
edit 2: To be clear, I know generally why it happens. I believe it has to do with the crystal lattice having 4 fold symmetry (lattice atoms arranged in square fashion), and thus, the atoms do not have an easy axis to align there magnetic moments. As a result, some weird vortices are created.
Best Answer
The Berezinskii-Kosterlitz-Thouless (BKT) scenario is one of the most beautiful transitions that is ubiquitous in 2d systems (though it can also occur in higher dimensions for particular kinds of models) that surprisingly requires non-perturbative effects (i.e. topological defects) to be realized. To understand all the fuss (and the nobel) around this transition, perhaps a bit of context might be useful.
There is a celebrated theorem in equilibrium statistical mechanics, the Mermin-Wagner-Hohenberg-Coleman theorem, that essentially tells us that a continuous symmetry cannot be broken spontaneously at any finite temperature in dimensions two or lower. This is because the goldstone modes generated upon breaking breaking a continuous symmetry have strong fluctuations in $d=1,2$ leading to the symmetry being restored at long distances (for $T>0$).
Now for a 2d superfluid or superconductor the relevant order parameter is a complex scalar field $\psi=|\psi|e^{i\phi}$ with a phase shift $U(1)$ symmetry. So immediately, one would imagine that the 2d superconducting or superfluidity transition would never occur at finite temperature (and hence these states would never exist in the thermodynamic limit). The same conclusion is reached for the XY ferromagnet ($O(2)$ classical spins on a 2d lattice) or a 2d nematic liquid crystal. What Kosterlitz and Thouless went on to show was that the theorem was true in that no continuous symmetry is spontaneously broken at finite temperature, but there still was a continuous phase transition (with a diverging correlation length) at some finite temperature in these systems. This is an important discovery as until then the Landau-Ginzburg paradigm used to describe continuous phase transitions and critical phenomena always associated spontaneous symmetry breaking with the transition (note that however it was rather well known that a first order transition didn't require any such symmetry breaking, c.f the regular Liquid-Gas transition). Later on, Polyakov extended this scenario to gauge theories (in the hope of describing confinement in QCD), resulting in some very nice work showing, for example, 2+1 "compact" QED has a gapped spectrum in the IR due to topological excitations (Phys. Lett. B 59, 1975, Nucl. Phys. B 120, 1977) and the SU($N$) Thirring model has a fermions condensing with finite mass in the IR without breaking the chiral symmetry of the theory (E. Witten, Nucl. Phys. B 145, 1978). It was also further extended by D. Nelson and B. Halperin in the context of 2d melting of crystalline solids (Phys. Rev. B 19, 1979) leading to the prediction of a new liquid crystalline hexatic phase.
After this very long preamble, let us now look at what the transition is really all about. The simplest model that exhibits the BKT transition is the XY model. Consider a 2d lattice with unit 2d vectors at each site. Each vector $\vec{S}_i$ (at site '$i$') being in the plane is specified by a single angle $\theta_i$ $$ \vec{S}_i=(\cos\theta_i,\sin\theta_i) $$ The model is now specified by the Hamiltonian of the system, which includes nearest neighbour interactions that prefer to align near by spins. In the absence of an external field, we have \begin{align} \beta\mathcal{H}&=-\dfrac{J}{k_B T}\sum_{\langle i,j\rangle}\vec{S}_i\cdot\vec{S}_j\\ &=-\dfrac{J}{k_B T}\sum_{\langle i,j\rangle}\cos(\theta_i-\theta_j) \end{align} where $J>0$ is the interaction coupling constant. Now at low temperatures, as the fluctuations in the angles are going to be small, at long distances, we take the continuum limit of the lattice model assuming the angle field is slowly varying. Therefore writing $\theta_i-\theta_j=a\nabla\theta(x)\cdot\hat{e}_{ij}+O(a^2)$, where $a\rightarrow0$ is the lattice spacing and $\hat{e}_{ij}$ is the unit vector along the lattice bond joining sites $i$ and $j$, we get $$ \beta\mathcal{H}_{\mathrm{cont.}}=\dfrac{\beta J}{2}\int\mathrm{d}^2x\ |\nabla\theta|^2 $$ For small fluctuations ($\theta\ll 1$), the fact that $\theta(x)$ is an angular variable is irrelevant, allowing us to compute the two point correlation function as the partition function is given by a gaussian integral. $$ \langle|\theta(q)|^2\rangle=\dfrac{k_B T}{J\ q^2} $$ Inverse Fourier transforming this, gives us \begin{gather} \langle\theta(x)^2\rangle=\dfrac{k_B T}{2\pi J}\ln\left(\dfrac{L}{a}\right)\\ \langle[\theta(x)-\theta(0)]^2\rangle=\dfrac{k_B T}{2\pi J}\ln\left(\dfrac{x}{a}\right) \end{gather} $L$ is the system size (IR cutoff) and $a$ the lattice spacing (UV cutoff). Hence as $L\rightarrow\infty$, $\langle\vec{S}(x)\rangle=0$ implying the absence of long-ranged order and $$ \langle\vec{S}(x)\cdot\vec{S}(0)\rangle=\left(\dfrac{x}{a}\right)^{-\frac{k_BT}{2\pi J}} $$ The two-point spin correlation goes to 0 as $x\rightarrow\infty$ denoting the absence of long ranged order (this is once again just the Mermin-Wagner theorem) though the decay is very slow. It is temperature dependent power law, instead of the usual exponential decay (with a finite correlation length) expected for a disordered phase. So the low temperature phase of the XY model has what is called quasi-long ranged order (QLRO) with infinite correlation length ($\xi=\infty$). As additional non-linearities coming from the gradient expansion can be shown to be irrelevant at long distances (in the RG sense), one is immediately led to believe that this power law decay persists for all temperatures. This is evidently wrong, as common sense (and also high temperature loop expansions of the lattice model) would tell us that at high temperatures, the interaction is irrelvant leaving each spin essentially independent and random, leading to decorrelations over a few lattice spacings.
The resolution is then obtained by noting that in forgetting the angular nature of $\theta(x)$, the continuum gaussian "spin wave" theory does not account for windings of the angular phase field from $0$ to $2\pi$. These are called vortices (and anti-vortices), and they correspond to topological defects in the $\theta(x)$ field (which is then not defined at the core of the defect). They are perfectly reasonable configurations on the lattice whose continuum limit corresponds to point singularities in the angle field.
Note that these configurations never appear in a perturbative gradient expansion and are hence non-perturbative in nature. At the continuum level, the vortex is a singular solution of the euler-lagrange equation.
$$
\nabla^2\theta=0\\
\oint_{\Gamma}\mathrm{d}s\cdot\nabla\theta=2\pi q
$$
where $\Gamma$ is closed loop encircling the origin and $q$ is the integer "charge" of the vortex. This basically says that as you go once around the origin, the phase field $\theta$ goes from $0$ to $2\pi q$ (which is the same as $0$ for a periodic function as $q$ is an integer). Now for a single such defect, we have $|\nabla\theta|=q/r$ ($r$ being the radial coordinate), we can compute its energy to be,
$$
E_q=\pi J q^2\ln\left(\dfrac{L}{a}\right)
$$
which diverges logarithmically in the thermodynamic limit. Hence single defects are never excited, but defect pairs with opposite charges (dipoles) have a finite energy and can be excited at finite temperature. Neglecting interactions for the time being, there is a very simple hand waving argument for the existence of a phase transition. The energy of a single free defect diverges, but at finite $T$ one must look at the free energy which includes entropic contributions too. The number of ways a single defect of size $\sim a^2$ can be placed in a area of $L^2$ is roughly $(L/a)^2$. taking the logarithm to get the entropy, we have for the free energy
$$
F=E_q-TS=(\pi J q^2-2 k_B T)\ln\left(\dfrac{L}{a}\right)
$$
As the lowest charge excitations correspond to $q=1$, we have for $T>T_c=\pi J/(2 k_B)$, the free energy becomes negative, which means that there is a proliferation of free defects in the system as entropy wins over the defect energetics. Including defect interactions doesn't change this picture (even $T_c$ remains the same). At $T=T_c$, one obtains a universal power law decay (upto log corrections)
$$
\langle\vec{S}(x)\cdot\vec{S}(0)\rangle=\left(\dfrac{x}{a}\right)^{-\eta}
$$
with $\eta(T_c)=1/4$. Above $T_c$, we have a finite correlation length ($\langle\vec{S}(x)\cdot\vec{S}(0)\rangle\sim e^{-x/\xi}$), which diverges exponentially fast as one approaches the transition form above.
So, here we have a model in which both low and high temperature phases are disordered, but there is a phase transition at finite $T$ that involves the proliferation and unbinding of pairs of topological defects. Thinking of the defects as electric charges, the transition is then from an insulating low temperature phase to a conducting plasma with freely moving ions at higher temperature.