The density operator of a quantum system in thermal equilibrium at temperature \(T\) is given by \[\rho_T = \frac{e^{-H / k_{\text{B}} T}}{Z} \,,\] where \(H\) is the system Hamiltonian, \(k_{\text{B}}\) is the Boltzmann constant, and \[Z = \operatorname{Tr} e^{-H / k_{\text{B}} T}\] is the partition function. This expression for the density matrix can be obtained by extremizing the von Neumann entropy \[S = -k_{\text{B}} \operatorname{Tr} \rho \ln \rho\] over all states \(\rho\) of the same mean energy \[E = \operatorname{Tr} \rho H \,.\] Proof: We want to extremize \(S\) subject to (i) the normalization constraint, \(\operatorname{Tr} \rho = 1\), and (ii) the mean energy constraint, \(\operatorname{Tr} \rho H = E\). To this end, we introduce two Lagrange multipliers, \(\alpha k_{\text{B}}\) and \(\beta k_{\text{B}}\), and perform unconstrained extremization of \[\mathcal{F} = -k_{\text{B}} \operatorname{Tr} \rho \ln \rho - \alpha k_{\text{B}} \left( \operatorname{Tr}