Why do molecules behave classically at high temperatures? In chapter 6 of "Charge and Energy Transfer Dynamics in Molecular Systems", the authors consider diatomic molecules in contact with a heat bath. We have established before that the nuclear coordinates $R$ are quantum mechanical and can be described in the adiabatic (Born-Oppenheimer) approximation. In this approximation, the nuclear motion is governed by the Hamiltonian
$$ H_{nuc} = T_{nuc}(R)+V_{nuc-nuc}(R)+E_a(R),$$
where $E_a$ is the ($R$-dependent) eigenvalue of the electronic Hamiltonian,
$$ \big(T_{el}+V_{el-el}+V_{nuc-el}(R) \big)\phi_a(r,R)=E_a(R)\phi_a(r,R).$$
Usually, the potential term in $H_{nuc}$ has a local minimum and can be Taylor expanded to yield a harmonic oscillator problem. The nuclear eigenfunctions can be approximated by the eigenfunctions of the harmonic oscillator.
Now the authors claim that if the temperature is much larger than the oscillator energy $k_B T >>\hbar \omega$, the molecule coordinates can be treated classically, which means in this case that they are localized at a certain value $R$, following the statistical distribution
$$P(R)=\frac{1}{Z} e^{-(V_{nuc-nuc}(R)+E_a(R))/k_B T} \tag{1}$$
(This also happens in a lot of molecular dynamics simulations: The molecular coordinates are localized and follow classical trajectories.)
Question: Why are the nuclear coordinates localized like this?
Naively, I would have expected that the nuclear system is described by the reduced density operator (Gibbs ensemble) $\rho=\frac{1}{Z} e^{-H_{nuc}/k_B T}$. There, the nuclei would occupy eigenstates of $H_{nuc}$, with the probability depending on the energy of that state. Clearly the eigenstates are delocalized across the whole harmonic potential, especially for high energy states.
Edit: I would like to point out another disparity between the QM and classical descriptions: The statistical operator $\rho$ does not evolve in time since $[H_{nuc},\rho]=0$ (neglecting non-adiabatic processes). The classical coordinates however oscillate with frequency $\omega$.
 A: The statement that "the molecules behave (approximately) classically" is really just saying that the classical statistical ensemble gives (approximately) the same results as the quantum Gibbs ensemble for any physically measurable quantity. So the claim is not that the quantum system doesn't approach the quantum Gibbs ensemble, only that it isn't necessary to do the full quantum analysis of the system because the classical ensemble will give the same result.
It may not be immediately obvious why the classical ensemble gives the same results as the quantum ensemble when they are constructed in very different ways (probability distribution in phase space, vs. probability distribution over eigenstates). I don't have a particularly satisfying answer to that; I think that the origins of classical-quantum correspondence in statistical mechanics is quite a subtle issue to really understand.
A: In a simplified diatomic molecule model approximation electronic and nuclear motions can essentially be treated independently because of the small ratio of the electron mass with respect to the nuclear mass $m_e/M_{nucl}\simeq 10^{-3}-10^{-5}$. The translational motion can be separated by introducing the centre of mass, which moves as a free particle. The other motions involved are the vibrational and the rotational ones.
In the high temperature limit, it can be shown that the quantum (q) vibrational and rotational partition functions tend to the classical ones (c).
$$Q_{vib}^{q}=\sum_{n=0}^{\infty}e^{-\beta \hbar \omega(n+1/2)}=\frac{e^{-\beta \hbar \omega/2}}{1-e^{-\beta \hbar \omega}}\rightarrow Q_{vib}^{c}=\frac{kT}{\hbar\omega}$$
The quantum rotational partition function can be written using the allowed values of the angular momentum:
$$Q_{rot}^{q}=\sum_{l=0}(2l+1)e^{-\beta \hbar^2 l(l+1)/2I}$$
For sufficiently high values of $T$, the terms in this equation vary slowly and the sum can be replaced by an integral of the form
$$Q_{rot}^{q}\rightarrow \int_{0}^{\infty}dx \ (2x+1)e^{-\beta \hbar^2 x(x+1)/2I}=\frac{2I kT}{\hbar^2}=Q_{rot}^{c}$$
Determining the temperature at which molecules begin to behave classically requires some complicated calculations, in general it depends on the specific volume $v_s=V/N$ available to the molecule and on the thermal volume $v_t=\lambda^3$ (or also De Broglie thermal wavelenght cubed). We get the classical regime when $v_s>v_t$.
A: This answer won't be very technical, but I'll try to be helpful with a rather elementary explanation that is intuitive. If the gasses particles all have high kinetic energy (also just known as heat/temperature), then they will more or less collied with each other perfectly elastically. If the gas particles had lower kinetic energy, then the particles' Van der Waals forces (Dipole attractions and/or dispersion forces) could cause the particles to stick together, making the collision somewhat inelastic. Liquids consist of inelastic collisions.
When particles are colliding perfectly elastically, we can predict them quite well with statistical models, such as the Maxwell-Boltzmann distribution.
I hope this helps, let me know if you need any clarification.
