‹ Abrams Group

2.4 Classical Statistical Mechanics

Analogous to the quasi-classical microcanonical paritition function of Eq. 7, here is the quasi-classical representation of the canonical partition function: \begin{equation} \label {eq:q_classical} Q_{\rm classical} = \frac {1}{h^{dN}N!}\int \int d{\bf r}^N d{\bf p}^N \exp \left [-\beta \mathscr {H}\left ({\bf r}^N,{\bf p}^N\right )\right ] \end{equation} \(\mathscr {H}\left ({\bf r}^N,{\bf p}^N\right )\) is the Hamiltonian function which computes the energy of a point in phase space. The probability of a point in phase space is represented as \begin{equation} P\left ({\bf r}^N,{\bf p}^N\right ) = \left (Q_{\rm classical}\right )^{-1} \exp \left [-\beta \mathscr {H}\left ({\bf r}^N,{\bf p}^N\right )\right ]. \end{equation} So, the general “sum-over-states” ensemble average of quantum statistical mechanics, first presented in Eq. 1, becomes an integral over phase space in classical statistical mechanics: \begin{equation} \label {eq:classical_ensemble_average} \left <G\right > = \frac {\displaystyle \int \int d{\bf r}^N d{\bf p}^N \exp \left [-\beta \mathscr {H}\left ({\bf r}^N,{\bf p}^N\right )\right ] G\left ({\bf r}^N,{\bf p}^N\right ) }{\displaystyle \int \int d{\bf r}^N d{\bf p}^N \exp \left [-\beta \mathscr {H}\left ({\bf r}^N,{\bf p}^N\right )\right ] }, \end{equation} where \(G\left ({\bf r}^N,{\bf p}^N\right )\) is the value of the observable \(G\) at phase space point \(\left ({\bf r}^N,{\bf p}^N\right )\). Before moving on, it is useful to recognize that we normally simplify this ensemble average by noting that, for a system of classical particles, the usual choice for the Hamiltonian has the form \begin{equation} \label {eq:classical_hamiltonian} \mathscr {H}\left ({\bf r}^N, {\bf p}^N\right ) = \mathscr {K}\left ({\bf p}^N\right ) + U\left ({\bf r}^N\right ) \end{equation} where \(\mathscr {K}\) is the kinetic energy, which is only a function of momenta, and \(U\) is the potential energy, which is only a function of position. The canonical partition function, \(Q\), can in this case be factorized: \begin{equation}\label {eq:q_fact} \begin {aligned} Q\left (N,V,T\right ) & = \frac {1}{h^{3N}N!} \left \{\int d{\bf p}^N \exp \left [-\beta \mathscr {K}\left ({\bf p}^N\right )\right ]\right \} \left \{\int d{\bf r}^N \exp \left [-\beta U\left ({\bf r}^N\right )\right ]\right \} \\ & = \left \{\frac {V^N}{h^{3N}N!} \int d{\bf p}^N \exp \left [-\beta \mathscr {K}\left ({\bf p}^N\right )\right ]\right \} \left \{V^{-N}\int d{\bf r}^N \exp \left [-\beta U\left ({\bf r}^N\right )\right ]\right \} \\ & = Q_{\rm ideal} Z \end {aligned} \end{equation} The quantity in the left-hand braces is the ideal gas partition function, because it corresponds to the case when the potential \(U\) is 0. (Note that we have multiplied and divided by \(V^N\); this is the equivalent of scaling the positions in the integration over positions.) The quantity in the right-hand braces is called the configurational partition function, \(Z\).

Because the kinetic energy \(\mathscr {K}\) has the simple form, \begin{equation} \mathscr {K}\left ({\bf p}^N\right ) = \sum _i \frac {{\bf p}_i^2}{2m_i}, \end{equation} where \(m_i\) is the mass of particle \(i\), the integral over particle momenta can be evaluated analytically: \begin{equation} \begin {aligned} \int d{\bf p}^N \exp \left [-\beta \mathscr {K}\left ({\bf p}^N\right )\right ] & = \prod _{i=1}^{3N}\int dp_i \exp \left (-\frac {p_i^2}{2mk_BT}\right ) \\ & = \left (2\pi m k_B T\right )^{3N/2}. \end {aligned} \end{equation} (We have assumed all particles have the same mass, \(m\); in the case of distinct masses, this is just a product of similar factors.)

\(Q_{\rm ideal}\) becomes \begin{equation} \begin {aligned} \frac {V^N}{N!h^{3N}} \int d{\bf p}^N \exp \left [-\beta \mathscr {K}\left ({\bf p}^N\right )\right ] & = \frac {V^N}{N!}\left (\sqrt {\frac {2\pi m k_B T}{h^2}}\right )^{3N} \\ = Q_{\rm ideal}\left (N,V,T\right ) & = \frac {V^N}{N!\Lambda ^{3N}} \end {aligned} \end{equation} where \(\Lambda \) is the de Broglie wavelength, a quantum-mechanical property of a particle inversely proportional to its momentum (and thus inversely proportional to the square root of temperature): \begin{equation} \label {eq:debroglie} \Lambda = \sqrt {\frac {h^2}{2\pi m k_BT}} \end{equation} As an example, for a hydrogen atom with mass 1 amu and at room temperature (298 K), \(\Lambda \approx \) 1 Å. The de Broglie wavelength limits the precision by which a particle’s position can be determined; for H atoms at room temperature, one is not permitted to specify their positions with a precision finer than about 1 ångstrom without violating the Heisenberg uncertainty principle of quantum mechanics. However, as we will see, in classical molecular simulations, we must lift this restriction, while never forgetting that this makes a classical representation of a molecule somewhat less realistic.

With the momentum degrees of freedom handled at finite temperature, when the observable \(G\) is a function of positions only, the ensemble average becomes a configurational average: \begin{equation} \label {eq:config_avg} \left <G\right > = Z^{-1}\int d{\bf r}^N \exp \left [-\beta U\left ({\bf r}^N\right )\right ] G\left ({\bf r}^N\right ). \end{equation} Note that the integration over momentum yields a factor \(Q_{\rm ideal}\) in both the numerator and denominator, and thus divides out. We can write this configurational average using a probability distribution, \(\rho _{NVT}\), as \begin{equation} \left <G\right > = \int d{\bf r}^N G\left ({\bf r}^N\right )\rho _{NVT}\left ({\bf r}^N\right ) \end{equation} where \begin{equation} \rho _{NVT}\left ({\bf r}^N\right ) \equiv Z^{-1}e^{-\beta U\left ({\bf r}^N\right )} \end{equation} is called the “canonical probability distribution.” As pointed out on p. 15 of Frenkel & Smit [1], Eq. 47 is “the starting point for virtually all classical simulations of many-body systems”; that is, it is the starting point for almost all simulations discussed in this course.