I was on vacation last week, and during some quiet time, I randomly happened to be thinking about explanations for damping in physical systems. I remember learning in ELE 456 — Quantum Optics, from last spring, that the phenomenological linear damping of a classical oscillator could be derived by coupling a quantum oscillator to a thermal bath of quantum oscillators; each linear oscillator is microscopically undamped, but by treating the bath through statistical thermodynamics, the coupling of the oscillator in question to a bath essentially produces a linear damping coefficient dependent on the spectrum of the bath (and the coupling too). Microscopically, the quantization of energy levels in a linear oscillator makes it easy to interpret how discrete excitations can move from one oscillator to another coupled oscillator, but I was wondering if quantum mechanics is really necessary to explain damping. Follow the jump to see an extremely rough sketch of ideas that may (or may not) justify the use of classical mechanics by itself. (Added after finishing: this turns out to be a rambling and possibly ultimately pointless post with a much clearer and more self-consistent explanation linked at the end, so for the time being, humor me.)
A counterexample (to the notion of quantization being required) that I thought of was that of a nonrelativistic ideal gas under no external potential and under the influence of only short-ranged two-body interactions, such that collisions between particles in the gas could be approximated by hard-sphere collisions (with the details specified by the particular interaction potential). Looking through the notes for 8.333 — Statistical Mechanics I at MIT, the BBGKY hierarchy that relates the $s$-particle density evolution with the collisions involving the $s + 1$-particle density, the microscopic interactions are fully reversible. Yet, by ignoring all $s$-body interactions for $s > 2$ (so considering only the 1- and 2-particle densities), by considering the 2-particle density to evolve in a manner similar to isolated 1-particle densities (due to the time between collisions being much larger than the time of a collision for a dilute gas, so the relaxation to independent 1-particle densities happens much faster than anything else), and by considering the 2-particle density (being a product of 1-particle densities) to be correlated after but not before a collision, information about the full collisional behavior is being thrown away, and the resulting Boltzmann equation displays irreversible macroscopic relaxation arising from reversible microscopic dynamics. Going a little further with this analysis by considering a Gaussian 1-particle density function and the corresponding first-order expansion in the Boltzmann equation yields the Euler equations for fluid flow, replete with viscous and thermal conduction terms that describe the damping/relaxation to equilibrium dependent on the timescale of the collision itself. At least superficially (and probably a little deeper than that too), the arguments in this analysis, based on timescales and correlations of energy transfer, are similar to those used to describe the macroscopic damping of a quantum oscillator coupled to a thermal bath of oscillators, yet this analysis requires neither quantization (as it is purely classical) nor even an explicit form of the 2-particle interaction interaction (only that it is short-range, so that the collisions can be cast in the form of hard spheres).
(As an aside, I'd like to point out how amazing the books by Mehran Kardar are. As I have progressed in my physics education, it's cool to see the connections between classical and quantum statistical mechanics becoming clearer to me. A particular example would be something like $f_{1} (\vec{q}, \vec{p}, t) = \langle \sum_{n = 1}^{N} \delta^{3} (\vec{q} - \vec{q}_{n}) \delta^{3} (\vec{p} - \vec{p}_{n}) \rangle$. The expectation value of a Dirac delta function suggests that this is essentially projecting onto a subspace of a smaller number of particles, tracing over the density of all of the other particles; indeed, that is exactly what is happening, and in quantum mechanics, the integral over the $N$-particle density multiplied by the Dirac delta function, which traces over classical states, would simply be replaced by a quantum trace over the other particles of the full density operator multiplied by an appropriate projection operator. It really helps that the notation in those books is consistent and clear in itself, and it helps to clarify some of these connections too.)
With that out of the way, what else can be done here? What I'd really like to see is whether I can recover linear damping of a harmonic oscillator coupled to a bath of linear oscillators by treating everything classically. Hopefully it works!
Let's start with the Hamiltonian for a linear harmonic oscillator at a given spring constant coupled to other oscillators of identical mass but with different spring constants (but for which the other oscillators are independent from each other). This looks like \[ H = \frac{p_{0}^2}{2m} + \frac{m\omega_{0}^2 x_{0}^2}{2} + \sum_{n} \left(\frac{p_{n}^2}{2m} + \frac{m\omega_{n}^2 x_{n}^2}{2} + r_{n} x_{0} x_{n} \right) \] where $r_{n}$ is the coupling of the given oscillator to each of the other oscillators in the bath. (Strictly speaking, the coupling would be of the form $\frac{r_{n}}{2} (x_{n} - x_{0})^2$, but expanding the square simply renormalizes the spring constants/frequencies of both the given oscillator as well as of all of the oscillators in the bath, leaving only $-r_{n} x_{0} x_{n}$, from which the sign can be absorbed into the definition of $r_{n}$.) I've mentioned in previous posts that it is a perfectly valid canonical transformation to go from $(x_{0}, p_{0})$ and $(x_{n}, p_{n})$ to $(\alpha_{0}, \alpha_{0}^{\star})$ and $(\beta_{n}, \beta_{n}^{\star})$ where $\alpha_{0} = \sqrt{\frac{m\omega_{0}}{2\hbar}} \left(x_{0} + \frac{\mathrm{i}}{m\omega_{0}} p_{0}\right)$ (with similar expressions for $\beta_{n}$), though strictly speaking, $\hbar$ classically is a somewhat arbitrary constant with units of action to link the energy of the oscillator to its frequency. These satisfy the Poisson bracket relations $\{\alpha_{0}, \alpha_{0}^{\star}\} = -\frac{\mathrm{i}}{\hbar}$ and $\{\beta_{m}, \beta_{n}\} = -\frac{\mathrm{i}}{\hbar} \delta_{mn}$, with all other Poisson brackets vanishing. The resulting Hamiltonian is then \[ H = \hbar\omega_{0} |\alpha_{0}|^2 + \sum_{n} \left(\hbar\omega_{n} |\beta_{n}|^2 + g_{n} (\alpha_{0} + \alpha_{0}^{\star})(\beta_{n} + \beta_{n}^{\star})\right) \] where $g_{n} = \frac{\hbar^{2} r_{n}}{2m \sqrt{\omega_{0} \omega_{n}}}$.
At this point, it is helpful to note that the natural oscillation frequency of $\alpha_{0}$ is $\omega_{0}$, and of $\beta_{n}$ is $\omega_{n}$, so new variables may be defined through $\alpha_{0} = a_{0} e^{-\mathrm{i}\omega_{0} t}$ and $\beta_{n} = b_{n} e^{-\mathrm{i} \omega_{n} t}$, which still satisfy $\{a_{0}, a_{0}^{\star}\} = -\frac{\mathrm{i}}{\hbar}$ and $\{b_{m}, b_{n}^{\star}\} = -\frac{\mathrm{i}}{\hbar} \delta_{mn}$, with all other Poisson brackets vanishing. Plugging this into the $x_{0} x_{n}$ term gives terms that go as $a_{0} b_{n} e^{-\mathrm{i} (\omega_{0} + \omega_{n})t} + a_{0} b_{n}^{\star} e^{-\mathrm{i} (\omega_{0} - \omega_{n})t} + \mathrm{c.c.}$; here, the terms involving $e^{-\mathrm{i} (\omega_{0} - \omega{n})t}$ and its conjugate rotate slowly compared to the terms involving $e^{-\mathrm{i} (\omega_{0} + \omega_{n})t}$ and its conjugate, so the latter terms, averaging to zero much quicker through fast oscillation, can be neglected. The resulting equations of motion are \[ \dot{a}_{0} = -\mathrm{i} \sum_{n} g_{n} b_{n} e^{-\mathrm{i} (\omega_{n} - \omega_{0})t} \] for the system and \[ \dot{b}_{n} = -\mathrm{i} g_{n} a_{0} e^{\mathrm{i} (\omega_{n} - \omega_{0})t} \] for the bath. After formally integrating the equation for the bath variables, these can be combined to yield the integrodifferential equation \[ \dot{a}_{0} = \tilde{f}(t) - \int_{0}^{t} K(\tau) a_{0} (t - \tau)~\mathrm{d}\tau \] for the system, where \[ \tilde{f}(t) = -\mathrm{i} \sum_{n} g_{n} b_{n} (0) e^{-\mathrm{i}(\omega_{n} - \omega_{0})t} \] expresses fluctuations in the bath oscillators in terms of their initial conditions and the couplings to the system, and \[ K(\tau) = \sum_{n} g_{n}^2 e^{-\mathrm{i}(\omega_{n} - \omega_{0})\tau} \] is the convolution kernel that only depends on the coupling of the system to the bath. Here, if the bath is thermodynamically large, then the summation over bath oscillators $n$ can be replaced by a continuous integral modulated by the density of states $\rho(\omega)$ describing the number of oscillators per unit frequency characterized by a given frequency. This allows writing \[ K(\tau) = \int \rho(\omega) g^2 (\omega) e^{-\mathrm{i}(\omega - \omega_{0})\tau}~\mathrm{d}\omega \] from which it is helpful to denote the spectral function of the bath as $J(\omega) = \rho(\omega) g^2 (\omega)$ for brevity.
Now, the approximations begin. The first is the Markov approximation, saying that $K(\tau)$ is characterized by a very short correlation time, meaning that it only contributes to the convolution with $a_{0}$ in the integrodifferential equation if it is sharply peaked around $t - t' = \tau = 0$. This means that $\int_{0}^{t} K(\tau) a_{0} (t-\tau)~\mathrm{d}\tau \approx a_{0} (t) \int_{0}^{t} K(\tau)~\mathrm{d}\tau$. The second goes further, saying that because $K(\tau)$ decays so much for any $\tau > 0$, then $\int_{0}^{t} K(\tau)~\mathrm{d}\tau \approx \int_{0}^{\infty} K(\tau)~\mathrm{d}\tau \equiv \Gamma$ is a constant. Plugging in $K(\tau) = \int J(\omega) e^{-\mathrm{i}(\omega - \omega_{0})\tau}~\mathrm{d}\omega$ and using the integral of a Dirac delta function over half of the real line gives a decay rate $\kappa = 2\pi J(\omega_{0})$, along with a frequency shift (Lamb shift) denoted for simplicity as $\Delta\omega$. This yields \[ \dot{a}_{0} + \left(\mathrm{i}\Delta \omega + \frac{\kappa}{2}\right)a_{0} = \tilde{f}(t) \] as the equation of motion for the system oscillator.
Assuming that the initial state of the bath is thermal with $\langle \beta_{n} (0) \rangle = 0$, classical statistical mechanics would say that the correlations of the bath variables obey $\langle \beta_{m}^{\star} (0) \beta_{n} (0) \rangle = \frac{k_{\mathrm{B}} T}{\hbar \omega_{n}} \delta_{mn}$. Plugging this into the definition of $\tilde{f}(t)$ and applying the Markov approximation consistently yields $\langle \tilde{f}^{\star} (t) \tilde{f}(t') \rangle \approx \frac{k_{\mathrm{B}} T\kappa}{\hbar\omega_{0}} \delta(t - t')$. The resulting equation of motion for the system oscillator energy $|a_{0}|^2$ (having been normalized to $\hbar\omega_{0}$) is $\left(\frac{\mathrm{d}}{\mathrm{d}t} + \kappa\right)|a_{0}|^2 = \frac{k_{\mathrm{B}} T\kappa}{\hbar\omega_{0}}$, which in the steady-state becomes $\hbar\omega_{0} |a_{0}|^2 \to k_{\mathrm{B}} T$ exactly as desired for a classical harmonic oscillator connected to a thermal bath.
Thus, the damping of an oscillator can be derived classically too. That said, the damping coefficient depends on the properties of the bath, for which $\hbar$ will sneak in one way or another; moreover, the use of $\hbar$ in the first place is essentially a dimensional crutch that doesn't make as much sense except with a knowledge of quantum mechanics. Mathematically, though, the point is that a notion of damping, as in the case of irreversible dynamics and viscous relaxation in a classical gas, does not necessarily require quantum mechanics from the start. (Of course, I could have saved all of this trouble by simply doing the true quantum operator derivation and then taking the classical limit $k_{\mathrm{B}} T \gg \hbar\omega$ everywhere, but my goal was to see if I could do the same using commuting classical variables and Poisson brackets everywhere; I can, and the point is that the lack of commutation of observables becomes irrelevant in the classical limit, because for large enough $T$, the Bose-Einstein occupation number is so large that the added contribution from lack of commutation won't make a significant difference anywhere. The reason this post seems a bit anticlimactic is because after deriving the equations of motion, whether for classical observables or quantum operators, the other derivations are all the same; moreover, $\hbar$ itself doesn't appear in the general equation for a damped oscillator, but either for dimensional reasons (classically) or for more fundamental reasons (quantum mechanically), it has to appear in expressions related to the bath or the coupling constants. With all that said, I found a much nicer derivation by Prof. Dmitry Garanin at CUNY that only considers the mechanical properties of the bath, without appealing to $\hbar$ or thermal considerations, to more cleanly end up with the same result.)
A counterexample (to the notion of quantization being required) that I thought of was that of a nonrelativistic ideal gas under no external potential and under the influence of only short-ranged two-body interactions, such that collisions between particles in the gas could be approximated by hard-sphere collisions (with the details specified by the particular interaction potential). Looking through the notes for 8.333 — Statistical Mechanics I at MIT, the BBGKY hierarchy that relates the $s$-particle density evolution with the collisions involving the $s + 1$-particle density, the microscopic interactions are fully reversible. Yet, by ignoring all $s$-body interactions for $s > 2$ (so considering only the 1- and 2-particle densities), by considering the 2-particle density to evolve in a manner similar to isolated 1-particle densities (due to the time between collisions being much larger than the time of a collision for a dilute gas, so the relaxation to independent 1-particle densities happens much faster than anything else), and by considering the 2-particle density (being a product of 1-particle densities) to be correlated after but not before a collision, information about the full collisional behavior is being thrown away, and the resulting Boltzmann equation displays irreversible macroscopic relaxation arising from reversible microscopic dynamics. Going a little further with this analysis by considering a Gaussian 1-particle density function and the corresponding first-order expansion in the Boltzmann equation yields the Euler equations for fluid flow, replete with viscous and thermal conduction terms that describe the damping/relaxation to equilibrium dependent on the timescale of the collision itself. At least superficially (and probably a little deeper than that too), the arguments in this analysis, based on timescales and correlations of energy transfer, are similar to those used to describe the macroscopic damping of a quantum oscillator coupled to a thermal bath of oscillators, yet this analysis requires neither quantization (as it is purely classical) nor even an explicit form of the 2-particle interaction interaction (only that it is short-range, so that the collisions can be cast in the form of hard spheres).
(As an aside, I'd like to point out how amazing the books by Mehran Kardar are. As I have progressed in my physics education, it's cool to see the connections between classical and quantum statistical mechanics becoming clearer to me. A particular example would be something like $f_{1} (\vec{q}, \vec{p}, t) = \langle \sum_{n = 1}^{N} \delta^{3} (\vec{q} - \vec{q}_{n}) \delta^{3} (\vec{p} - \vec{p}_{n}) \rangle$. The expectation value of a Dirac delta function suggests that this is essentially projecting onto a subspace of a smaller number of particles, tracing over the density of all of the other particles; indeed, that is exactly what is happening, and in quantum mechanics, the integral over the $N$-particle density multiplied by the Dirac delta function, which traces over classical states, would simply be replaced by a quantum trace over the other particles of the full density operator multiplied by an appropriate projection operator. It really helps that the notation in those books is consistent and clear in itself, and it helps to clarify some of these connections too.)
With that out of the way, what else can be done here? What I'd really like to see is whether I can recover linear damping of a harmonic oscillator coupled to a bath of linear oscillators by treating everything classically. Hopefully it works!
Let's start with the Hamiltonian for a linear harmonic oscillator at a given spring constant coupled to other oscillators of identical mass but with different spring constants (but for which the other oscillators are independent from each other). This looks like \[ H = \frac{p_{0}^2}{2m} + \frac{m\omega_{0}^2 x_{0}^2}{2} + \sum_{n} \left(\frac{p_{n}^2}{2m} + \frac{m\omega_{n}^2 x_{n}^2}{2} + r_{n} x_{0} x_{n} \right) \] where $r_{n}$ is the coupling of the given oscillator to each of the other oscillators in the bath. (Strictly speaking, the coupling would be of the form $\frac{r_{n}}{2} (x_{n} - x_{0})^2$, but expanding the square simply renormalizes the spring constants/frequencies of both the given oscillator as well as of all of the oscillators in the bath, leaving only $-r_{n} x_{0} x_{n}$, from which the sign can be absorbed into the definition of $r_{n}$.) I've mentioned in previous posts that it is a perfectly valid canonical transformation to go from $(x_{0}, p_{0})$ and $(x_{n}, p_{n})$ to $(\alpha_{0}, \alpha_{0}^{\star})$ and $(\beta_{n}, \beta_{n}^{\star})$ where $\alpha_{0} = \sqrt{\frac{m\omega_{0}}{2\hbar}} \left(x_{0} + \frac{\mathrm{i}}{m\omega_{0}} p_{0}\right)$ (with similar expressions for $\beta_{n}$), though strictly speaking, $\hbar$ classically is a somewhat arbitrary constant with units of action to link the energy of the oscillator to its frequency. These satisfy the Poisson bracket relations $\{\alpha_{0}, \alpha_{0}^{\star}\} = -\frac{\mathrm{i}}{\hbar}$ and $\{\beta_{m}, \beta_{n}\} = -\frac{\mathrm{i}}{\hbar} \delta_{mn}$, with all other Poisson brackets vanishing. The resulting Hamiltonian is then \[ H = \hbar\omega_{0} |\alpha_{0}|^2 + \sum_{n} \left(\hbar\omega_{n} |\beta_{n}|^2 + g_{n} (\alpha_{0} + \alpha_{0}^{\star})(\beta_{n} + \beta_{n}^{\star})\right) \] where $g_{n} = \frac{\hbar^{2} r_{n}}{2m \sqrt{\omega_{0} \omega_{n}}}$.
At this point, it is helpful to note that the natural oscillation frequency of $\alpha_{0}$ is $\omega_{0}$, and of $\beta_{n}$ is $\omega_{n}$, so new variables may be defined through $\alpha_{0} = a_{0} e^{-\mathrm{i}\omega_{0} t}$ and $\beta_{n} = b_{n} e^{-\mathrm{i} \omega_{n} t}$, which still satisfy $\{a_{0}, a_{0}^{\star}\} = -\frac{\mathrm{i}}{\hbar}$ and $\{b_{m}, b_{n}^{\star}\} = -\frac{\mathrm{i}}{\hbar} \delta_{mn}$, with all other Poisson brackets vanishing. Plugging this into the $x_{0} x_{n}$ term gives terms that go as $a_{0} b_{n} e^{-\mathrm{i} (\omega_{0} + \omega_{n})t} + a_{0} b_{n}^{\star} e^{-\mathrm{i} (\omega_{0} - \omega_{n})t} + \mathrm{c.c.}$; here, the terms involving $e^{-\mathrm{i} (\omega_{0} - \omega{n})t}$ and its conjugate rotate slowly compared to the terms involving $e^{-\mathrm{i} (\omega_{0} + \omega_{n})t}$ and its conjugate, so the latter terms, averaging to zero much quicker through fast oscillation, can be neglected. The resulting equations of motion are \[ \dot{a}_{0} = -\mathrm{i} \sum_{n} g_{n} b_{n} e^{-\mathrm{i} (\omega_{n} - \omega_{0})t} \] for the system and \[ \dot{b}_{n} = -\mathrm{i} g_{n} a_{0} e^{\mathrm{i} (\omega_{n} - \omega_{0})t} \] for the bath. After formally integrating the equation for the bath variables, these can be combined to yield the integrodifferential equation \[ \dot{a}_{0} = \tilde{f}(t) - \int_{0}^{t} K(\tau) a_{0} (t - \tau)~\mathrm{d}\tau \] for the system, where \[ \tilde{f}(t) = -\mathrm{i} \sum_{n} g_{n} b_{n} (0) e^{-\mathrm{i}(\omega_{n} - \omega_{0})t} \] expresses fluctuations in the bath oscillators in terms of their initial conditions and the couplings to the system, and \[ K(\tau) = \sum_{n} g_{n}^2 e^{-\mathrm{i}(\omega_{n} - \omega_{0})\tau} \] is the convolution kernel that only depends on the coupling of the system to the bath. Here, if the bath is thermodynamically large, then the summation over bath oscillators $n$ can be replaced by a continuous integral modulated by the density of states $\rho(\omega)$ describing the number of oscillators per unit frequency characterized by a given frequency. This allows writing \[ K(\tau) = \int \rho(\omega) g^2 (\omega) e^{-\mathrm{i}(\omega - \omega_{0})\tau}~\mathrm{d}\omega \] from which it is helpful to denote the spectral function of the bath as $J(\omega) = \rho(\omega) g^2 (\omega)$ for brevity.
Now, the approximations begin. The first is the Markov approximation, saying that $K(\tau)$ is characterized by a very short correlation time, meaning that it only contributes to the convolution with $a_{0}$ in the integrodifferential equation if it is sharply peaked around $t - t' = \tau = 0$. This means that $\int_{0}^{t} K(\tau) a_{0} (t-\tau)~\mathrm{d}\tau \approx a_{0} (t) \int_{0}^{t} K(\tau)~\mathrm{d}\tau$. The second goes further, saying that because $K(\tau)$ decays so much for any $\tau > 0$, then $\int_{0}^{t} K(\tau)~\mathrm{d}\tau \approx \int_{0}^{\infty} K(\tau)~\mathrm{d}\tau \equiv \Gamma$ is a constant. Plugging in $K(\tau) = \int J(\omega) e^{-\mathrm{i}(\omega - \omega_{0})\tau}~\mathrm{d}\omega$ and using the integral of a Dirac delta function over half of the real line gives a decay rate $\kappa = 2\pi J(\omega_{0})$, along with a frequency shift (Lamb shift) denoted for simplicity as $\Delta\omega$. This yields \[ \dot{a}_{0} + \left(\mathrm{i}\Delta \omega + \frac{\kappa}{2}\right)a_{0} = \tilde{f}(t) \] as the equation of motion for the system oscillator.
Assuming that the initial state of the bath is thermal with $\langle \beta_{n} (0) \rangle = 0$, classical statistical mechanics would say that the correlations of the bath variables obey $\langle \beta_{m}^{\star} (0) \beta_{n} (0) \rangle = \frac{k_{\mathrm{B}} T}{\hbar \omega_{n}} \delta_{mn}$. Plugging this into the definition of $\tilde{f}(t)$ and applying the Markov approximation consistently yields $\langle \tilde{f}^{\star} (t) \tilde{f}(t') \rangle \approx \frac{k_{\mathrm{B}} T\kappa}{\hbar\omega_{0}} \delta(t - t')$. The resulting equation of motion for the system oscillator energy $|a_{0}|^2$ (having been normalized to $\hbar\omega_{0}$) is $\left(\frac{\mathrm{d}}{\mathrm{d}t} + \kappa\right)|a_{0}|^2 = \frac{k_{\mathrm{B}} T\kappa}{\hbar\omega_{0}}$, which in the steady-state becomes $\hbar\omega_{0} |a_{0}|^2 \to k_{\mathrm{B}} T$ exactly as desired for a classical harmonic oscillator connected to a thermal bath.
Thus, the damping of an oscillator can be derived classically too. That said, the damping coefficient depends on the properties of the bath, for which $\hbar$ will sneak in one way or another; moreover, the use of $\hbar$ in the first place is essentially a dimensional crutch that doesn't make as much sense except with a knowledge of quantum mechanics. Mathematically, though, the point is that a notion of damping, as in the case of irreversible dynamics and viscous relaxation in a classical gas, does not necessarily require quantum mechanics from the start. (Of course, I could have saved all of this trouble by simply doing the true quantum operator derivation and then taking the classical limit $k_{\mathrm{B}} T \gg \hbar\omega$ everywhere, but my goal was to see if I could do the same using commuting classical variables and Poisson brackets everywhere; I can, and the point is that the lack of commutation of observables becomes irrelevant in the classical limit, because for large enough $T$, the Bose-Einstein occupation number is so large that the added contribution from lack of commutation won't make a significant difference anywhere. The reason this post seems a bit anticlimactic is because after deriving the equations of motion, whether for classical observables or quantum operators, the other derivations are all the same; moreover, $\hbar$ itself doesn't appear in the general equation for a damped oscillator, but either for dimensional reasons (classically) or for more fundamental reasons (quantum mechanically), it has to appear in expressions related to the bath or the coupling constants. With all that said, I found a much nicer derivation by Prof. Dmitry Garanin at CUNY that only considers the mechanical properties of the bath, without appealing to $\hbar$ or thermal considerations, to more cleanly end up with the same result.)