Thermal Density Operators

For a quantum system with Hamiltonian $H$ in equilibrium at absolute temperature $T$, the canonical thermal state is the Gibbs density operator

$$\rho_\beta = \frac{e^{-\beta H}}{Z(\beta)}, \qquad Z(\beta) = \operatorname{Tr} e^{-\beta H},$$

where

$$\beta = \frac{1}{k_{\mathrm B}T}.$$

The partition function $Z$ normalizes the operator:

$$\operatorname{Tr}\rho_\beta=1.$$

The central task is to understand this as an operator: its eigenvalues are probabilities, degeneracies count distinct states, and its normalization must exist. For fixed-number thermodynamic derivatives see the canonical ensemble; for particle exchange see the grand-canonical ensemble.

Helpful background. Density operators provide positivity and trace normalization; the ensemble comparison explains the choice of trace space.

Gibbs operators and trace convergence

In a finite-dimensional Hilbert space, $e^{-\beta H}$ is well defined for every self-adjoint $H$ and every finite real $\beta$. For the usual positive-temperature state, $\beta>0$.

In an infinite-dimensional Hilbert space, more is required:

1. $H$ must be self-adjoint.
2. The exponential $e^{-\beta H}$ must be trace class.
3. The partition function must satisfy

$$0<Z(\beta)<\infty.$$

A Hamiltonian bounded below with a sufficiently fast-growing discrete spectrum often satisfies these conditions for $\beta>0$. They can fail for an unconfined continuum system in infinite volume or for spectra with excessive degeneracy.

When $Z$ diverges, the formal ratio is not a density operator. One can first normalize in finite volume and then study limits of local observables, or use an infinite-system equilibrium formalism. Dividing a divergent operator trace by volume does not itself construct a normalized density operator on the original infinite-volume Hilbert space.

Operator Exponential

The notation $e^{-\beta H}$ is defined by functional calculus. In finite dimension it can be represented by the convergent series

$$e^{-\beta H} = \sum_{m=0}^{\infty} \frac{(-\beta H)^m}{m!}.$$

The spectral theorem gives the more useful form. If

$$H = \sum_a E_a P_a,$$

where $P_a$ projects onto the eigenspace with energy $E_a$, then

$$e^{-\beta H} = \sum_a e^{-\beta E_a}P_a.$$

Therefore

$$\rho_\beta = \frac{1}{Z} \sum_a e^{-\beta E_a}P_a.$$

The formula is basis independent. The energy eigenbasis merely makes its probabilities visible.

Energy-Basis Representation

For a nondegenerate discrete spectrum,

$$H|E_n\rangle = E_n|E_n\rangle,$$

the thermal state is

$$\rho_\beta = \sum_n p_n |E_n\rangle\langle E_n|,$$

with

$$p_n = \frac{e^{-\beta E_n}}{Z}, \qquad Z = \sum_n e^{-\beta E_n}.$$

Thus an energy measurement returns $E_n$ with Boltzmann probability $p_n$.

The state is generally mixed:

$$\operatorname{Tr}\rho_\beta^2 < 1$$

whenever more than one energy eigenstate has nonzero weight. At finite positive temperature in finite dimension, every energy eigenstate has nonzero weight and $\rho_\beta$ has full rank.

Degeneracy

Let the energy $E_a$ have degeneracy

$$g_a = \operatorname{Tr}P_a.$$

Every normalized vector in that eigenspace receives the same weight per basis state:

$$\rho_\beta P_a = \frac{e^{-\beta E_a}}{Z} P_a.$$

The probability of measuring the energy $E_a$ is

$$\Pr(E_a) = \frac{ g_a e^{-\beta E_a} }{ Z }.$$

The factor $g_a$ belongs to the probability of the energy value, not to the weight of each orthonormal microstate.

Because $\rho_\beta$ is proportional to the identity within a degenerate eigenspace, it does not choose a preferred basis there. A basis rotation within that eigenspace cannot create physical thermal coherence.

Positivity and stationarity

The positive eigenvalues $e^{-\beta E_n}/Z$ show that $\rho_\beta$ is Hermitian and positive. Normalization gives unit trace. Since it is a function of $H$, $[H,\rho_\beta]=0$, so it is stationary under the same time-independent Hamiltonian.

The converse fails: an arbitrary mixture of energy eigenstates is stationary without having Boltzmann weights. Even coherences inside an exactly degenerate eigenspace can be stationary, whereas the Gibbs state is proportional to the identity in that space.

Invariance Under an Energy-Zero Shift

Replace the Hamiltonian by

$$H' = H+C\mathbb I.$$

Then

$$\begin{aligned} e^{-\beta H'} &= e^{-\beta C} e^{-\beta H}, \\ Z' &= e^{-\beta C}Z. \end{aligned}$$

The scalar factor cancels:

$$\rho_\beta' = \rho_\beta.$$

Thermal probabilities do not depend on the arbitrary zero of energy. Thermodynamic potentials such as $F=-k_{\mathrm B}T\ln Z$ do shift by the same constant $C$, so one must distinguish state invariance from potential conventions.

Thermal Expectation Values

For any observable $A$,

$$\langle A\rangle_\beta = \operatorname{Tr} \left( \rho_\beta A \right).$$

In a discrete energy basis,

$$\langle A\rangle_\beta = \frac{1}{Z} \sum_n e^{-\beta E_n} \langle E_n|A|E_n\rangle.$$

Only the energy-basis diagonal matrix elements of $A$ enter this one-time expectation, even when $[A,H]\neq0$. Off-diagonal matrix elements remain essential for time-dependent correlations and response functions.

For a function $f(H)$,

$$\langle f(H)\rangle_\beta = \frac{1}{Z} \sum_n f(E_n)e^{-\beta E_n}.$$

The quantum trace has reduced the calculation to a weighted spectral sum.

Entropy of a Gibbs State

The von Neumann entropy with thermodynamic units is

$$S(\rho) = -k_{\mathrm B} \operatorname{Tr} \left( \rho\ln\rho \right).$$

For the Gibbs state,

$$\ln\rho_\beta = -\beta H-\ln Z,$$

so

$$S(\rho_\beta) = k_{\mathrm B} \left( \beta U+\ln Z \right).$$

Using $T=(k_{\mathrm B}\beta)^{-1}$ gives

$$F = U-TS = -k_{\mathrm B}T\ln Z.$$

This entropy is the entropy of the full thermal mixed state. It should not be confused with the entanglement entropy of a subsystem in a pure state, though reduced thermal states can also carry entanglement and classical correlations.

Relative-Entropy Characterization

The Gibbs state also minimizes the nonequilibrium free-energy functional

$$\mathcal F_T(\rho) = \operatorname{Tr}(\rho H) - T S(\rho)$$

at fixed $H$ and positive $T$. Using quantum relative entropy,

$$D(\rho\Vert\rho_\beta) = \operatorname{Tr} \left[ \rho \left( \ln\rho-\ln\rho_\beta \right) \right],$$

one finds

$$\mathcal F_T(\rho) - \mathcal F_T(\rho_\beta) = k_{\mathrm B}T D(\rho\Vert\rho_\beta) \geq 0.$$

For a normalizable negative-temperature state the identity still holds, but its sign reverses: the Gibbs state maximizes this free-energy functional. At either sign, fixing the same mean energy makes it the entropy maximum. See maximum entropy for that constrained proof and quantum thermodynamics for operational work questions.

High-Temperature Limit

For a finite $D$-dimensional Hilbert space,

$$\lim_{\beta\to0} \rho_\beta = \frac{\mathbb I_D}{D}.$$

The state becomes maximally mixed. To first order,

$$\rho_\beta ={} \frac{1}{D} \left[ \mathbb I_D - \beta \left( H-\overline E\,\mathbb I_D \right) \right] + O(\beta^2),$$

where

$$\overline E = \frac{\operatorname{Tr}H}{D}.$$

The subtraction of $\overline E$ preserves normalization at first order.

In an infinite-dimensional space, the formal maximally mixed state $\mathbb I/D$ does not exist. The limit $\beta\to0^+$ can make $Z$ diverge, so high-temperature statements must be made for regulated systems or appropriate observables.

Zero-Temperature Limit

Assume a finite-dimensional system, or a discrete Hamiltonian with a finite ground-space degeneracy and a finite Gibbs trace at some positive inverse temperature. Let the ground energy be $E_0$ and the ground-space projector be $P_0$, with $g_0=\operatorname{Tr}P_0$.

As $\beta\to\infty$,

$$\rho_\beta \longrightarrow \frac{P_0}{g_0}.$$

For a nondegenerate ground state,

$$\rho_\beta \longrightarrow |E_0\rangle\langle E_0|.$$

For a degenerate ground space, the Gibbs limit without additional fields or boundary conditions is the equal mixture on that space, not a selected pure ground state.

If the first excited energy is separated by a gap

$$\Delta = E_1-E_0>0,$$

low-temperature corrections are suppressed by factors of order

$$e^{-\beta\Delta}.$$

In a thermodynamic limit the gap can close, degeneracies can grow, and the order of $\beta\to\infty$, volume, and source limits can matter.

Negative Temperatures

The formula permits $\beta<0$ only when

$$Z(\beta) = \operatorname{Tr} e^{-\beta H}$$

remains finite. An upper spectral bound is necessary, but is not sufficient in infinite dimension: the sum over all states must still converge. Finite-dimensional spin systems automatically satisfy this trace condition.

For $\beta<0$, higher-energy states receive larger Boltzmann weights. Such a state is hotter than every positive-temperature Gibbs state in the thermodynamic sense associated with population inversion; it is not colder than zero.

An unbounded-above Hamiltonian such as the harmonic oscillator does not admit a normalizable negative-temperature Gibbs state. Population inversion in a subsystem is not by itself enough; the equilibrium constraints and bounded spectrum must be specified.

Imaginary-Time Interpretation

Real-time evolution is generated by

$$U(t) = e^{-iHt/\hbar}.$$

After the formal substitution

$$t=-i\tau,$$

the propagator becomes

$$e^{-H\tau/\hbar}.$$

Setting

$$\tau = \beta\hbar$$

produces the unnormalized Gibbs operator $e^{-\beta H}$. The partition function is its trace:

$$Z(\beta) = \operatorname{Tr} e^{-H\tau/\hbar}.$$

In path-integral language, the trace identifies the endpoints and produces an imaginary-time circle of circumference $\beta\hbar$. For full grand-canonical Fock traces, number-changing bosonic and fermionic fields lead respectively to periodic and antiperiodic thermal Green functions. A fixed-number projection changes that argument: moving a creation or annihilation operator through the sector projector changes the sector. One cannot assume the same simple boundary condition for that projected trace.

This relation is structural, not literal laboratory evolution. See imaginary time for the generator and ensemble conventions, and Euclidean path integrals for the path construction.

Thermal State Versus Thermalization

The Gibbs operator specifies an equilibrium state. It does not specify how that state is prepared or reached.

A closed isolated system evolves unitarily:

$$\rho(t) = e^{-iHt/\hbar} \rho(0) e^{iHt/\hbar}.$$

Its full spectrum as a density operator is conserved, so generic unitary evolution cannot turn an arbitrary pure global state into the mixed state $\rho_\beta$. Thermalization of an isolated many-body system concerns local observables, reduced states, dephasing, and ensembles, not literal convergence of the full pure state in trace norm.

An open system can relax toward a Gibbs state under appropriate bath, weak-coupling, Markovian, and detailed-balance assumptions. A thermal fixed point is a dynamical conclusion, not part of the definition of $\rho_\beta$.

Reduced Equilibrium at Strong Coupling

Let a system $S$ interact with a bath $B$:

$$H_{\mathrm{tot}} = H_S+H_B+H_I.$$

The joint equilibrium state is

$$\rho_{SB}^{\mathrm{eq}} = \frac{ e^{-\beta H_{\mathrm{tot}}} }{ Z_{\mathrm{tot}} }.$$

The reduced system state is

$$\rho_S^{\mathrm{eq}} = \operatorname{Tr}_B \rho_{SB}^{\mathrm{eq}}.$$

When $H_I$ is not negligible, this need not equal

$$\frac{e^{-\beta H_S}}{ \operatorname{Tr}_S e^{-\beta H_S} }.$$

For finite trace-class systems define the Hamiltonian of mean force by

$$e^{-\beta H_S^*}=\frac{\operatorname{Tr}_B e^{-\beta H_{\mathrm{tot}}}}{Z_B}, \qquad Z_B=\operatorname{Tr}_B e^{-\beta H_B}.$$

Then $\rho_S^{\mathrm{eq}}=e^{-\beta H_S^*}/Z_S^*$ with $Z_S^*=Z_{\mathrm{tot}}/Z_B$. Usually $H_S^*$ depends on temperature, so derivatives of $Z_S^*$ must include that dependence. The bare Gibbs state of $H_S$ is an approximation, recovered exactly when $H_I=0$.

For the physical interpretation and thermodynamic derivatives see Strong-Coupling Open-System Effects.

Generalized Exponential States

If several commuting conserved quantities $Q_j$ are constrained, maximum entropy gives

$$\rho \propto \exp \left( -\sum_j\lambda_j Q_j \right).$$

The grand-canonical state is the important example

$$\rho_{\beta,\mu} \propto e^{-\beta(H-\mu N)}.$$

Integrable systems may require a generalized Gibbs ensemble involving many conserved charges. These states are not ordinary canonical Gibbs states unless the additional multipliers vanish or the charges are redundant.

For noncommuting constraints, the variational problem is still meaningful but its operational interpretation and parameterization require care. One should not import classical joint-distribution intuition without checking the operator structure. The general construction and Kubo–Mori response matrix are developed in Maximum Entropy Principle.

Example: Two-Level System

Choose the ground energy to be zero and let

$$H = \Delta |1\rangle\langle1|, \qquad \Delta>0.$$

Then

$$Z = 1+e^{-\beta\Delta}$$

and

$$\rho_\beta = \frac{ |0\rangle\langle0| + e^{-\beta\Delta} |1\rangle\langle1| }{ 1+e^{-\beta\Delta} }.$$

The excited-state population is

$$p_1 = \frac{1}{ e^{\beta\Delta}+1 }.$$

The mean energy is

$$U = \frac{\Delta}{ e^{\beta\Delta}+1 }.$$

As $T\to0^+$, the state approaches $|0\rangle\langle0|$. As $T\to\infty$, both levels have probability $1/2$.

The plus sign in this two-level expression does not make it a Fermi–Dirac gas. It follows from a partition function with exactly two available levels for one finite system.

Example: Harmonic Oscillator

For

$$H = \hbar\omega \left( a^\dagger a+\frac{1}{2} \right),$$

the partition function is

$$\begin{aligned} Z &= \sum_{n=0}^{\infty} e^{-\beta\hbar\omega(n+1/2)} \\ &= \frac{ e^{-\beta\hbar\omega/2} }{ 1-e^{-\beta\hbar\omega} }. \end{aligned}$$

The thermal state is diagonal in number states:

$$\rho_\beta = \left( 1-e^{-\beta\hbar\omega} \right) \sum_{n=0}^{\infty} e^{-n\beta\hbar\omega} |n\rangle\langle n|.$$

The zero-point factor cancels from the normalized state. The mean occupation is

$$\langle a^\dagger a\rangle_\beta = \frac{1}{ e^{\beta\hbar\omega}-1 }.$$

For $\beta>0$, the geometric sum converges. For $\beta<0$, it diverges, illustrating why an unbounded-above spectrum cannot support a negative-temperature Gibbs state.

Example: Free Particle and Volume Regularization

For a free particle on all of $\mathbb R^d$,

$$H = \frac{\mathbf p^2}{2m},$$

the position volume is infinite and the canonical trace diverges. Place the particle in a finite box of volume $V$ first. In the large-box continuum approximation,

$$Z_1 \sim \frac{V}{ \lambda_T^d },$$

where

$$\lambda_T = \sqrt{ \frac{ 2\pi\hbar^2 }{ mk_{\mathrm B}T } }$$

is the thermal de Broglie wavelength.

The divergence proportional to $V$ is physical state counting. It is handled by finite-volume normalization and the thermodynamic limit, not by pretending that $e^{-\beta H}$ is trace class on infinite-volume one-particle space.

Exercises

Degenerate ground space

A finite Hamiltonian has a threefold-degenerate ground energy $E_0$ and a gap to all excited states. Find the positive-temperature Gibbs limit as $T\to0^+$.

Solution

Let $P_0$ project onto the three-dimensional ground space. Factoring $e^{-\beta E_0}$ from numerator and denominator shows that excited-state terms vanish exponentially relative to the ground terms. Hence

$$\lim_{\beta\to\infty} \rho_\beta = \frac{P_0}{3}.$$

The limit is an equal mixture on the ground space. A pure symmetry-broken or boundary-selected ground state requires an additional selection procedure and a specified order of limits.

Two-level purity

For the two-level Hamiltonian

$$H = \Delta|1\rangle\langle1|,$$

compute $\operatorname{Tr}\rho_\beta^2$ and its limits as $T\to0^+$ and $T\to\infty$.

Solution

The probabilities are

$$p_0 = \frac{1}{1+e^{-\beta\Delta}}, \qquad p_1 = \frac{e^{-\beta\Delta}}{ 1+e^{-\beta\Delta} }.$$

Therefore

$$\operatorname{Tr}\rho_\beta^2 = \frac{ 1+e^{-2\beta\Delta} }{ \left( 1+e^{-\beta\Delta} \right)^2 }.$$

As $\beta\to\infty$, the purity tends to $1$. As $\beta\to0$, it tends to $1/2$, the purity of the maximally mixed qubit state.

Purity from two temperatures

Show that $\operatorname{Tr}\rho_\beta^2=Z(2\beta)/Z(\beta)^2$ when the traces exist. Why does an overall energy shift cancel?

Solution

Squaring $\rho_\beta=e^{-\beta H}/Z(\beta)$ gives $e^{-2\beta H}/Z(\beta)^2$; taking its trace proves the identity. Shifting $H$ by $C I$ multiplies the numerator by $e^{-2\beta C}$ and the denominator by the same factor. The purity is unchanged because the normalized state is unchanged.

A bath changes the effective field

Let two spins have $H_S=-h\sigma_z$, $H_B=-b\tau_z$, and $H_I=-J\sigma_z\tau_z$. Trace out the bath and calculate the ratio of the system’s $\sigma_z=+1$ and $-1$ probabilities. Is noncommutativity needed to change the bare Gibbs state?

Solution

All terms commute. For a system spin value $s=\pm1$, summing over $t=\pm1$ gives the unnormalized probability

$$w_s=\sum_{t=\pm1}e^{\beta(hs+bt+Jst)} =2e^{\beta hs}\cosh[\beta(b+Js)].$$

Consequently

$$\frac{p_+}{p_-}=e^{2\beta h} \frac{\cosh[\beta(b+J)]}{\cosh[\beta(b-J)]}.$$

The last ratio differs from one for generic nonzero $b,J$. A changed reduced equilibrium state therefore does not require noncommuting interactions. Special symmetries can still make the bare populations exact, as happens here when $b=0$.

References

• H.-P. Breuer and F. Petruccione, The Theory of Open Quantum Systems, Oxford University Press (2002).
• K. Huang, Statistical Mechanics, 2nd ed., Wiley (1987).
• M. Kardar, Statistical Physics of Particles, Cambridge University Press (2007).
• R. K. Pathria and P. D. Beale, Statistical Mechanics, 3rd ed., Elsevier (2011).
• J. J. Sakurai and J. Napolitano, Modern Quantum Mechanics, 3rd ed., Cambridge University Press (2020), chapter on identical particles and statistical mechanics.
• P. Talkner and P. Hänggi, “Colloquium: Statistical mechanics and thermodynamics at strong coupling: Quantum and classical”, Reviews of Modern Physics 92, 041002 (2020).
• J. von Neumann, Mathematical Foundations of Quantum Mechanics, Princeton University Press (1955).