11.3.1: The harmonic Oscillator - Expansion about the Classical Path

It will be shown how to compute the density matrix for the harmonic oscillator:

\[ H = {P^2 \over 2m} + {1 \over 2}m\omega^2 X^2 \nonumber \]

using the functional integral representation. The density matrix is given by

\[ \rho(x,x';\beta) = \int_{x(0)=x}^{x(\beta\hbar)=x'}{\cal D} x (\tau ) exp \left [ - {1 \over \hbar } \int _0^{\beta \hbar } d \tau \left ({1 \over 2}m\dot{x}^2 + {1 \over 2}m\omega^2x^2\right)\right] \nonumber \]

As we saw in the last lecture, paths in the vicinity of the classical path on the inverted potential give rise to the dominant contribution to the functional integral. Thus, it proves useful to expand the path \( x (\tau )\) about the classical path. We introduce a change of path variables from \(x (\tau ) \) to \( y (\tau ) \), where

\[ x(\tau) = x_{\rm cl}(\tau) + y(\tau) \nonumber \]

where \(x_{\rm cl}(\tau)\) satisfies

\[ m\ddot{x}_{\rm cl} = m\omega^2 x_{\rm cl} \nonumber \]

subject to the conditions

\[ x_{\rm cl}(0) = x,\;\;\;\;\;\;\;\;\;\;x_{\rm cl}(\beta\hbar)=x' \nonumber \]

so that \(y (0) = y (\beta \hbar ) = 0 \).

Substituting this change of variables into the action integral yields

\[ \begin{align*} S &= \int_0^{\beta\hbar}d\tau\left[{1 \over 2}m\dot{x}^2 + {1 \over 2}m\omega^2 x^2\right] \\[4pt] &= \int_0^{\beta\hbar} d\tau\left[{1 \over 2}m(\dot{x}_{\rm cl} + \dot{y})^2 + {1 \over 2}m\omega^2(x_{\rm cl}+y)^2\right] \\[4pt] &= \int_0^{\beta\hbar} d\tau\left[{1 \over 2}m\dot{x}_{\rm cl}^2 + {1 \over 2} m \omega^2 x^2_{cl} \right ] + \int _0^{\beta \hbar}d\tau\left[{1 \over 2}m\dot{y}^2 + {1 \over 2}m\omega^2y^2\right]+ \int_0^{\beta\hbar}d\tau \left[m\dot{x}_{cl}\dot{y} +m\omega^2 x_{\rm cl}y\right] \end{align*}\]

An integration by parts makes the cross terms vanish:

\[\int_0^{\beta\hbar}d\tau\left[m\dot{x}_{\rm cl}\dot{y} +m\omega^2x_{cl} = m \dot {x}_{cl} y \right ] = m \dot {x}_{cl}y \vert _0^{\beta \hbar} + \int _0^{\beta \hbar } dtau\left[-m\ddot{x}_{\rm cl} +m\omega^2 x{\rm cl}\right]y = 0 \nonumber \]

where the surface term vanishes because \(y(0) = y ( \beta \hbar ) = 0 \) and the second term vanishes because \(xcl\) satisfies the classical equation of motion.

The first term in the expression for \(S\) is the classical action, which we have seen is given by

\[ \int_0^{\beta\hbar}d\tau\left[{1 \over 2}m\dot{x}_{\rm cl}^2 + {1 \over 2} m\omega ^2 x^2_{cl} \right] = {m\omega \over 2 \sinh (2\beta \hbar \omega} \left[(x^2 + x^{'2}){\rm cosh}(\beta\hbar\omega) - 2xx'\right] \nonumber \]

Therefore, the density matrix for the harmonic oscillator becomes

\[ \rho(x,x';\beta) = I[y]\exp\left[-{m\omega \over 2{\rm sinh} (\beta \hbar \omega )} \left ( (x^2 + x^2 + x^{'2}){\rm cosh}(\beta\hbar\omega) - 2xx'\right)\right] \nonumber \]

where \(I (y)\) is the path integral

\[ I[y] = \int_{y(0)=0}^{y(\beta\hbar)=0}{\cal D}y(\tau) exp\left [ - { 1 \over \hbar } \int _0^{\beta \hbar} \left ({m \over 2}\dot{y}^2 + {m \omega^2 \over 2}y^2\right)\right] \nonumber \]

Note that \(I [y] \) does not depend on the points \(x\) and \(x'\) and therefore can only contribute an overall (temperature dependent) constant to the density matrix. This will affect the thermodynamics but not any averages of physical observables. Nevertheless, it is important to see how such a path integral is done.

To compute \(I [y] \), we note that it is a functional integral over functions \(y (\tau) \) that vanish at \(\tau = 0 \) and \(\tau = \beta \hbar \). Thus, they are a special class of periodic functions and can be expanded in a Fourier sine series:

\[ y(\tau) = \sum_{n=1}^{\infty}c_n \sin(\omega_n\tau) \nonumber \]

where

\[ \omega_n = {n\pi\over\ \beta\hbar} \nonumber \]

Thus, we wish to change from an integral over the functions \(y (\tau ) \) to an integral over the Fourier expansion coefficients \( {c_n} \). The two integrations should be equivalent, as the coefficients uniquely determine the functions \( y (\tau ) \). Note that

\[ \dot{y}(\tau) = \sum_{n=1}^{\infty} \omega_n c_n \cos(\omega_n \tau) \nonumber \]

Thus, terms in the action are:

\[\int_0^{\beta\hbar} d\tau {1 \over m}\dot{y}^2 ={m \over 2} \sum_{n=1}^{\infty} \sum_{n'=1}^{\infty} c_nc_{n'}\omega _n\omega _{n'} \int._0^{\beta\hbar} d\tau \cos(\omega_n \tau)\cos(\omega_{n'}\tau) \nonumber \]

Since the cosines are orthogonal between \(\tau = 0\) and \(\tau = \beta \hbar \), the integral becomes

\[ \int_0^{\beta\hbar} d\tau {1 \over m}\dot{y}^2 ={m \over 2}\sum_{n=1}^{\infty} c_n^2 \omega _n^2 \int _0^{\beta \hbar} d\tau \cos ^2 (\omega _n \tau ) = {m \over 2} \sum _{n=1}^{\infty} c_n^2 \omega _n^2 \int _0^{\beta \hbar } d \tau \left [ {1 \over 2} + {1 \over 2} \cos (2 \omega _n \tau ) \right] = {m\beta\hbar \over 4}\sum_{n=1}^{\infty}c_n^2 \omega_n^2 \nonumber \]

similarly,

\[\int_0^{\beta\hbar} {1 \over 2}m\omega^2 y^2 = {m\beta\hbar \over 4}\omega^2\sum_{n=1}^{\infty}c_n^2 \nonumber \]

The measure becomes

\[{\cal D}y(t) \rightarrow \prod_{n=1}^{\infty} {dc_n \over \sqrt{4\pi/m\beta\omega_n^2}} \nonumber \]

which, is not an equivalent measure (since it is not derived from a determination of the Jacobian), but is chosen to give the correct free-particle ( \(\omega = 0\)) limit, which can ultimately be corrected by attaching an overall factor of \(\sqrt{m/2\pi\beta\hbar^2}\).

With this change of variables, \(I [y] \) becomes

\[ I[y] = \prod_{n=1}^{\infty} \int_{-\infty}^{\infty}{dc_n \over \sqrt {4\pi / m \beta \omega ^2_n}} exp \left [ - {m \beta \over 4} (\omega ^2 + \omega ^2_n )c^2_n \right ] = \prod _{n=1}^{\infty} \left[{\omega_n^2 \over \omega^2 + \omega_n^2}\right]^{1/2} \nonumber \]

The infinite product can be written as

\[\prod_{n=1}^{\infty} \left[{\pi^2 n^2/\beta^2\hbar^2 \over \omega^2 + \pi^2 n^2/ \beta^2 \hbar^2} \right ] = \left [\prod _{n=1}^{\infty} \left ( {\beta^2 \hbar^2 \omega^2 \over \pi^2 n^2}\right)\right]^{-1} \nonumber \]

the product in the square brackets is just the infinite product formula for \({\rm sinh}(\beta\hbar\omega)/(\beta\hbar\omega)\), so that \(I [y]\) is just

\[I[y] = \sqrt {\beta\hbar\omega \over {\rm sinh}(\beta\hbar\omega)} \nonumber \]

Finally, attaching the free-particle factor \(\sqrt {m/2\pi\beta\hbar^2}\), the harmonic oscillator density matrix becomes:

\[ \rho(x,x';\beta) = \sqrt {m\omega \over 2\pi\hbar{\rm sinh}(\beta \hbar \omega )} exp \left [ - {m \omega \over 2 \sinh (\beta \hbar \omega ) } \left ( (x^2 + x'^2 ) \cosh (\beta \hbar \omega ) - 2xx' \right ) \right] \nonumber \]

Notice that in the free-particle limit \( (\omega \rightarrow 0 ) \), \(\sinh (\beta \hbar \omega ) \approx \beta \hbar \omega \) and \(\cosh (\beta \hbar \omega ) \approx 1 \), so that

\[ \rho(x,x';\beta) \rightarrow\sqrt {m \over 2\pi\beta\hbar^2} \exp\left[-{m \over 2\beta\hbar^2}(x-x')^2\right] \nonumber \]

which is the expected free-particle density matrix.