$\newcommand{\I}{\mathrm{i}} \newcommand{\E}{\mathrm{e}} \newcommand{\D}{\mathop{}\!\mathrm{d}} \newcommand{\bra}[1]{\langle{#1}|} \newcommand{\ket}[1]{|{#1}\rangle} \newcommand{\braket}[1]{\langle{#1}\rangle} \newcommand{\bbraket}[1]{\langle\!\langle{#1}\rangle\!\rangle} \newcommand{\bm}{\boldsymbol}$

Lectures on statistical mechanics.

Homework assignments

Noninteracting systems

Use the maximum entropy principle to find the momentum distribution $f(p)$ of a particle of mass $m$, taking into account the constraint on its energy
$$E= \braket{\frac{p^2}{2m}}$$
Remember that the entropy is given by,
$$S[f] = - \int_{-\infty}^\infty \D p \, f(p) \ln f(p)\,.$$
Two level system at fixed temperature $T$ (canonical). We consider a gas of two level atoms whose (one particle) hamiltonian is
$$H = \epsilon \sigma_z\,,$$
where $\sigma_z = \mathrm{diag}(1,-1)$ is the $z$-pauli matrix. (We neglect the other degrees of freedom.)

Calculate the partition function and the free energy. Deduce the energy and the specific heat. Discuss the system’s behavior at temperatures close to zero. Calculate the entropy and discuss the significance of negative temperatures.
System of oscillators. Calculate the partition function in the canonical ensemble of a set of $N$ oscillators. The hamiltonian of each oscillator is,
$$H = \hbar \omega (n_x + 1/2)\,, \quad n_x = 0, 1,\ldots\,.$$
where $x = 1,2,,\ldots,N$ labels the oscillators, and $n_x$ is the energy quantum number (energy level). Calculate the specific heat (at constant volume), and discuss its properties as a function of the temperature.

Fermions

A degenerate gas of free fermions of mass $m$ and momentum $p = \hbar k$, at low temperature, has a finite density related to the fermi energy, $\epsilon_F$. We are interested in the finite temperature $T$ correction of the chemical potential

$$\mu = \epsilon_F \left[ 1 - \frac{\pi^2}{12}\left(\frac{T}{\epsilon_F} \right)^2\right]$$

Using,

$$n_k = \frac{1}{z^{-1}\E^{\epsilon_k/T} + 1},\quad z= \E^{\mu/T}, \; \epsilon_k =\frac{\hbar^2 k^2}{2m}\,,$$

for the number of fermions as a function of their energy $\epsilon_k$, demonstrate that the density is at $T=0$, given by,

$$\frac{N}{V} = 2 \int \frac{d^3 k}{(2\pi)^3} n_k = \frac{1}{3\pi^2} \left( \frac{2m}{\hbar^2}\right)^{3/2}\,\epsilon_F^{3/2}\,.$$

In order to find the temperature correction, use the previous expression to show that,

$$\frac{N}{V}= \frac{2^{1/2}}{\pi^2} \left( \frac{mT}{\hbar^2} \right)^{3/2} I = \frac{3}{2} \frac{N}{V} \left(\frac{T}{\epsilon_F}\right)^{3/2} I \,,$$

where the integral $I$, can be approximated by,

$$I = \int_0^\infty \D x \frac{x^{1/2}}{z^{-1}\E^x + 1}\approx \frac{2}{3} (\ln z)^{3/2} + \frac{\pi^2}{12}(\ln z)^{-1/2}\,.$$

Virial expansion

The interactions of argon atoms in a classical gas is roughly approximated by the potential,

$$w(r) = \begin{cases} \infty & \text{if } r < r_0 \\ -\epsilon & \text{if } r \in (r_0,r_1) \\ 0 & \text{if } r > r_1 \end{cases}\,,$$

where $r$ is the interparticle distance, $r_0$ is the repulsion radius, and $r_1$ the range of the attractive energy, which is of the order of $\epsilon$. Compute the second virial coefficient $B_2(T)$,

$$B_2(T) = 2\pi \int_0^\infty r^2 \D r\, \left(1 - \E^{-w(r)/T} \right)\,,$$

and discuss its behavior as a function of the temperature $T$.

Calculate the equation of state

$$PV=NT(1+B_2 N/V)$$

in the van der Waals form,

$$(P+n^2a)(1-nb)=nT, \quad n = N/V$$

and deduce the expression of the coefficients $a$ and $b$. Make a plot of $P = P(V)$ for different temperatures and discuss the result.

Ising 2D: Monte Carlo computation

The goal of this homework is to investigate the phase transition of the Ising model in two dimensions, and to compare it with the Onsager solution.

We follow the presentation of the lecture. The lattice size is $L^2$, and we take $a=1$. Use the Monte Carlo algorithm to

compute the magnetization per site
$$m = \frac{1}{L} \sum_x^{L^2} s_x$$
as a function of the temperature and fit the monte carlo data with the onsager formula
$$m(T) = \left[1 - \sinh^{-4} \frac{2}{T} \right]^{1/8}$$
determine the transition temperature
$$T_c = 2/\log(1+\sqrt{2})$$
show some typical configurations for different temperatures ($T < T_c$, $T \approx T_c$ and $T > T_c$)
calculate the energy and heat capacity, and discuss their behavior near the transition.

Your report should contain a brief description of the model and a scheme of the monte carlo algorithm (you may add in an appendix your code); the figures must be of publication quality (as generated by matplotlib with tex labels).

You may implement, in addition to the usual metropolis method, the cluster algorithm, which is faster, but it is not mandatory.

Electron plasma Debye-Hückel screening theory¹

The goal of this assignment is to study, using the virial cumulant expansion, the Debye-Hückel theory of an electron plasma.

A simple diluted plasma is a neutral gas of electrostatically interacting electrons in a background of ions. The electron mass and charge are noted $m$ and $-e$, respectively; the system hamiltonian is

\begin{equation}\label{e:H} H_N = H_0(\bm P) + W(\bm X) = \sum_{n=1}^N \frac{p_n^2}{2m} + \sum_{nm=1}^N w(|\bm x_n - \bm x_m|) \end{equation}

where the first term $H_0$ is the kinetic energy of $N$ electrons, and the second term $W$, the interaction energy, contains the two particles coulomb potential

\begin{equation} \label{e:w} w(r) = \frac{\alpha}{r} - w_0 \,, \quad \alpha = \frac{e^2}{4\pi \varepsilon_0} \end{equation}

where $r = |\bm x|$ is the distance between the two electrons, and we added a constant $w_0$ which takes into account the neutralizing ion background. It ensures that the integral

\begin{equation} \label{e:w0} \int_0^\infty 4\pi r^2 \D r \, w(r) = 0\,, \end{equation}

vanishes; we also defined the phase space coordinates $\bm X = (\bm x_1, \ldots, \bm x_N)$ and $\bm P = (\bm p_1, \ldots, \bm p_N)$.

The partition function,

\begin{equation} \label{e:Z} Z(T,V,N) = \frac{1}{N!} \int_{\mathbb{R}^N} \frac{ \D \bm P }{ (2\pi \hbar)^{3N} } \int_{V^N} \D \bm X \E^{ -H_0(\bm P)/T } \E^{ - W(\bm X)/T } \end{equation}

where $V$ is the system’s volume.

Show that
\begin{equation} \label{e:ZZ0} Z = Z_0 \braket{\E^{ - W(\bm x)/T }}\,,\quad \braket{O(\bm X)} = \int_{V^N} \D \bm X \, O(\bm X) \end{equation}
with $Z_0$ the ideal gas partition function ($W=0$) given by:
\begin{equation} \label{e:Z0} Z_0 = \frac{1}{N!} \left( \frac{V}{\lambda} \right)^{3N}\,, \quad \lambda = \sqrt{\frac{2\pi\hbar^2}{mT}}\,. \end{equation}
The brackets $\braket{\cdots}$ are for the mean over the non-interacting system probability distribution (which is the uniform probability in space):
$$\braket{\cdots} = \frac{1}{N!} \int_{\mathbb{R}^N} \frac{ \D \bm p }{ (2\pi \hbar)^{3N} } \int_{V^N} \D \bm X \,\frac{ \E^{ -H_0(\bm p)/T } }{Z_0} \, (\cdots) \,.$$

Ring graphs partition function

In the case of the coulomb interaction, because of its long range, the usual cumulant expansion breaks down. We need to compute terms with arbitrary powers of the density. We focus here on the ring graphs, which give, as we will demonstrate, a physically interesting picture of the diluted plasma high temperature plasma.

Therefore, we calculate now the contribution of these ring graphs to the logarithm of the partition function, neglecting the other possible graphs. This corresponds to the diluted limit.

Compute the fourier transform of the coulombian two particles interaction $w(\bm k)$, where $\bm k$ is the wavenumber, conjugate to $\bm x$. The neutrality condition ensures $w(0) = 0$. Show that for $\bm k \ne 0$,
\begin{equation} \label{e:fw} w(\bm k) = \frac{4\pi \alpha}{k^2} \end{equation}
Compute, using the convolution theorem, the ring graph of order $n \ge 2$:
\begin{equation} \label{e:Rn} R_n = \int_{V^n} \prod_{i=1}^n \frac{\D \bm x_i}{V}w(\bm x_1 - \bm x_2) w(\bm x_2 - \bm x_3) \ldots w(\bm x_n - \bm x_1)\,. \end{equation}
Show that the result of this integration is
\begin{equation} \label{e:Rn1} R_n = \frac{1}{V^{n-1}} \int \frac{\D \bm k}{(2\pi)^3} w(\bm k)^n\,. \end{equation}

The logarithm of the partition function

Explain the formula
\begin{equation} \label{e:Nn} N_R(n) = \frac{N!}{(N-n)!} \times (n-1)! \times \frac{1}{2} \end{equation}
of the number of ring graphs with $n$ points (the coordinates appearing in the binary interaction factors). Use the stirling formula to show that
\begin{equation} \label{e:Nn1} N_R(n) = \frac{(n-1)!}{2} N^n \end{equation}
The rings contribution to $\ln Z$ is
\begin{equation} \label{e:ZR} \ln Z_R = \ln Z_0 +\sum_{n=2}^\infty \frac{(-1)^n}{n!} N_R(n) \frac{R_n}{T^n} \end{equation}
Replacing the previous expression \eqref{e:Rn1} and \eqref{e:Nn1} demonstrate the relation
\begin{equation} \label{e:lnZR} \ln Z_R = \ln Z_0 + \frac{V}{12\pi} \kappa^2, \quad \kappa = \sqrt{\frac{4\pi \alpha n}{T}} \end{equation}
where $\kappa^{-1}$ is the Debye length and $n$ the density.
Compute the equation of state:
\begin{equation} \label{e:PVT} P = T \frac{\partial }{\partial V} \ln Z_R \end{equation}
and show that,
\begin{equation} \label{e:P} P = P_0 - \frac{T}{24 \pi} \kappa^3\,. \end{equation}
Discuss the pressure behavior as a function of the temperature. What happens when $\kappa n^{1/3} > 1$? (dense state limit).

Debye effective potential

The previous discussion shows that the collective effect of the coulomb interaction modifies the pressure, introducing a \emph{negative} correction of the order $\kappa n^{1/3}$. This is related to the long range of the coulomb potential. We ask now what is the effective interaction of two charges, in the same approximation (related to the ring graphs).

We define the effective potential by

\begin{equation} \label{e:veff} \overline{w}(\bm x - \bm y) = w(\bm x - \bm y) + \sum_{n=1}^\infty \left(\frac{-N}{TV}\right)^n \int_{V^N} \D \bm x_1 \ldots \D \bm x_n w(\bm x - \bm x_1) \ldots w(\bm x_n - \bm y) \end{equation}

which leads to the screened Debye potential:

\begin{equation} \label{e:debye} \overline{w}(\bm x) = \alpha \frac{\E^{-\kappa |\bm x|}}{|\bm x|}\,. \end{equation}

This formula takes into account the paths form $\bm x$ to $\bm y$ without loops; identifying the extreme points gives the ring graph.

To demonstrate this formula, show that the second order term (in $w^2$), can be written as,
\begin{equation} \label{e:weff2} \overline{w}_2(\bm x - \bm y) = w(\bm x - \bm y) - \frac{N}{VT} \int \frac{\D \bm k}{(2\pi)^3} w(\bm k)^2 \E^{\I \bm k \cdot (\bm x - \bm y) }\, \end{equation}
The generalization of \eqref{e:weff2} is straightforward:
\begin{equation} \label{e:weffn} \overline{w}(\bm x) = w(\bm x - \bm y) + \sum_{n=1}^\infty \left(\frac{N}{VT}\right)^n \int \frac{\D \bm k}{(2\pi)^3} w(\bm k)^{n+1} \E^{\I \bm k \cdot \bm x }\, \end{equation}
Compute the integral via the residue theorem to obtain the Debye potential \eqref{e:debye}.

(Optional question) Use the Poisson formula for the electrostatic potential, and consider the electrons distributed according to the Boltzmann distribution in a uniform ion background, to find the Debye length and the Debye potential.

The Debye-Hückel theory is discussed in the book by Kardar (Statistical Physics of Particles, problems section of chapter 5).

The file PS-hw3.pdf contains the complete text of the Debye-Hückel assignment, including some useful formulas. ↩