Statistical mechanics of plasma

Thanks go to Salvatore Pace (MIT) for the notes that I used to prepare this post.

Firstly, notation. Statistical mechanics is notoriously challenging to wrap one’s head around due to the need to index all of the particles in a macroscopic system, which inevitably leads to some lengthy expressions and tortured notations. I am conciously deviating from the canonical notation here in the vain hope of making the presentation marginally cleaner and less reliant on knowledge of the canon.

Consider the grand canonical partition function of a gas of several species of interacting, specieswise indistinguishable, classical particles. Let the index set for all particles in the gas be \(I(\underline{N}) = \{1, ..., N_1\}\times ... \times \{1, ..., N_S\}\). $S$ is simply the number of species present. The index \(\underline{N} = (N_a)_{a \in \text{Species}}\) is a vector in $\mathbb{N}^S$ representing the number of particles in each species. $\mathbb{N}$ includes 0.

\[\mathcal{Z}(\beta, (\mu_a)_{a \in \text{Species}}) = \sum_{\underline{N} \in \mathbb{N}^S} \int \prod_{a\in\text{Species}} \left\{ \frac{d^{N_a}\Omega_{a}}{N_a!} \right\} \exp\left(-\beta H[ (r, p)_{j \in I(\underline{N})} ] + \sum_{b\in\text{Species}} \beta \mu_b N_b \right)\]

I’ve notated a volume of phase space by $d\Omega = \frac{dr dp}{h^d}$. Throughout this text I will use $dr$ to mean a top form, i.e. a $d$ dimensional volume element. Classically, $h$ is needed (with units [xp] = [Energy][Time], considering $p$ to be canonically conjugate momenta). Note that this constant must be present if we want the summand to have consistent dimensions. The factorials account for overcounting due to indistinguishability.

Now, let $H$ describe classical particles with Coulomb interactions,

\[H[ (r, p)_{j \in I(\underline{N})} ] = \sum_{i \in I(\underline{N})} \frac{p_i^2}{2 m_i} + \frac{1}{2}\sum_{i\neq j \in I(\underline{N})} e_i e_j V(r_i - r_j)\]

This problem is too hard, but we can make progress using perturbation theory in the weak coupling limit. We will consider the interactions as a perturbation to the ideal gas partition function, i.e. the free theory:

\[\mathcal{Z}_0(\beta, \underline{\mu}) = \sum_{\underline{N} \in \mathbb{N}^S} \prod_{a\in\text{Species}} \left\{ \int \frac{d^{N_a}\Omega_{a}}{N_a!} \exp \left( \sum_{j=1}^{N_a} \left[ -\frac{\beta p_j^2}{2 m_a}\right] + \beta \mu_a N_a \right) \right\}\] \[= \prod_{a\in\text{Species}} \exp[V \lambda_a^{-d} z_a]\]

where the final expression is written in terms of the de Broglie wavelength $\lambda_a^2 = 2\pi \hbar^2 \beta/m_a$ and fugacity $z_a = \exp(\beta \mu_a)$.

The expected number desnsities of the free theory are then given by

\[\frac{1}{V} \frac{\partial F}{\partial \mu_a} = \int \mathcal{D}\Omega\ \frac{N_a}{V} \frac{\exp(-\beta H + \beta \mu_b N_b) }{\mathcal{Z}_0} = \langle n_a \rangle_0 = \lambda_a^{-3}z_a\]

The full partition function may then be simplified to

\[\mathcal{Z}(\beta, \underline{\mu}) = \sum_{\underline{N} \in \mathbb{N}^S} \int \prod_{a\in\text{Species}} \left\{ \frac{\left[dr_a n_{0,a} \right]^{N_a} }{N_a!} \right\} \exp \left( -\frac{\beta}{2 } \sum_{i\neq j} \int_{M} dr dr' \rho_i(r)\rho_j(r') V(r - r')\right)\]

In this expression, we have also introduced an electric charge density distribution $\rho_j(r) = e_j \delta^d(r-r_j)$. To avoid the unpleasant question of what happens when $r=r’$, generally speaking a singularity of $V$, we tactfully remove all points from the domain $M$ in which $r=r’$. (It may become necessary to actually remove an $\epsilon$ ball around these singularities, but we’ll cross that bridge when we come to it.) This trickery allows the sum in the integrand to factorise cleanly:

\[\mathcal{Z}(\beta, \underline{\mu}) = \sum_{\underline{N} \in \mathbb{N}^S} \int \prod_{a\in\text{Species}} \left\{ \frac{\left[dr_a n_{0,a} \right]^{N_a} }{N_a!} \right\} \exp \left( -\frac{\beta}{2 } \int dr dr' \left(\sum_i\rho_i(r)\right)\left(\sum_j\rho_j(r')\right) V(r - r')\right)\]

Hubbard-Stratonovich transformation

Advanced multiplication by 1

Reference: Atland and Simons pg.196

Using the might of linear algebra, it’s possible to show that for a positive definite symmetric real matrix $K$,

\[\int d^n x \exp\left( -\frac{1}{2} x^T K x + h \cdot x \right) = \sqrt{\frac{(2\pi)^n}{\det K}} \exp\left[ \frac{1}{2}h^TK^{-1} h\right]\]

Now suppose that each of the $x’s$ represents the value of a field $\phi$ evaluated on some kind of finite grid $\Gamma$ (suppose we have periodic boundaries). It’s perfectly legitimate to then rearrange these degrees of freedom into a form resembling the dreaded path integral,

\[1 = \frac{ \int \prod_{x \in \Gamma} \left\{d\phi_x\right\} \exp\left( -\frac{1}{2} \phi \cdot K[\phi] + h \cdot \phi \right) }{ \sqrt{\frac{(2\pi)^|\Gamma|}{\det K}} \exp\left[ \frac{1}{2}h \cdot K^{-1} [h]\right]} = \exp\left[ -\frac{1}{2} h \cdot K^{-1}[h]\right] \int \mathcal{D} \phi \exp\left( -\frac{1}{2} \phi \cdot K[\phi] + h \cdot \phi \right)\]

\(K[\psi](x) = \sum_{y \in \Gamma} K_{xy} \psi_y\) is a symmetric, positive definite operator (that can be used to define an inner product), and likewise $h\cdot \phi = \sum_{x\in \Gamma} h_x \phi_x$. $K^{-1}$ is the inverse of the correlation matrix $K_{xy}$.

The temptation is now to densen up the grid, ‘let $N\to infty$’ if you like. This is problematic for several reasons, the main one being that the denominator diverges, restricting us from making any kind of statement about the path integral in isolation. For this reason, we have absorbed this normalisation into the definition of $\mathcal{D}$. The essential point is that the denominator has no physical meaning - it depends only on the grid shape and the operator determinant of K, both of which diverge, which is fine so long as they diverge uniformly with respect to the physical parameters $\beta, \mu$. This being the case, the quantity $\mathcal{Z}|\Gamma|^{1/2} \det^{-1/2}K$ remains finite in the limit of infinite lattices.

Decoupling the classical plasma

The expression for the full plasma partition function, including interactions, is

\[\mathcal{Z}(\beta, \underline{\mu}) = \sum_{\underline{N} \in \mathbb{N}^S} \int \prod_{a\in\text{Species}} \left\{ \frac{\left[dr_a n_{0,a} \right]^{N_a} }{N_a!} \right\} \exp \left( -\frac{\beta}{2 } \int dr dr' \left(\sum_i\rho_i(r)\right)\left(\sum_j\rho_j(r')\right) V(r - r')\right)\]

The name of the game is to perform a decoupling of the $\rho \rho V$ term. The current integrand appearing in the partition function as we have so far defined it may be viewed as $\exp(-\frac{1}{2} h \cdot \hat{V}[h])$, where $\hat{V}$ convolves with $V(x)$ and $\cdot$ is the standard functional inner product. We need to find an operator for $H = V^{-1}$, something that satisfies

\[\hat{H}\left[\int dr V(r-r') \phi(r)\right](r) = \int dr V(r-r') \hat{H}[\phi](r) = \phi(r)\]

We are simply looking for the (integro-)differential operator $K^{-1}$ that $K$ is the Green’s function for. Since $K(r-r’)=\frac{\beta}{|r-r’|}$ is the Coulomb potential, it is easy to see that $K^{-1} = -\frac{1}{\beta} \nabla^2$. Regarded as an operator, $\nabla^2$ is local, while $V(r,r’)$ was not - we will only have one integral in the final argument of $\exp$.

Using a Hubbard-Stratonovich transformation with $h = i\rho = i \sum_j \rho_j$,

\[\exp \left( -\frac{\beta}{2 } \int dr dr' \left(\sum_i\rho_i(r)\right)\left(\sum_j\rho_j(r')\right) V(r - r')\right) = \int \mathcal{D} \phi \exp \left[ \int d^3 x \left\{-\frac{1}{2\beta} \phi \nabla^2 \phi + i \rho \phi \right\}\right]\]