# Coupled quantum harmonic oscillators (exact $\neq$ perturbative)

Suppose we are given the Hamiltonian $$\hat H = \hat H_0 + \hat H_p(\varepsilon) = \frac 1 {2m}(\hat p_1^2 + \hat p_2^2) +\frac 1 2 m \omega^2(\hat x_1^2 + \hat x_2^2) + \varepsilon m\omega^2\hat x_1\hat x_2$$ and, after switching to CM coordinates, $$\begin{split} \hat X= \frac {\hat x_1 + \hat x_2}{2}, &\qquad \hat \rho = \hat x_1 -\hat x_2, \\ \hat P = \hat p_1 + \hat p_2, &\qquad \hat \pi= \frac {\hat p_1-\hat p_2}{2}, \end{split}$$ we are asked to solve the corresponding eigenvalue problem $$\hat H \psi = E \psi$$ first exactly and then perturbatively. I managed to find the following exact (normalized) solution: $$\begin{split} \psi_{00}(X,\rho) &= \Gamma_{00} \exp\left\{-\frac{m\omega}{2\hbar}\left(X^2\sqrt{1+\varepsilon} + \rho^2\sqrt{1-\varepsilon} \right) \right\}, \\ \Gamma_{00} &= \sqrt{\frac{m\omega}{\pi\hbar}\sqrt[4]{(1+\varepsilon)(1-\varepsilon)}}, \\ \psi_{n\ell}(X,\rho) &= \frac{1}{\sqrt{2^{n+\ell}n!\ell!}}\operatorname{He}_n\left(X\sqrt{\frac{2m\omega}{\hbar}\sqrt{1+\varepsilon}} \right) \operatorname{He}_\ell\left(\rho \sqrt{\frac{m\omega}{2\hbar}\sqrt{1-\varepsilon}} \right) \psi_{00}(X,\rho), \\ E_{n\ell} &= \hbar \omega \left\{\left(n+\frac 1 2 \right) \sqrt{1+\varepsilon} + \left(\ell + \frac 12 \right)\sqrt{1-\varepsilon}\right\}, \end{split}$$ where $$\operatorname{He}_q(\xi)$$ is the usual $$q$$-th Hermite polynomial. Treating the problem perturbatively, again in the new coordinates, the analysis leads to the first-order energy corrections $$\delta^{(1)}E_{00} = 0, \qquad \delta^{(1)}E_{01} = -\frac 1 2 \varepsilon\hbar \omega, \qquad \delta^{(1)}E_{10} = \frac12 \varepsilon \hbar \omega,$$ which are perfectly compatible with the first-order expansion of the exact energies. However, the first-order correction to the unperturbed eigenstate $$\psi_{00}^{(0)}$$ (I haven't tackled the others) does not match what I get when expanding $$\psi_{00}$$ above to first order in $$\varepsilon$$. Why is this?

Calculations. The theory of stationary perturbations gives $$\delta^{(1)}\psi_{00} = \sum_{n,\ell\neq 0} \frac{\langle \psi_{n\ell}^{(0)}| \hat H_p \psi_{00}^{(0)} \rangle}{E_{00}^{(0)}-E_{n\ell}^{(0)}} \psi_{n\ell}^{(0)};$$ the matrix element in the numerator can be calculated through appropriate ladder operators $$\hat A$$ and $$\hat \alpha$$ (respectively for the unperturbed oscillator in $$\hat X$$ and $$\hat \rho$$) to be, knowing $$\hat H_p = \varepsilon m \omega^2(\hat X^2 - \frac 1 4 \hat \rho^2)$$ and remembering $$(\hat a^\dagger + \hat a)^2 = (\hat a^\dagger)^2 + \hat a^2 + 2\hat a^\dagger \hat a + \hat 1$$ when $$[\hat a,\hat a^\dagger] = \hat 1$$, $$\begin{split} \langle \psi_{n\ell}^{(0)}| \hat H_p \psi_{00}^{(0)} \rangle &= \varepsilon m \omega^2 \Bigg\{ \frac{\hbar}{4m\omega} \left( \langle \psi_{n\ell}^{(0)}| (\hat A^\dagger)^2 \psi_{00}^{(0)} \rangle + \underbrace{\langle \psi_{n\ell}^{(0)}| (\hat A)^2 \psi_{00}^{(0)} \rangle}_{=0} + 2\underbrace{\langle \psi_{n\ell}^{(0)}| \hat A^\dagger \hat A \psi_{00}^{(0)}\rangle}_{=0} + \underbrace{\langle \psi_{n\ell}^{(0)}| \psi_{00}^{(0)} \rangle}_{=0} \right) \\ &\qquad\qquad - \frac{\hbar}{4m\omega} \left( \langle \psi_{n\ell}^{(0)}| (\hat \alpha^\dagger)^2 \psi_{00}^{(0)} \rangle + \underbrace{\langle \psi_{n\ell}^{(0)}| (\hat \alpha)^2 \psi_{00}^{(0)} \rangle}_{=0} + 2\underbrace{\langle \psi_{n\ell}^{(0)}| \hat \alpha^\dagger \hat \alpha \psi_{00}^{(0)}\rangle}_{=0} + \underbrace{\langle \psi_{n\ell}^{(0)}| \psi_{00}^{(0)} \rangle}_{=0} \right) \Bigg\} \\ &= \frac{1}{4}\varepsilon \hbar \omega \left\{\sqrt 2 \langle \psi_{n\ell}^{(0)}| \psi_{20}^{(0)} \rangle - \sqrt 2 \langle \psi_{n\ell}^{(0)}| \psi_{02}^{(0)} \rangle \right\}, \end{split}$$ so that $$\delta^{(1)}\psi_{00} = \frac{\varepsilon \hbar \omega}{2\sqrt 2} \left[ \frac{1}{-\hbar\omega(2+0)}\psi_{20}^{(0)} - \frac{1}{-\hbar\omega(0+2)} \psi_{02}^{(0)} \right] = -\frac 1 8 \varepsilon \sqrt 2 (\psi_{20}^{(0)} - \psi_{02}^{(0)}).$$ However, $$\begin{split} \psi_{00}(X,\rho) &= \left( \sqrt{\frac{m\omega}{\pi\hbar}} + \mathscr O(\varepsilon^2) \right) \exp \left\{ - \frac{m\omega}{2\hbar} \left[X^2\left(1+\frac 1 2 \varepsilon\right) + \rho^2\left(1-\frac 1 2 \varepsilon\right) + \mathscr O(\varepsilon^2) \right] \right\} \\ &= \sqrt{\frac{m\omega}{\pi\hbar}} \exp \left\{ - \frac{m\omega}{2\hbar} \left(X^2 + \rho^2 \right) \right\} \exp \left\{ - \frac{m\omega}{4\hbar} \varepsilon (X^2-\rho^2) \right\} + \mathscr O(\varepsilon^2) \\ &= \psi_{00}^{(0)}(X,\rho) \left\{1 - \frac{m\omega}{4\hbar} \varepsilon (X^2-\rho^2) + \mathscr O(\varepsilon^2) \right\} + \mathscr O(\varepsilon^2); \end{split}$$ hence the correction must be a linear combination of the unperturbed states $$\alpha_{20}\psi_{20}^{(0)}+ \alpha_{02}\psi_{02}^{(0)}$$ (which agrees with the perturbative result), in such a way that (substituting those states explicitly) $$- \frac{m\omega}{4\hbar} \varepsilon (X^2-\rho^2) = \alpha_{20} X^2 \frac{m\omega}{\hbar} \sqrt{8} - \frac{1}{\sqrt 2}\alpha_{20} + \alpha_{02} \rho^2 \frac{m\omega}{\hbar} \frac{1}{\sqrt 2} - \frac{1}{\sqrt 2}\alpha_{02}.$$ But comparing coefficients leads to an (overdetermined) system of equations for the $$\alpha$$s, giving $$\alpha_{20}= - \frac{1}{8\sqrt 2}\varepsilon, \qquad \alpha_{02} = \frac{1}{2\sqrt 2} \varepsilon,$$ which disagrees with the perturbative prediction (apart from being inconsistent with the requirement $$\alpha_{20} = -\alpha_{02}$$ coming from the constant terms above). Where have I gone wrong?

The problem is that you haven't used a consistent normalization of $$X$$ and $$\rho$$. To see this explicitly, set $$\epsilon = 0$$. Then by symmetry, the ground state should only depend on the combination $$x_1^2 + x_2^2$$. But your ground state depends only on $$X^2 + \rho^2$$. This would have been correct, had you defined $$X = \frac{x_1 + x_2}{2}, \quad \rho = \frac{x_1 - x_2}{2}$$ but you instead defined one of them with a $$1/2$$ and one of them without it. This missing factor of $$2$$ squares into a factor of $$4$$, which is exactly the factor of $$4$$ you're missing at the end.
As a aside, it's generally a good idea to set $$m = \omega = \hbar = 1$$ at the beginning of calculations like these. The simpler form makes it much easier to spot issues, and you can always restore the factors with dimensional analysis.
• I didn't define them this way, they were handed to me in the problem statement. Moreover, with that definition, $X$ is the center of mass of the two particles and $\rho$ is their relative displacement, whereas your definition (or the same with $\sqrt 2$ in the denominator, which simplifies calculations even further) does not have an equally direct physical meaning. Jun 29, 2020 at 7:13