I have two coupled quantum harmonic oscillators given by the following Hamiltonian:

$$H=\frac{p_{x}^{2}}{2}+\frac{\omega^{2} x^{2}}{2}+\frac{p_{y}^{2}}{2}+\frac{\Omega^{2} y^{2}}{2}+\frac{C p_{x} y}{2}.$$

As you can see, the coupling is done over the position variable of one oscillator and the momentum of the other. I need to find the wavefunction of states where either of the two oscillators (or both) is excited.

What I tried to do: In general when I have couple harmonic oscillators where the coupling term is of the form $C (x_1-x_2)^2$ I start by diagonalising the Hamiltonian, then define normal coordinates and write down the Hamiltonian using those. I.e., I write the Hamiltonian under the form

$$H=\frac{1}{2} \sum_{i=1}^{2} p_{i}^{2}+\frac{1}{2} \sum_{i, j=1}^{2} x_{i} K_{i j} x_{j}.$$

Then I diagonalise $K$,

$$K_D = UKU^T.$$

And define normal coordinates

$$\left(\begin{array}{l} x_{+} \\ x_{-} \end{array}\right) \equiv U\left(\begin{array}{l} x_{1} \\ x_{2} \end{array}\right),$$

leading to

$$H=\frac{1}{2}\left[p_{+}^{2}+p_{-}^{2}+\omega_{+}^{2} x_{+}^{2}+\omega_{-}^{2} x_{-}^{2}\right].$$

I am unable to use this here since I do not know how to diagonalise this Hamiltonian with the term $C p_x y$ present.

  • $\begingroup$ Since you have posted a bounty asking for a more detailed answer, could you clarify what about Mike Stone's answer is insufficient? It looks rather comprehensive to me. $\endgroup$ – Richard Myers Jun 27 at 4:43
  • $\begingroup$ @RichardMyers MS’s answer is great but there are a few details in it that I do not understand how to implement exactly. For instance how to I get $R$? How exactly do I compute $H^{-1/2}$? $\endgroup$ – Y2H Jun 27 at 7:22
  • $\begingroup$ @RichardMyers another issue is that, when I define $D$ as Mike Stones notes and assume no constraints on $S$ at all, Mathematica simply shows that the equation $S^T H S$ cannot be satisfied. $\endgroup$ – Y2H Jun 27 at 12:10

Your problem is basically that of a charged particle in a magnetic field in the Landau gauge, but I'll explain the general theory for such systems as it's both interesting and suprisingly complicated --- and seldom found in textbooks. It's essentially the theory of Bose Bogoluibov transformations.

Suppose we are given classical (or quantum) quadratic Hamiltonian $$ {\mathcal H}[{\bf p},{\bf x}]= \frac 12 M_{ij} p_ip_j+ \frac 12 V_{ij} x_ix_j + K_{ij}\,p_i x_j \nonumber\\ =\frac 12 \left[\matrix{{\bf p}^T& {\bf x}^T}\right] H \left[\matrix{ {\bf p}\cr {\bf x}}\right], \nonumber $$ where and $H$ is the real, positive definite, symmetric $2N$-by-$2N$ matrix $$ H= \left[\matrix{M& K\cr K^T &V}\right]. $$ (If $H$ is not positive definite the diagonalization is not possible as the energy of the system is not bounded below and the eigenfrequencies are pure imaginary.)

We seek a transformation $$ \left[\matrix{ {\bf p}\cr {\bf x}}\right]= S \left[\matrix{ {\bf P}\cr {\bf X}}\right] $$ that diagonalizes $H$ via the congruence transformation ${H}\to S^T HS $ and preserves the classical Poisson brackets and/or the quantum commutation relations. This preservation requires $$ S^T JS=J, \quad J= \left[\matrix{0&- {\mathbb I}_n \cr {\mathbb I}_n&\phantom{-} 0}\right], $$ so $S\in {\rm Sp}[2N, {\mathbb R}]$, the non-compact symplectic group of linear canonical transformations. That positive definiteness of $H$ is sufficient to ensure that we can find such an $S$ is the statement of Williamson's theorem. The resulting frequencies are called the symplectic spectrum of $H$. They are unrelated to the eigenvalues of $H$. Further matrix $S$ is neither orthogonal nor unitary, and so does not obey $S^{-1}=S^\dagger$

We observe that $H$ being positive definite ensures that the matrices $H^{\pm 1/2}$ are well defined. We then see that
the matrix $$ \tilde J= H^{-1/2}JH^{-1/2} $$ is skew symmetric and hence an element of $\mathfrak{so}[2N]$. Therefore there exists an $R\in {\rm O}(2N)$ that conjugates $\tilde J$ into the Cartan algebra of $\mathfrak{so}(2n)$ --- i.e. $$ R^TH^{-1/2} J H^{-1/2} R= \left[\matrix{0&- \Omega^{-1} \cr \Omega^{-1} &0}\right], $$ with $\Omega^{-1}$ diagonal. As we are allowing $R\in {\rm O}[2N]$ rather than demanding $ R\in {\rm SO}[2N]$, we can ensure that $\Omega^{-1}$ possesses strictly positive entries. This being so, we define
$$ D= \left[\matrix{\Omega^{1/2}&0 \cr 0& \Omega^{1/2}}\right], $$ so that $$ DR^T H^{-1/2} J H^{-1/2} RD=J,\nonumber\\ DR^T H^{-1/2}H H^{-1/2} R D= D^2.\nonumber $$ Thus $$ S= H^{-1/2} RD \in {\rm Sp}[2N, {\mathbb R}], $$ and $$ S^T HS= D^2 = \left[\matrix{\Omega&0 \cr 0& \Omega}\right] $$ where $\Omega= {\rm diag}(\omega_1,\ldots \omega_N)$. We have $$ {\mathcal H}[{\bf P}, {\bf Q}]=\frac 12 \left[\matrix{ {\bf P}^T& {\bf X}^T}\right]\left[\matrix{ \Omega &0\cr 0&\Omega}\right]\left[\matrix{ {\bf P}\cr {\bf X}}\right]= \sum_i \frac {\omega_i}{2} (P_i^2+X_i^2), $$ which is a set of decoupled harmonic oscillators with frequency $\omega_i$.

We can of course scale $X_i\to \sqrt{\omega} X_i$ and $P_i\to P_i/\sqrt{\omega}$ without changing the commutation relations, and so
$$ \sum_i \frac {\omega_i}{2} (P_i^2+X_i^2)\to \sum_i \frac {1}{2} (P_i^2+\omega_i^2X_i^2), $$ which probably looks more familiar.

The eigenvalues $\omega_i$ are most easily found from those of $JH$ which are $\pm i\omega_i$.

Note added: I see that this same strategy is discussed in a math stack exchange answer: https://math.stackexchange.com/questions/1171842/finding-the-symplectic-matrix-in-williamsons-theorem

  • $\begingroup$ Never seen such computations. Very cool! +1 $\endgroup$ – Davide Morgante Jun 13 at 13:43
  • $\begingroup$ I am unsure If one could find a symmetric K matrix such that the proposed form of the Hamiltonian matches the question's Hamiltonian since there is only a $p_x y$ term. Also, is it 100% sure that the $p_k$ in the first equation should not be a $p_j$ ? Aside from that, thank you very much for posting this since it looks like it works beautifully when there are xipj and xjpi terms. $\endgroup$ – elscan Jun 13 at 17:31
  • $\begingroup$ Thanks for catching typo (now fixed)! I think that you have $H=(1/2)(p_x+\lambda y)^2- (1/2)\lambda^2 y^2 +\ldots $ which is why I was suggesting that you regard $\lambda y=A_x$ as a gauge field. YOu can also add a $xp_y-yp_x$ angular momentum term to make symmetric but may nt help much. $\endgroup$ – mike stone Jun 13 at 18:19
  • 1
    $\begingroup$ I don't know of one. I chased this stuff down in papers. In particular Simon, Chaturvedi, Srinivasan "Congruences and Canonical Forms for a Positive Matrix: Application to the Schweinler-Wigner Extremum Principle" arXiv:math-ph/9811003 $\endgroup$ – mike stone Jun 14 at 11:49
  • 1
    $\begingroup$ Juts look at my notes oin this: people.physics.illinois.edu/stone/bose_bogoliubov.pdf $\endgroup$ – mike stone Jun 26 at 13:10

I am liking the answers which have been posted, but figure to offer an answer which uses more primitive means. This system has a classical analog with Lagrangian of the form:

$$ L=\frac{1}{2}\cdot \dot{x}\cdot M\cdot \dot{x}-\frac{1}{2}\cdot x\cdot K\cdot x\mp \frac{1}{2}\cdot \dot{x}\cdot C\cdot x $$ I say $\mp \frac{1}{2}\cdot \dot{x}\cdot C\cdot x$ since I am not 100 % certain that the $p_x y$ term is to be considered momentum - potential energy as is the case in electrodynamics.

Lagrange' s equation $\frac{d}{dt}\cdot \frac{\partial L}{\partial \dot{x}}=\frac{\partial L}{\partial x}$ produces:

$$ \left( M^T+M \right) \cdot \ddot{x} \mp C\cdot \dot{x}=\mp C^T\cdot \dot{x}-\left(K^T+K\right)\cdot x $$

Which produces an eigenvalue equation and I assume that all who are reading this know how to find the normal modes of this system such that the Lagrangian becomes:

$$ L=\frac{1}{2} \underset{\alpha }{\Sigma } \left(\dot{Q}_{\alpha }^2-\omega _{\alpha }^2 Q_{\alpha }^2\right) $$

where the $\left\{\overset{\to }{Q}_{\alpha }\right\}$ and $\left\{\omega _{\alpha }\right\}$ are the eigenmodes and eigenfrequencies. And so the classical Hamiltonian is:

$$ H=\frac{1}{2} \underset{\alpha }{\Sigma } \left(p_{\alpha }^2+\omega _{\alpha }^2 Q_{\alpha }^2\right) $$

And the only step left is to use canonical quantization :

$$ Q_{\alpha }\to \hat{Q}_{\alpha },p_{\alpha }\to \hat{p}_{\alpha }=-i\hbar \cdot \frac{\partial }{\partial Q_{\alpha }} $$

where I use the hats to indicate that $\hat{Q}$ and $\hat{p}$ are quantum operators.


I started off by assuming a general transformation of the form:

$$\left\{\begin{array}{l} p_{x}=v p_{-}+c p_{+}+d x_{+}, \\ p_{y}=k p_{-}+l p_{+}, \\ x=n x_{+}+q x_{-}, \\ y=u x_{+}+ax_{-}+b p_{+} , \end{array}\right.$$

with ${p_+, p_-,x_+, x_-}$ the new variables. To avoid confusion, I am also renaming the coupling coefficient $C$ to $\gamma$. Next I plugged the above transformations in the Hamiltonian and demanded that no off-diagonal terms survive. This constraints the transformation coefficients. I also demand that the coefficients next to $p_+^2 $ and $p_-^2$ in the Hamiltonian be $1$. This leads to the set of transformations

$$\left\{\begin{array}{l} p_{x}=-\sqrt{1-l^{2}} p_{+}+l p_{m}-\frac{a \gamma}{2} x_{-}, \\ p_{y}=k p_{+}+\sqrt{1-l^{2}} p_{-}, \\ x=n x_{+}, \\ y=a x_{-} . \end{array}\right.$$

We are now left with the set of coefficients ${a, l, n}$ to be determined. This can be done by enforcing the canonical commutation relations on the new and old sets of coordinates. We then find

$$[x, p_y]=0 \Rightarrow l=0, \\ [x, p_x]=i \hbar \Rightarrow n=-1, \\ [y, p_y] = i \hbar \Rightarrow a=1.$$

The resulting set of transformations is hence

$$\left\{\begin{array}{l} p_{x}=-p_{+}- \frac{\gamma}{ 2} x . \\ p_{y}=p_{-}, \\ x=-x_{+}, \\ y=x . \end{array}\right.$$

leading to the Hamiltonian

$$H = p_{+}^{2}+p_{-}^{2}+\omega^{2} x_{+}^{2}+\left(\Omega^{2}-\gamma^{2} / 4\right) x_{-}^{2}.$$

The normal modes are obviously

$$\omega_+ = \pm \omega, \,\,\,\,\,\,\,\, \omega_- = \sqrt{\left(\Omega^{2}-\gamma^{2} / 4\right)}.$$


Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.