Quantum scattering problem in anisotropic medium Problem

For simplicity, let's consider a two dimensional problem. A plane wave
  is incident from $x$ axis: $$ \psi_{inc}=e^{ikx} $$ The 2D Hamiltonian
  reads: $$ H=\hat{p}_x^2+\alpha \hat{p}_y^2 +V(r) $$ where $\alpha>0$
  but $\alpha\neq1$. $V(r)=+\infty$  when $x^2+y^2<r_0^2$ and $V(r)=0$
  otherwise.

Discussion
In the standard quantum mechanics book, it seems we always have a spherically symmetric Hamiltonian, therefore we use the partial wave decomposition and different partial waves scatter independently. In the two dimensional problem, the radial wave function is the Hankel function of the first kind because we require the scattering wave to have the asymptotic form $\sim e^{ikr}/r$.
In the present case, should we also require the scattering wave to have the asymptotic form:
$$
\psi_{ref}\sim \frac{e^{ikr}}{r}
$$
If we make a change of variables $x\to x', y\to \sqrt{\alpha}y'$, then the scattering problem becomes:
$$
H=\hat{p}_{x'}^2+ \hat{p}_{y'}^2 +V(x',y')
$$
where $V(x',y')=+\infty$  when $x'^2+\alpha y'^2<r_0^2$ and $V(x',y')=0$ otherwise. From this, maybe we should require the boundary condition to be:
$$
\psi_{ref}\sim \frac{e^{ikr'}}{r'}
$$
Even if we set the boundary condition, then what is the next step towards the wave function in whole space?
Question
In summary: how to solve this seemingly simple scattering problem, what is the boundary condition? How to get the near field ($r\sim r_0$) wave function? Analytical and numerical methods are both welcomed.
 A: After the change of variables, the problem is equivalent to the solution of the 2D Helmholtz equation. $$\nabla'^2\psi + A^2\psi = 0$$
If the boundary had remained circular, the radial part of the solution would be a Hankel function of the first kind which asymptotically goes as $${H_{n}^{(1)}(k r')\sim {\sqrt {\frac {2}{\pi k r'}}}\exp \left(i\left(k r'-{\frac {n \pi }{2}}-{\frac {\pi }{4}}\right)\right)},$$ so, in the primed coordinates, the solution would be of the form $$ \psi_{ref}\sim \frac{e^{ikr'}}{\sqrt{r'}}.$$ In 2D, the solution decays radially as ${r'}^{-1/2}$, unlike the 3D case, which decays as ${r'}^{-1}$. 
Unfortunately, in your case, the boundary becomes elliptic. Although the asymptotic behavior should not change, the actual solution is a bit more complicated. From now on I will be considering the problem expressed in the new coordinates, so I am dropping the primes, for the sake of brevity. I will also assume, with no loss of generality, that $\alpha>1$, so that the boundary forms an ellipse with focal points at $(-c,0) $ and $(c,0)$, where
$$ c = r_0 \sqrt{\frac{\alpha-1}{\alpha}} $$ and the ellipticity $$ e = \sqrt\frac{\alpha-1}{\alpha}.$$ 
Such problems are best solved in the elliptic coordinates $(\mu,\theta)$, where
$$
\begin{aligned}
x &= c \cosh \mu \cos \theta \newline
y &= c \sinh \mu \sin \theta.
\end{aligned}
$$
Line $\mu=const$ is a closed ellipse with the same focal points as the boundary, with $c \cosh \mu = r_0$ at the boundary. Line $\cos\theta = cost $ is a family of concentric parabolas. In the new variables, Helmholtz equation takes the form $$
\frac{\partial^2\psi}{\partial \mu^2} + \frac{\partial^2\psi}{\partial\theta^2} + c^2A^2 \left[\cosh^2\mu - cos^2 \theta \right]\psi = 0, 
$$
which takes solutions of the form
$$
\psi_n = M_n^{(1)}(c A,\mu)\left[C c_n(c A,\theta) +D s_n(c A,\theta)\right].
$$
Functions $c_n(c A,\theta)$ and $s_n(c A,\theta)$ are periodic solutions to Mathieu equation. They form an orthogonal basis and for $c\rightarrow 0$ they go to $\cos(n \theta)$ and $\sin(n \theta)$ respectively. The family $M_n^{(1)}(c^2 A,\mu)$ are solutions of the modified Mathieu equation that asymptotically go to $H_n^{(1)}(\sqrt 2 c A e^\mu)$.
The last step is to express the boundary value condition at the ellipse $$\psi_{ref} = - \psi_{inc} = - e^{i k r_0 \cos\theta}$$ as a linear combination of $s_n$ and $c_n$.
The study of Mathieu functions is rather technical and a complete solution to the problem would blow the size of this answer way out of proportion. Instead, I will point to some resources that can provide all tools you need to carry it on yourself.


*

*For a general introduction to Mathieu functions and their relation to Helmholtz equation, see Morse & Feshbach "Methods of theoretical physics",vol. 1 and vol. 2, chapter 5 and 11.

*For more details, but less physics, see Abramowitz and Stegun,
Handbook of Mathematical Functions, chapter 20.

*For the integrals necessary to perform the boundary condition expansion, see Gradshteyn and Ryzhik, "Table of integrals Series and Products", 6.924.


Since all three of them are classic reference books, you shouldn't find it difficult to locate a copy. 
