An important case of a statistical mixture is that of a system in thermodynamic equilibrium with a heat reservoir at temperature T.  The various possible dynamical states are the eigenstates of the Hamiltonian H.  The statistical weight of a given eigenstate depends upon the corresponding eigenvalue of H.  It is proportional to the Boltzmann factor exp(-E/(kT)), in which E is the eigenvalue of H and k is the Boltzmann constant.  (k = 1.38*10-23J/K)

A system in thermodynamic equilibrium can be represented by the density operator
ρ = ∑npnρn = ∑npn|n><n|,
where {|n>} are the eigenstates of H,
H|n> = En|n>.
Then pn = Nexp(-En/(kT)), where N is a normalization constant to make the total probability equal to one.
ρ = ∑nNexp(-En/(kT))|n><n| = ∑nNexp(-H/(kT))|n><n| = Nexp(-H/(kT)),
since ∑n|n><n| = I.

We need Tr(ρ) = 1,   Tr(Nexp(-H/(kT)) = 1,   Tr(exp(-H/(kT)) = N-1,
N-1 is called the partition function Z, and we write  ρ = Z-1exp(-H/(kT)).

Let us calculate the partition function for the harmonic oscillator.
Z = ∑n<n|exp(-En/(kT)|)n> =  ∑nexp(-(n+½)ħω/(kT)) = exp(-ħω/(2kT)) ∑nexp(-(nħω/kT).
(1 - x)-1 = 1 + x + x2 + x3 + ...  =  ∑nxn.
Therefore ∑nexp(-(nħω/(kT) ) = ∑nexp(-(ħω/(kT))n = (1 - exp(-(ħω/(kT))-1,
and
Z = exp(-ħω/(2kT))/[1 - exp(-ħω/(kT))].

Now let us calculate the average energy.
<H> = Tr(ρH) = Z-1Tr(exp(-H/(kT)) H) = Z-1n(n+½)ħω exp(-(n+½)ħω/(kT)) = kT2(1/Z)dZ/dT,
since
dZ/dT = d/dT nexp(-(n+½)ħω/(kT)) = (1/(kT2))∑n(n+½)ħω exp(-(n+½)ħω/(kT)).

Using Z = exp(-ħω/(2kT))/[1 - exp(-ħω/(kT))] we find a simpler expression for dZ/dT.
dZ/dT = (ħω/(2kT2)) exp(-ħω/(2kT))/[1 - exp(-ħω/(kT))] +
(ħω/(kT2))exp(-ħω/(kT))exp(-ħω/(2kT))/[1 - exp(-ħω/(kT))]2
= Z(ħω/(2kT2)) + Z(ħω/(kT2)))exp(-ħω/(kT))/[1 - exp(-ħω/(kT))].

We now have <H> = ½ħω +  ħω/[exp(ħω/(kT)) - 1].
This is Planck’s formula (to within a constant ½ħω) for the average energy of a quantized oscillator.

The energy of a classical one dimensional oscillator is E(x,p) = ½p2/m + ½mω2x2.  The mean energy of such an oscillator in thermodynamic equilibrium at temperature T is
<E> = ∫∞E(x,p)exp(-E(x,p)/(kT))dxdp/∫∞exp(-E(x,p)/(kT))dxdp = kT.

 Temperature QM oscillator Classical oscillator T --> 0 ½ħω 0 kT << ħω ½ħω + ħω exp(-ħω/(kT)) kT kT >> ħω ½ħω + ħω/(ħω/(kT)) ≈ kT kT

Note: For the three-dimensional harmonic oscillator in thermodynamic equilibrium the mean energy <H> is three times that of a one-dimensional oscillator with the same frequency.