Solve analytically the time-independent Schrödinger equation for the Hamiltonian $\hat H$ The Hamiltonian $\hat H$ is given by:
$$\hat H=\frac{1}{2}(\hat p^2+\hat q^2)+K(\hat p\hat q+\hat q\hat p)$$
where $K$ is a real constant, in the coordinate space:
$$\hat p=-iℏ\frac{\partial}{\partial q}$$
$$\hat q=q$$
I try to solve it by expanding those operater, and end up with:
$$\hat H=-\frac{ℏ^2}{2}\frac{\partial^2}{\partial q^2}+\frac{q^2}{2}+K[(-iℏ)\frac{\partial}{\partial q}q+(-iℏ)q\frac{\partial}{\partial q}]$$
$$\hat H=-\frac{ℏ^2}{2}\frac{\partial^2}{\partial q^2}+\frac{q^2}{2}+K(-iℏ)[(1+q\frac{\partial}{\partial q})+q\frac{\partial}{\partial q}]$$
$$ \hat H=-\frac{ℏ^2}{2}\frac{\partial^2}{\partial q^2}-i(2Kℏ)q\frac{\partial}{\partial q}+(\frac{q^2}{2}-iKℏ)$$
However I don't know how to solve the second order differiential equation:
$$\hat H\psi=E\psi$$
$$-\frac{ℏ^2}{2}\frac{\partial^2\psi}{\partial q^2}-i(2Kℏ)q\frac{\partial\psi}{\partial q}+(\frac{q^2}{2}-iKℏ)\psi=E\psi$$
$$-\frac{ℏ^2}{2}\frac{\partial^2\psi}{\partial q^2}-i(2Kℏ)q\frac{\partial\psi}{\partial q}+(\frac{q^2}{2}-iKℏ-E)\psi=0$$
$$\frac{\partial^2\psi}{\partial q^2}+i\frac{4K}{ℏ}q\frac{\partial\psi}{\partial q}+\biggl(\frac{2E}{ℏ^2}+i\frac{2K}{ℏ}-\frac{q^2}{ℏ^2}\biggl)\psi=0$$
Can someone give me some hints? Or correct me if I am using a wrong approach.
 A: It seems that it is quite useful to notice that the last term in the Hamiltonian reduces,
$$K(pq+qp)=K\left[pq+(i+pq)\right]=2Kpq+iK,$$
with help of $[q,p]=i$ with $\hbar=1$. Then, we have
$$H=\frac{1}{2}(p^2+q^2)+2Kpq+iK.$$
The first terms can be rewritten with help of ladder operators,
$$p=i\sqrt{\frac{\omega}{2}}(a^{\dagger}-a),\quad q=\sqrt{\frac{1}{2\omega}}(a^{\dagger}+a).$$
Next, it seems that in terms of ladder operators the first terms simply gives
$$\frac{1}{2}\left(\omega+\frac{1}{\omega}\right)\left(a^{\dagger}a+\frac{1}{2}\right),$$
whereas the second term gives
$$\frac{iK}{2}(-aa^{\dagger}+a^{\dagger}a)=-\frac{iK}{2},$$
where I tacitly assume that you are interested in averages $\langle n|a^{\dagger}a|n\rangle$. Finally, we have
$$\frac{1}{2}\left(\omega+\frac{1}{\omega}\right)\left(a^{\dagger}a+\frac{1}{2}\right)+\frac{iK}{2}$$
and it looks like the QHO with constant energy shift.
However, I am not sure that everything is correct. The differential operator that you have write down is so close to the differential operator in the Landau levels problem and it seems that this question can help
A: A Bogoliubov transformation to new operators $\hat b$ and $\hat b^\dagger$
\begin{align}
\hat b= u\hat a+v\hat a^\dagger\, ,\qquad \hat b^\dagger u^*\hat a^\dagger +v^*\hat a
\end{align}
with $\vert u\vert^2-\vert v\vert^2=1$ and suitably chosen $u$ and $v$ will bring $H$ to a diagonal form.
Alternatively, this Hamiltonian is expressible in terms of generators of the $\mathfrak{su}(1,1)$:
\begin{align}
 \hat K_+&=\frac{1}{2}\hat a^\dagger\hat a^\dagger \qquad 
        \hat K_-=\frac{1}{2}\hat a\hat a\qquad 
        \hat K_0=\frac{1}{4}(\hat a\hat a^\dagger+\hat a^\dagger\hat a)
\end{align}
and can be diagonalized $UHU^\dagger=H_d\sim K_0$ by a transformation of the form $U=\exp(-i\Lambda K_x)$ for suitable $\Lambda$.  If you need to find the eigenkets, it is then a matter of finding the action of $U$ on the kets $\vert n\rangle$: this can be found in

Ui, H., 1970. Clebsch-Gordan formulas of the SU (1, 1) group. Progress of Theoretical Physics, 44(3), pp.689-702,

or in

Perelomov, A., 2012. Generalized coherent states and their applications. Springer Science & Business Media.

A: This is not a solution, but multiple suggestions for attacking the problem:
Brutal force
The equation is likely reducible to a hypergeometric equation, perhaps even something simpler, like parabolic cylinder functions or confluent hypergeometric. The way to start is to try to reduce it to one of the canonical forms by writing ψ(q)=f(q)ϕ(q) and choosing f(q) in such a way as to have a constant coefficient in the last term.
Analogy with magnetic field
Another option is re-writing it as a particle in a constant magnetic field and parabolic potential... then you change the gauge, reducing it to a harmonic oscillator... and then you transform back to the original gauge.
Actually, it is an oscillator...
One could even try diagonalize the quadratic form in the Hamiltonian, introducing new position and momentum operators.
