Statistical Mechanics: A Brief Introduction

2 Statistical Mechanics: A Brief Introduction

This course is centered upon a mathematical statement called an “ensemble average”: \begin{equation} \label {eq:ensemble_average} \left <G\right > = \sum _\nu P_\nu G_\nu \end{equation} That is, the expectation value, $\left <G\right >$, of some observable property $G$ is an average over all possible microstates available to a system, indexed by $\nu $, where $P_\nu $ is the probability of observing the system in microstate $\nu $, and $G_\nu $ is the value of the measured property $G$ when the system is in microstate $\nu $. Before even considering how to use computer simulation to make such a measurement of a particular property for a particular system, there are three main issues to consider:

1.: What is a microstate?
2.: What is meant by observing the system?
3.: How do we calculate probabilities?

In the following subsections, we give a cursory treatment of elmentary statistical mechanics aimed at answering these questions. The aim is to give the student an appreciation of the basic physics that underlies a majority of current molecular simulation.

2.1 Microstates and Degeneracy

A microstate is a full specification of all degrees of freedom of a system. In quantum mechanics, degrees of freedom are quantum numbers. The index $\nu $ in Eq. 1 runs over all unique combinations of quantum number values. Equilibrium (eigen)solutions of the Schrödinger equation define the energy $E_\nu $ of any state $\nu $: \begin{equation} \label {eq:energy_eigenstates} \mathscr {H}\left |\nu \right > = E_\nu \left |\nu \right >, \end{equation} where $\mathscr {H}$ is the Hamiltonian operator, $\left |\nu \right >$ is shorthand for the system wavefunction in state $\nu $, and $E_\nu $ is the energy of state $\nu $ ($E_\nu $ is an eigenvalue of the Schrödinger equation). In contrast to model systems usually considered in elementary quantum mechanics, the number of distinct microstates of systems with energy $E$ and comprising $\sim 10^{23}$ particles is very large, and this set of eigenstates is in practice impossible to obtain explicitly. This is indeed why we must instead treat this set statistically. We refer to the number of states that satisfy a given energy as the degeneracy of the energy level $E$, denoted $\Omega (E,N,V)$: \begin{equation} \label {eq:microcanonical} \Omega (N,V,E) = \begin {array}{l} \text {number of microstates with } N \text { and } V \text { and}\\ \text {energy between } E \text { and } E+\delta E. \end {array} \end{equation} The many states contributing to the count $\Omega (N,V,E)$ is called a microcanonical ensemble. Because one in principle can partition state space into non-overlapping sets of states where each set represents a unique value of $E$, $\Omega $ is also called the “microcanonical partition function”.

The Ising spin lattice is a simple statistical mechanical model with discrete energy levels which we can now introduce to gain some understanding of what it means to say $\Omega $ is “large.” Imagine a linear array of $N$ spins, each pointing either “up” or “down.”

Let us suppose that the Hamiltonian of this system is given by \begin{equation} \label {eq:ising_fieldhamil} \mathscr {H} = h\sum _{i=1}^{N}s_i \end{equation} where $s_i$ is -1 if spin $i$ is “down” and +1 if spin $i$ is “up,” and $h$ is some unit of energy. Let $m$ denote the number of spins that are up out of the $N$ spins, and let $\Omega (N,m)$ be the number of realizations of the $N$ spin states for which $m$ spins are up. The ground state, the state with the lowest energy, has all spins down, so $\Omega (N,0) = 1$. The next state up has one spin up, but there are $N$ possible microstates that have this energy: $\Omega (N,1) = N$. The next state up has two spins, and there are $N(N-1)/2$ such microstates: $\Omega (N,2) = N(N-1)/2$. For $m$ spins flipped, there are \begin{equation} \label {eq:choose} \Omega _m = \frac {N!}{(N-m)!m!} \end{equation} distinct microstates. You can now easily see that working with $\Omega $ for statistical mechanical systems means working with enormous numbers. For the rather small system of 100 spins, if we ask how many states are there with 50 spins up, we see $\Omega (100,50) \approx 10^{29}$.

Although quantum mechanics tells us that atomic systems have discrete energy levels, when systems contain very large numbers of atoms, these energy levels become so closely spaced relative to their span that they may effectively be considered a continuum. We can thus pass into a classical (as opposed to quantum mechanical) representation, where the microstate for a system of $N$ particles is specified by a point in a 6$N$-dimensional phase space: \begin{equation} \left (r^N,p^N\right ) \equiv \left (\rb _1,\rb _2,...,\rb _N;\pb _1,\pb _2,...,p_N\right ). \end{equation} where $\rb _i$ is the 3-space position of particle $i$ and $\pb _i$ is its momentum. We can denote the number of states in a classical microcanonical ensemble by integrating over all of phase space and plucking out those states for which the energy is $E$ using the Dirac delta function: \begin{align} \label {eq:microcanon_class} \overline \Omega (N,V,E) = & \frac {1}{h^{3N}N!}\int _V d\rb _1\int _Vd\rb _2\cdots \int _Vd\rb _N \int _{R^3}d\pb _1\int _{R^3}d\pb _2\cdots \int _{R^3}d\pb _N \delta \left [\mathscr {H}\left ({\bf r}^N,{\bf p}^N\right )-E\right ]\\ & \equiv \int \int d\rb ^Nd\pb ^N \delta \left [\mathscr {H}\left ({\bf r}^N,{\bf p}^N\right )-E\right ], \end{align}

(The second line introduces some shorthand notation.) The delta function integrand has units that are the reciprocal of its argument, so $\overline \Omega $ is more precisely termed the “density of states” and is directly related to $\Omega $: \begin{equation} \label {eq:dos} \overline \Omega \delta E = \Omega \end{equation} This definition satisfies the idea that the integral of $\overline \Omega $ over the entire continuous domain of energy $E$ should equal a complete integral over all of phase space.

The microcanonical ensemble represents a hyperdimensional surface in the phase space dimensioned by $N$ particles with positions limited by the extent of $V$. The factorial in Eq. 7, $N!$, takes into account that the particles are indistinguishable; that is, all orderings of particle indices 1, 2, $\dots $ N, are treated identically. $h$ is Planck’s constant; note that it has units of length$\bullet $momentum. Think of it as a quantum-mechanically-required “mesh discretization” for continuous space (it arises due to the Heisenberg uncertainty relation). It also nondimensionalizes the partition function. We will encounter it again in the next section, but we will also see why these “prefactors” are not essential ingredients of most molecular simulations.

You may wonder why there seem to be two viewpoints of statistical mechanics, quantum and classical. First, there really aren’t two viewpoints: the classical picture is an approximation of the more general quantum mechanical picture. But statistical mechanics as a discipline was first formalized by Gibbs and Boltzmann before quantum mechanics was widely accepted, so it dealt necessarily with systems of classical particles obeying Newtonian equations of motion; that is, on classical mechanics. There appears to be a general consensus that it is easier to introduce statistical mechanical concepts using the “sum-over-states” notation of quantum statistical mechanics, rather than the apparently more cumbersome (and anyway approximate) “integral-over-phase-space” notation of classical statistical mechanics.

2.2 Making Observations: The Ergodic Hypothesis

Scientists are taught early on that when conducting measurements, one must perform repeated experiments and average the results. If one makes $\mathscr {N}$ independent measurements of some observable $G$, one computes the mean value as \begin{equation} G_{\rm obs} = \frac {1}{\mathscr {N}}\sum _{i=1}^{\mathscr {N}}G_i. \end{equation} We now imagine that the time of a measurement is so short that we know that the system is in only one of its many possible microstates. This means we can write \begin{equation} \label {eq:exp_obs} G_{\rm obs} = \sum _{\nu }\left [\frac {1}{\mathscr {N}}\left (\begin {array}{c} \mbox {number of times state $\nu $ is} \\ \mbox {observed in the $\mathscr {N}$ observations} \end {array}\right )\right ]G_\nu , \end{equation} where $G_\nu $ was introduced previously (Eq. 1) as the value of observable $G$ when the system is in state $\nu $.

Now we have to imagine that our system is evolving in time. As it evolves, its degrees of freedom change values, and the system is thought of as tracing out a trajectory in state space. (“State space” is a Hilbert space spanned by all states $\left |\nu \right >$ in the quantum mechanical case, or phase space in the classical case.) How is the system evolving? The system wavefunction evolves according to Schrödinger’s equation, while particles in a classical system follow Newtonian mechanics. As the experimenters, we control the system by specifying a handful of variables, such as its total energy, $E$, the number of particles, $N$, and the volume, $V$. These constraints force the system’s trajectory to remain in a designated partition of state space.

The key assumption we make at this point is that, if we wait long enough, our system will visit every possible state; that is, the trajectory will eventually pass through every available point in state space consistent with our constraints (that is, all states in the partition). If this is true, and we make $\mathscr {N}$ independent observations, then the number of times we observe the system in state $\nu $ divided by the number of observations, $\mathscr {N}$, is the probability of observing state $\nu $, $P_\nu $, if we happen to make a random observation. So, Eq. 10 above becomes the ensemble average first presented in Eq. 1: \begin{equation} G_{\rm obs} = \sum _\nu P_\nu G_\nu = \left <G\right > \end{equation} This assumption is important: it is referred to as the ergodic hypothesis. A system is “ergodic” if, after a sufficiently long time, it visits all possible state space points consistent with whatever constraints are put on it. We cannot in general prove that any system is ergodic; it is something we are comfortable assuming for most systems based on our physical intuition. There are, however, many systems which are non-ergodic on the time-scales they are measured. The ergodic hypothesis is only strictly true when the number of measurements taken is sufficiently long that the computed probabilities $P_\nu $ no longer change, which can only be guaranteed for an infinite number of observations: \begin{equation} \left <G\right > = \lim _{\mathscr {N}\rightarrow \infty }\frac {1}{\mathscr {N}}\sum _{i=1}^\mathscr {N} G_{\nu (i)} \end{equation} where $\nu (i)$ refers to the state of the system at measurement $i$.

Another important consideration is the following: How far apart must the $\mathscr {N}$ independent measurements be from one another in time to be considered truly “independent”? To answer this question, we must introduce the notion of a relaxation time, $\tau _{\rm relax}$, which arises naturally due to the presumably chaotic nature of the microscopic system. Given some initial conditions, after a time $\tau _{\rm relax}$ has elapsed, the system has “lost memory” of the initial condition. We measure this loss of memory in terms of correlation functions, which will be discussed in more detail later. If we wait at least $\tau _{\rm relax}$ between successive observations, we know they are more likely to be independent. It turns out that one can use simulation methods to estimate relaxation times (and their spectra; many systems display a broad spectrum of relaxation times, each element cooresponding to a particular type of molecular motion). We will pay particularly close attention to $\tau _{\rm relax}$ in upcoming sections.

2.3 Entropy and Temperature

Eq. 3 introduced the quantity $\Omega (N,V,E)$ as the number of states available to a system under the constraints of constant number of particles, $N$, volume, $V$, and energy $E$. The fundamental postulate of statistical mechanics, also called the “rational basis”, is the following:

In statistical equilibrium, all states consistent with the constraints of $N$, $V$, and $E$ are equally probable.

or \begin{equation} P_\nu = 1/\Omega (N,V,E). \end{equation} This relation is often referred to as a statement of the “equal a priori probabilities in state space.” Another way of saying the same thing: The probability distribution for states in the microcanonical ensemble is uniform.

This postulate reflects the fact that we are maximally uncertain with regards to the probabilities of any particular arrangements of degrees of freedom inside a closed system. A closed system is one which cannot exchange energy, volume, or particles with the environment. As such, there is quite literally no way for us to learn anything at all about how particles are arranged, so we must assume all arrangments that satisfy the given energy, volume and number of particles are equiprobable.

One link between statistical mechanics and classical thermodynamics is given by a definition of entropy: \begin{equation} \label {eq:entropy} S \equiv k_B\ln \Omega \end{equation} Note two important properties of $S$. First, it is extensive: if we consider a compound system made of subsystems $A$ and $B$ with $\Omega _A$ and $\Omega _B$ as the respective number of states, the total number of states is $\Omega _A\Omega _B$, and therefore $S = S_A + S_B$. Second, it is consistent with the second law of thermodynamics: putting any constraint on the system lowers its entropy because the constraint lowers the number of accessible states.

Temperature is defined using entropy: $1/T \equiv \left (\partial S/\partial E\right )_{N,V}$, or \begin{equation} \beta = \left (k_BT\right )^{-1} = \left (\partial \ln \Omega /\partial E\right ) \end{equation} Now we will consider constraining our system not with constant $E$, but with constant $T$. The set of all possible states satisfying constraints of $N$, $V$, and $T$ is called the canonical ensemble. We now ask, what is the probability of any microstate in this ensemble? Consider a closed system divided into a small subsystem $A$ surrounded by a large “bath” $B$. We imagine that these two subsystems exchange only thermal energy, but no particles, and their volumes remain fixed. We seek to compute the probability of finding the total system in a state such that subsystem $A$ has energy $E_A$. The entire system is microcanonical, so the total energy, $E$ is constant, as is the total number of states available to the system, $\Omega (E)$ (we omit the $N$ and $V$ for simplicity).

When $A$ has energy $E_A$, the total system energy is $E = E_A + E_B$, where $E_B$ is the energy of the bath. By constraining system $A$’s energy, we have reduced the number of states available to the whole system to $\Omega (E-E_A)$. So, using the fundamental postulate, the probability of observing the closed system in a state in which subsystem $A$ has energy $E_A$ is \begin{equation} P_A = \frac {\left [\begin {array}{l} \mbox {Number of states for which}\\ \mbox {subsystem $A$ has energy $E_A$} \end {array}\right ]} {\left [\begin {array}{l} \mbox {Number of states available}\\ \mbox {to entire system} \end {array}\right ]} = \frac {\Omega \left (E - E_A\right )}{\Omega \left (E\right )} \end{equation} We can expand $\Omega \left (E - E_A\right )$ in a Taylor series around $E_A = 0$:

\begin{eqnarray} P_A \propto \Omega \left (E - E_A\right ) & = & \exp \left [\ln \Omega \left (E-E_A\right )\right ]\\ & = & \exp \left [\ln \Omega \left (E\right ) - E_A\frac {\partial \ln \Omega }{\partial E} + \cdots \right ], \end{eqnarray}

where the partial derivative implies we are holding $N$ and $V$ fixed. We can truncate the Taylor expansion at the first-order term, because higher order terms become less and less important as the size of subsystem $B$ becomes larger and larger. What results is the Boltzmann distribution law for energies of a system at constant temperature: \begin{equation} P_A \propto \exp \left (-\beta E_A\right ) \end{equation} The normalization condition requires that for all energies of subsystem $A$, $E_A$, \begin{equation} \sum _A P_A = 1 = \sum _A \exp \left (-\beta E_A\right ) \equiv Q\left (N,V,T\right ), \end{equation} which defines the canonical partition function, $Q$. Therefore, \begin{equation} P_A = Q^{-1}\exp \left (-\beta E_A\right ) \end{equation} Because some energies can correspond to more than one microstate, we should distinguish between “states” and “energy levels.” We can express the canonical partition function as \begin{equation} \label {eq:canonical_1} Q = \sum _{\begin {array}{cc}\nu \\\mbox {states}\end {array}} e^{-\beta E_\nu } = \sum _{\begin {array}{cc}l\\\mbox {levels}\end {array}} \Omega \left (E_l\right ) e^{-\beta E_l} \end{equation} where, as we have seen, $\Omega \left (E\right )$ is the number of microstates with energy $E$. Moving to the continuum limit, and assuming a reference energy of $E_{ref} = 0$, \begin{equation} \label {eq:dos2} Q \rightarrow \int _0^\infty dE \overline \Omega \left (E\right )e^{-\beta E} \end{equation} where $\overline \Omega \left (E\right )$ is the density of states. What is this equation telling us? It is telling us that $Q$ is the Laplace transform of $\overline \Omega $. We know that transform pairs are unique, and hence, both $Q$ and $\overline \Omega $ contain the same information.

We recognize that for a system described by a canonical ensemble, the energy is a fluctuating quantity. And we now have the probability of observing a state with a given energy, so we can use Eq. 1 to compute the average energy, $\left <E\right >$. Consider

\begin{eqnarray} \left <E\right > & = &\left <E_\nu \right > = \sum _\nu P_\nu E_\nu \\ & = & \left [\sum _\nu E_\nu \exp \left (-\beta E_\nu \right )\right ] \left / \left [\sum _\nu \exp \left (-\beta E_\nu \right )\right ] \right . \end{eqnarray}

Notice that \begin{equation} \sum _\nu E_\nu \exp \left (-\beta E_\nu \right ) = -\left (\partial Q \left / \partial \beta \right .\right ) \end{equation} Recalling that $d\ln f(x) / dx = 1/f df/dx$, we see that

\begin{eqnarray} \left <E\right > & = & -\left (\partial Q \left / \partial \beta \right .\right )\left /Q\right .\\ \label {eq:ave_ener_canon} & = & -\left (\partial \ln Q \left / \partial \beta \right .\right )_{N,V} \end{eqnarray}

Now, let us consider the average magnitude of the fluctuations in energy in the canonical ensemble.

\begin{eqnarray} \left <\left (\delta E\right )^2\right > & = & \left <\left (E - \left <E\right >\right )^2\right >\\ & = & \left <E^2\right > - \left <E\right >^2\\ & = & \sum _\nu P_\nu E_\nu ^2 - \left (\sum _\nu P_\nu E_\nu \right )^2\\ & = & Q^{-1}\left (\frac {\partial ^2Q}{\partial \beta ^2}\right )_{N,V} - Q^{-2}\left (\frac {\partial Q}{\partial \beta }\right )^2_{N,V}\\ & = & \left (\frac {\partial \ln Q}{\partial \beta ^2}\right )_{N,V} = -\left (\frac {\partial \left <E\right >}{\partial \beta }\right )_{N,V} \end{eqnarray}

Now, noting that the definition of heat capacity at constant volume, $C_v$, is \begin{equation} C_v = \left (\frac {\partial E}{\partial T}\right ) \end{equation} we see that \begin{equation} \left <\left (\delta E\right )^2\right > = k_B T^2 C_v \end{equation} This is an interesting statement. It relates the magnitude of spontaneous fluctuations in the total energy of a system to that system’s capacity to store or release energy due to changing its temperature.

The fact (Eq. ??) that the average energy in the canonical ensemble is related to a derivative of the log of the partition function implies that $\ln Q$ is an important thermodynamic quantity. So, let’s go back to our undergraduate thermodynamics course(s) and recall the following statement of the 1st and 2nd Law: \begin{equation} dA = -SdT - pdV + \mu dN \end{equation} where $A$ is the Helmholtz free energy, defined in terms of internal energy and entropy as \begin{equation} \label {eq:helmholtz} A = \left <E\right > - TS \end{equation} Now, consider the following derivative of $A$:

\begin{eqnarray} \left (\frac {\partial \left (A/T\right )}{\partial \left (1/T\right )}\right )_{N,V} & = & A + \frac {1}{T}\left (\frac {\partial A}{\partial \left (1/T\right )}\right )_{N,V}\\ & = & A - T\left (\frac {\partial A}{\partial T}\right )_{N,V} = A + TS = \left <E\right >. \end{eqnarray}

Therefore, \begin{equation} \label {eq:dAdBeta_E} \left (\frac {\partial \left (\beta A\right )}{\partial \beta }\right )_{N,V} = \left <E\right >. \end{equation} Considering Eq. ??, we see that \begin{equation} \label {eq:C} \ln Q + C = -\beta A \end{equation} which does indeed suggest an important link between $\ln Q$ and the important thermodynamic quantity, the Helmholtz free energy. But what is the constant $C$? To evaluate it, consider the “boundary condition” as $T \rightarrow 0$: \begin{equation} Q = \sum _\nu e^{-\beta E_\nu } \stackrel {T\rightarrow 0}{\longrightarrow } e^{-\beta E_{\rm ground}}. \end{equation} Here, we have assumed that the degeneracy of the ground state, $\Omega \left (E_{\rm ground}\right )$ is 1. This tells us that \begin{equation} \lim _{T\rightarrow 0} \ln Q = -\beta E_{\rm ground} \end{equation} Using this fact, and combining Eqs. 26 and 28, as $T \rightarrow 0$, we see that \begin{equation} -\beta E_{\rm ground} + C = -\beta \underbrace {\left <E\right >}_{E_{\rm ground}} - \underbrace {\frac {S}{k_B}}_{\rightarrow 0\ (\Omega = 1)}, \end{equation} Hence, $C$ = 0. So, \begin{equation} \ln Q = -\beta A. \end{equation} The quantity $-\beta ^{-1}\ln Q$ is the Helmholtz free energy, $A$.

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{eqnarray} \label {eq:q_fact} 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{eqnarray}

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{eqnarray} \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{eqnarray}

(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{eqnarray} \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{eqnarray}

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. 39 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.

← Previous Contents Next →