Derivation of Gibbs Ensembles

Classical Hamiltonian Mechanics

We start with Hamiltons equations.

\dot{\vec{q}} = \frac{\partial H (\vec{q}, \vec{p})}{\partial \vec{p}}

\dot{\vec{p}} = - \frac{\partial H (\vec{q}, \vec{p})}{\partial \vec{q}}

Classical Microcanonical Ensemble

We start with some probability distribution in the phase space of the Hamiltonian $H (\vec{q}, \vec{p})$ .

ρ (\vec{q}, \vec{p}, t) = ? ? ?

The probability distribution evolves by the continuity equation.

\frac{\partial ρ (\vec{q}, \vec{p}, t)}{\partial t} = \nabla \cdot (ρ \vec{v})

Here $\vec{v} = ⟨ {\dot{q}}_{1}, {\dot{q}}_{2}, \dots, {\dot{q}}_{N}, {\dot{p}}_{1}, {\dot{p}}_{2}, \dots, {\dot{p}}_{N} ⟩$

In equilibrium the probability distribution at every point in phase space shouldn't change.

\nabla \cdot (ρ \vec{v}) = 0

\sum_{i = 1}^{2 N} \frac{\partial ρ (\vec{q}, \vec{p}) {\vec{v}}_{i}}{\partial a_{i}} = 0

\sum_{i = 1}^{N} \frac{\partial ρ (\vec{q}, \vec{p}) {\dot{q}}_{i}}{\partial q_{i}} + \frac{\partial ρ (\vec{q}, \vec{p}) {\dot{p}}_{i}}{\partial p_{i}} = 0

Simple chain rule.

\sum_{i = 1}^{N} ρ (\vec{q}, \vec{p}) \frac{\partial {\dot{q}}_{i}}{\partial q_{i}} + {\dot{q}}_{i} \frac{\partial ρ (\vec{q}, \vec{p})}{\partial q_{i}} + ρ (\vec{q}, \vec{p}) \frac{\partial {\dot{p}}_{i}}{\partial p_{i}} + {\dot{p}}_{i} \frac{\partial ρ (\vec{q}, \vec{p})}{\partial p_{i}} = 0

Substitute in Hamiltons equations.

\sum_{i = 1}^{N} ρ (\vec{q}, \vec{p}) \frac{\partial H (\vec{q}, \vec{p})}{\partial q_{i} \partial p_{i}} + {\dot{q}}_{i} \frac{\partial ρ (\vec{q}, \vec{p})}{\partial q_{i}} - ρ (\vec{q}, \vec{p}) \frac{\partial H (\vec{q}, \vec{p})}{\partial q_{i} \partial p_{i}} + {\dot{p}}_{i} \frac{\partial ρ (\vec{q}, \vec{p})}{\partial p_{i}} = 0

\sum_{i = 1}^{N} {\dot{q}}_{i} \frac{\partial ρ (\vec{q}, \vec{p})}{\partial q_{i}} + {\dot{p}}_{i} \frac{\partial ρ (\vec{q}, \vec{p})}{\partial p_{i}} = 0

\sum_{i = 1}^{N} \frac{\partial H (\vec{q}, \vec{p})}{\partial p_{i}} \frac{\partial ρ (\vec{q}, \vec{p})}{\partial q_{i}} - \frac{\partial H (\vec{q}, \vec{p})}{\partial q_{i}} \frac{\partial ρ (\vec{q}, \vec{p})}{\partial p_{i}} = 0

Next, we say that the probability distribution is a function only of the total energy, so $ρ (\vec{q}, \vec{p}, t_{0}) = f (H (\vec{q}, \vec{p}))$ . This means that $\frac{\partial ρ (\vec{q}, \vec{p}, t_{0})}{\partial a_{i}} = \frac{\partial f (H (\vec{q}, \vec{p}))}{\partial H} \frac{\partial H (\vec{q}, \vec{p})}{\partial a_{i}}$

\sum_{i = 1}^{N} \frac{\partial f (H)}{\partial H} (\frac{\partial H (\vec{q}, \vec{p})}{\partial q_{i}} \frac{\partial H (\vec{q}, \vec{p})}{\partial p_{i}} - \frac{\partial H (\vec{q}, \vec{p})}{\partial q_{i}} \frac{\partial H (\vec{q}, \vec{p})}{\partial p_{i}}) = 0

That assumption creates a solution to the continuity equation. Now in the microcanonical ensemble (NVE) we say the total energy is known, so the only possibility for $f$ is a uniform distribution over the selected energy level.

ρ (\vec{q}, \vec{p}) = \frac{δ (H (\vec{q}, \vec{p}) - E)}{\int_{V} δ (H ({\vec{q}}^{'}, {\vec{p}}^{'}) - E) d {\vec{q}}^{'} d {\vec{p}}^{'}}

This equation is the principle of equal a priori for the equilibrium in a closed system.

Classical Canonical Ensemble

Define the number of states accessable by a system with $N$ degrees of freedom at energy $E$ to be the below.

Ω (E) = \frac{1}{h^{N}} \int_{V} δ (H (\vec{q}, \vec{p}) - E) d \vec{q} d \vec{p}

Next split the degrees of freedom into the "bath" and "system". Let the degrees of freedom for the bath be $r_{bath} = ⟨ q_{M + 1}, \dots, q_{N}, p_{M + 1}, \dots, p_{N} ⟩$ and the degrees of freedom for the system $r_{sys} = ⟨ q_{1}, \dots, q_{M},, p_{1}, \dots, p_{M} ⟩$ .

In the Canonical ensemble the combined energy of the bath and system should be constant. The combination of the two makes a microcanonical ensemble.

H_{total} (r_{bath}, r_{sys}) = H_{bath} (r_{bath}) + H_{sys} (r_{sys}) + H_{coupling} (r_{bath}, r_{sys})

The total accessable states with an energy should be the total number combined number of states the bath and system can be in that sum to that energy.

Ω_{total} (E_{total}) = \frac{1}{h^{N}} \int_{V} Ω_{bath} (E_{total} - H_{sys} (r_{sys})) d r_{sys}

And in the microcanonical ensemble all microstates $(r_{bath}, r_{sys})$ with the correct energy have an equal probability density.

ρ (r_{bath}, r_{sys}) = \frac{1}{Ω_{total} (E_{total})} \frac{1}{h^{N}} δ (E_{total} - (H_{sys} (r_{sys}) + H_{bath} (r_{bath})))

We can calculate the probability density of just the microstate for the system.

ρ (r_{sys}) = \int_{V} ρ (r_{bath}, r_{sys}) d r_{bath}

= \int_{V} \frac{1}{Ω_{total} (E_{total})} \frac{1}{h^{N}} δ (E_{total} - (H_{sys} (r_{sys}) + H_{bath} (r_{bath}))) d r_{bath}

= \frac{Ω_{bath} (E_{total} - H_{sys} (r_{sys}))}{h^{M} Ω_{total} (E_{total})}

Take the natural log of both sides

\ln (ρ (r_{sys})) = \ln (\frac{Ω_{bath} (E_{total} - H_{sys} (r_{sys}))}{h^{M} Ω_{total} (E_{total})})

\ln (ρ (r_{sys})) = \ln (Ω_{bath} (E_{total} - H_{sys} (r_{sys}))) - \ln (h^{M} Ω_{total} (E_{total}))

Then we expand the logarithm to its Taylor series.

f (x + c) = f (x) + \frac{d f (x)}{d x} c + \frac{1}{2} \frac{d^{2} f (x)}{d x^{2}} c^{2} + \frac{1}{6} \frac{d^{3} f (x)}{d x^{3}} c^{3} + \dots

We set $f (x) = \ln (Ω_{bath} (x))$ , $x = E_{total}$ , and $c = - H_{sys} (r_{sys})$ we get the following expression.

\ln (ρ (r_{sys})) = \ln (Ω_{bath} (E_{total})) - \frac{d \ln (Ω_{bath} (E_{total}))}{d E_{total}} H_{sys} (r_{sys})

+ \frac{1}{2} \frac{d^{2} \ln (Ω_{bath} (E_{total}))}{d E_{total}^{2}} H_{sys} (r_{sys})^{2} - \frac{1}{6} \frac{d^{3} \ln (Ω_{bath} (E_{total}))}{d E_{total}^{3}} H_{sys} (r_{sys})^{3} + \dots

\dots - \ln (h^{M} Ω_{total} (E_{total}))

By definition if the bath was a microcanonical ensemble the entropy would be $S_{bath} (E_{total}) = k_{B} \ln (Ω_{bath} (E_{total}))$

= \ln (Ω_{bath} (E_{total})) - \frac{1}{k_{B}} \frac{d S_{bath} (E_{total})}{d E} H_{sys} (r_{sys}) + \frac{1}{2 k_{B}} \frac{d^{2} S_{bath} (E_{total})}{d E^{2}} H_{sys} (r_{sys})^{2}

- \frac{1}{6 k_{B}} \frac{d^{3} S_{bath} (E_{total})}{d E^{3}} H_{sys} (r_{sys})^{3} + \dots - \ln (h^{M} Ω_{total} (E_{total}))

And by definition $\frac{1}{T_{bath} (E)} = \frac{d S_{bath} (E)}{d E}$ .

= \ln (Ω_{bath} (E_{total})) - \frac{1}{k_{B} T_{bath} (E_{total})} H_{sys} (r_{sys}) + \frac{1}{2 k_{B} T_{bath} (E_{total})} \frac{d T_{bath} (E_{total})}{d E^{2}} H_{sys} (r_{sys})^{2}

- \frac{1}{6 k_{B} T_{bath} (E_{total})} \frac{d^{2} T_{bath} (E_{total})}{d E^{2}} H_{sys} (r_{sys})^{3} + \dots - \ln (h^{M} Ω_{total} (E_{total}))

Then if we assume the heat capacity of the bath is infinite so that $\frac{d T_{bath}}{d E} = 0$ we can remove the extra terms, and the temperature just becomes a constant so $T_{bath} = T_{bath} (E_{total}) = T_{bath} (E_{bath})$ .

\ln (ρ (r_{sys})) = - \frac{1}{k_{B} T_{bath}} H_{sys} (r_{sys}) + \ln (\frac{Ω_{bath} (E_{bath})}{h^{M} Ω_{total} (E_{total})})

Then exponentiate both sides.

ρ (r_{sys}) = \frac{e^{- \frac{H_{sys} (r_{sys})}{k_{B} T_{bath}}}}{(\frac{h^{M} Ω_{total} (E_{total})}{Ω_{bath} (E_{bath})})}

And the denominator is just a normalizing factor so the final expression is below.

ρ (r_{sys}) = \frac{e^{- \frac{H_{sys} (r_{sys})}{k_{B} T_{bath}}}}{\int_{V} e^{- \frac{H_{sys} (r_{sys}^{'})}{k_{B} T_{bath}}} d r_{sys}^{'}}

Quantum Microcanonical Ensemble

We assume a time evolution of a system of the Schrodinger equation for some Hamiltonian operator.

\frac{\partial Ψ (\vec{x}, t)}{\partial t} = - \frac{i}{ℏ} \hat{H} Ψ (\vec{x}, t)

To measure an observable when the stat is known is

A = ⟨ Ψ | \hat{A} | Ψ ⟩

but if there is a probability distrubition over some states then it is below.

= ⟨ A ⟩ = \sum_{i} p_{i} ⟨ Ψ_{i} | \hat{A} | Ψ_{i} ⟩

If we have an orthanormal basis $n$ so that $\sum_{k} | n_{k} ⟩ ⟨ n_{k} | = I$ we can insert it into the equation.

= \sum_{i} p_{i} ⟨ Ψ_{i} | \hat{A} (\sum_{k} | n_{k} ⟩ ⟨ n_{k} |) | Ψ_{i} ⟩

= \sum_{k} \sum_{i} p_{i} ⟨ Ψ_{i} | \hat{A} | n_{k} ⟩ ⟨ n_{k} | Ψ_{i} ⟩

= \sum_{k} \sum_{i} p_{i} ⟨ n_{k} | Ψ_{i} ⟩ ⟨ Ψ_{i} | \hat{A} | n_{k} ⟩

= \sum_{k} ⟨ n_{k} | (\sum_{i} p_{i} | Ψ_{i} ⟩ ⟨ Ψ_{i} | \hat{A}) | n_{k} ⟩

= Tr (\sum_{i} p_{i} | Ψ_{i} ⟩ ⟨ Ψ_{i} | \hat{A})

If we define an operator

ρ = \sum_{i} p_{i} | Ψ_{i} ⟩ ⟨ Ψ_{i} |

then the expected value of any observable is

⟨ A ⟩ = Tr (ρ \hat{A})

Next we find the time evolution of $ρ$ .

\frac{d ρ}{d t} = \sum_{i} p_{i} \frac{d (| Ψ_{i} ⟩ ⟨ Ψ_{i} |)}{d t} = \sum_{i} p_{i} (\frac{d | Ψ_{i} ⟩}{d t} ⟨ Ψ_{i} | + | Ψ_{i} ⟩ \frac{d ⟨ Ψ_{i} |}{d t})

Then substitute the Schrodinger equation

= \sum_{i} p_{i} (- \frac{i}{ℏ} \hat{H} | Ψ_{i} ⟩ ⟨ Ψ_{i} | + | Ψ_{i} ⟩ ⟨ Ψ_{i} | \frac{i}{ℏ} \hat{H})

= - \frac{i}{ℏ} \sum_{i} p_{i} (\hat{H} | Ψ_{i} ⟩ ⟨ Ψ_{i} | - | Ψ_{i} ⟩ ⟨ Ψ_{i} | \hat{H})

= - \frac{i}{ℏ} (\hat{H} ρ - ρ \hat{H})

So it must be that for equilibrium they commute

\hat{H} ρ = ρ \hat{H}

Similar to the classical case, if we have the $p_{i}$ be only a function of the energy for each eigenvector of the Hamiltonian, the system will be in equilibrium. Since observables are Hermetian, an orthogonal eigenbases must exist. Each operator will be defined in the eigenbasis.

\hat{H} = \sum_{a} E_{a} | n_{a} ⟩ ⟨ n_{a} |

ρ = \sum_{a} p (E_{a}) | n_{a} ⟩ ⟨ n_{a} |

We can show they commute

\hat{H} ρ = \sum_{a} E_{a} | n_{a} ⟩ ⟨ n_{a} | \sum_{b} p (E_{b}) | n_{b} ⟩ ⟨ n_{b} |

= \sum_{a} \sum_{b} E_{a} p (E_{b}) | n_{a} ⟩ ⟨ n_{a} | | n_{b} ⟩ ⟨ n_{b} |

= \sum_{a} \sum_{b} p (E_{b}) | n_{a} ⟩ ⟨ n_{a} | E_{a} | n_{b} ⟩ ⟨ n_{b} |

= \sum_{a} p (E_{a}) | n_{a} ⟩ ⟨ n_{a} | \sum_{b} E_{b} | n_{b} ⟩ ⟨ n_{b} |

ρ \hat{H}

So the quantum microcanonical ensemble must have a density matrix of the form

ρ = \sum_{a} p (E_{a}) | n_{a} ⟩ ⟨ n_{a} |

where the $n$ 's are eigenvectors of the Hamiltonian. Since the energy is known, it must be a uniform distribution over all eigenstates with that energy

ρ = \frac{1}{Ω_{total} (E)} \sum_{a = 1}^{Ω_{total} (E)} | n_{a}, E ⟩ ⟨ n_{a}, E |

Quantum Cannonical Ensemble

We assume the total Hamiltonian is the form

\hat{H} = {\hat{H}}_{sys} \otimes {\hat{0}}_{bath} + {\hat{0}}_{sys} \otimes {\hat{H}}_{bath} + {\hat{H}}_{coupling}

And we assume the coupling energy is none. From the microcannonical ensemble the total density matrix is

ρ_{total} = \frac{1}{Ω_{total} (E_{total})} \sum_{a = 1}^{Ω_{total} (E_{total})} | n_{a}, E_{total} ⟩ ⟨ n_{a}, E_{total} |

And the density matrix of just the system is the partial trace

ρ_{sys} = \sum_{b} \sum_{s_{1}} \sum_{s_{2}} | n_{s_{1}} ⟩ ⟨ n_{s_{2}} | (⟨ n_{s_{1}} \otimes n_{b} | ρ_{total} | n_{s_{2}} \otimes n_{b} ⟩)

Where $n_{b}, n_{s}$ are the eigenbasis of $H_{bath}, H_{sys}$ . Measuing a system microstate gives.

⟨ ψ | ρ_{sys} | ψ ⟩ = {Tr}_{bath} (ρ_{total}) = \frac{1}{Ω_{total} (E_{total})} \sum_{b} \sum_{s_{1}} \sum_{s_{2}} ⟨ ψ | n_{s_{1}} ⟩ ⟨ n_{s_{2}} | ψ ⟩ (⟨ n_{s_{1}} \otimes n_{b} | ρ_{total} | n_{s_{2}} \otimes n_{b} ⟩)

= \frac{1}{Ω_{total} (E_{total})} \sum_{b} ⟨ ψ \otimes n_{b} | ρ_{total} | ψ \otimes n_{b} ⟩

And the values of $ρ_{total}$ are only the eigenvectors of eigenvalue $E_{total}$ , so it would just be the ones where the bath has the complement energy.

= \frac{Ω_{bath} (E_{total} - E_{sys})}{Ω_{total} (E_{total})}

Take the ln of both sides

\ln (⟨ ψ | ρ_{sys} | ψ ⟩) = \ln (\frac{Ω_{bath} (E_{total} - ⟨ ψ | H_{sys} | ψ ⟩)}{Ω_{total} (E_{total})}) = \ln (Ω_{bath} (E_{total} - ⟨ ψ | H_{sys} | ψ ⟩)) - \ln (Ω_{total} (E_{total}))

Then Taylor series again

= \ln (Ω_{bath} (E_{total})) - \frac{d \ln (Ω_{bath} (E_{total}))}{d E} (⟨ ψ | H_{sys} | ψ ⟩) + \frac{d^{2} \ln (Ω_{bath} (E_{total}))}{d E^{2}} (⟨ ψ | H_{sys} | ψ ⟩)^{2} + \dots - \ln (Ω_{total} (E_{total}))

We are assuming that the number of states is dense enough that the derivative exists. Substitute the definition of entropy

= \ln (Ω_{bath} (E_{total})) - \frac{1}{k_{B}} \frac{d S_{bath} (E_{total})}{d E} (⟨ ψ | H_{sys} | ψ ⟩) + \frac{d^{2} S_{bath} (E_{total})}{d E^{2}} (⟨ ψ | H_{sys} | ψ ⟩)^{2} + \dots - \ln (Ω_{total} (E_{total}))

Then the definition of temperature

= \ln (Ω_{bath} (E_{total})) - \frac{1}{k_{B} T (E_{total})} (⟨ ψ | H_{sys} | ψ ⟩) + \frac{d^{2} T_{bath}^{- 1} (E_{total})}{d E^{2}} (⟨ ψ | H_{sys} | ψ ⟩)^{2} + \dots - \ln (Ω_{total} (E_{total}))

Assume heat capacity is infinte

= \ln (Ω_{bath} (E_{total})) - \frac{1}{k_{B} T} (⟨ ψ | H_{sys} | ψ ⟩) - \ln (Ω_{total} (E_{total}))

then exponentiate both sides to get

⟨ ψ | ρ_{sys} | ψ ⟩ = e^{- \frac{1}{k_{B} T} (⟨ ψ | H_{sys} | ψ ⟩)} \frac{Ω_{bath} (E_{total})}{Ω_{total} (E_{total})}

and the extra term is just a normalizing constant so

⟨ ψ | ρ_{sys} | ψ ⟩ = \frac{1}{Tr (H_{sys})} e^{- \frac{1}{k_{B} T} (⟨ ψ | H_{sys} | ψ ⟩)}

And since $ψ$ is any eigenstate it is just the operator exponential.

ρ_{sys} = \frac{1}{Tr (H_{sys})} e^{- \frac{1}{k_{B} T} H_{sys}}