Infinite annular potential well. Trouble with solving Bessel equation to get eigenstates and energy I have infinite annular potential well (scheme in the picture). 
 
Schrodinger equtation in the anullus (for $R_1 <r<R_2$ is $V=0$) with polar coordinates is
\begin{equation} - \frac{ \hbar }{2m} \left[ \frac{1}{r}\frac{\partial}{\partial r}\left(r   \frac{\partial \psi}{\partial r}\right) + \frac{1}{r^2}\frac{\partial^2 \psi }{\partial \phi^2} \right] = E \psi
\end{equation}
I assume separable solution $\psi = R(r) F(\phi)$. So after modification I get equation
$$
- \frac{ \hbar }{2m} \left[ r \frac{R'}{R} + r^2\frac{R''}{R} \right] -r^2 E = \frac{ \hbar }{2m} \frac{F''}{F} = - \epsilon_{\phi}
$$
From the side with $F$ i got $F(\phi) = e^{i l \phi}$, $\epsilon_{\phi} = \frac{\hbar^2 l^2}{2m}$ and $l$ is integer (from periodicity)
From the side with $R$ I get after modification Bessel equation:
$$
r^2 R'' + r R' +\left( k^2 r^2 - l^2 \right) = 0,
$$
where $k = \sqrt(2mE / \hbar^2)$.
The general solution should consists of Bessel's functions of first kind $J_l(kr)$ and second kind $Y_l(kr)$:
$$
R(r) = A_l J_l(kr) + B_l Y_l(kr)
$$
I have got two boundary conditions $R(R_1) =0$ and $R(R_2) =0$ However I need to find 3 constants: $A_l$, $B_l$ and $k$, to get eigenfunctions $\psi$ and eigenvalues $E$. 
Did I forget some other condition? I can't find out how to continue. 
How to continue with analytical solution?
 A: The problem isn't fully solvable in exact form, since it requires the solution of a transcendental equation to get $k$, and with it the energy eigenvalue. This shouldn't be surprising $-$ it's exactly the same situation as for a particle in a ring (i.e. without the internal cutout).
You are correct in all the steps that lead down to your radial eigenfunction in the form
$$
R(r) = A_l J_l(kr) + B_l Y_l(kr),
$$
and now the task is indeed to solve for the coefficients $A_l$ and $B_l$ and for the wavenumber $k$. From the first two, one is easy, and you can get it by setting $R(r)$ to zero at one of the two borders, say,
$$
R(R_1) = A_l J_l(kR_1) + B_l Y_l(kR_1)=0,
$$
or in other words
$$
B_l = -\frac{ J_l(kR_1)}{Y_l(kR_1)}A_l,
$$
so that's $B_l$ down and $A_l$ to go. Shifting notation slightly by defining $C_l=A_l/Y_l(kR_1)$, you get your full eigenfunction in the form 
$$
R(r) = C_l \bigg[Y_l(kR_1) J_l(kr) - J_l(kR_1) Y_l(kr)\bigg],
$$
and your second boundary condition in the form
$$
R(R_2) = C_l \bigg[Y_l(kR_1) J_l(kR_2) - J_l(kR_1) Y_l(kR_2)\bigg] = 0,
$$
and this is where things get tricky.
The reason that things get hairy at this point is that you are now solving for the wavenumber $k$ (itself a proxy for the energy eigenvalue $E=\hbar^2 k^2/2m$) in the transcendental equation
$$
Y_l(kR_1) J_l(kR_2) - J_l(kR_1) Y_l(kR_2) = 0,
$$
and this simply does not accommodate any exact-form solutions (much in the same way e.g. the finite square well will reduce to a transcendental equation of the form $\tan(ka) = \sqrt{C+k^2b^2}$ with no exact solutions). This is analogous to the situation for the no-inner-cutout case, where you'll just be left with
$$
J_l(kR)=0,
$$
where you know that the product $kR$ needs to be a Bessel zero, but that is as far as the analytical approaches can tell you.
In the no-inner-cutout case, of course, you're in luck: the Bessel zeros are extremely common and well-studied objects, and they are included in most tabulations of special functions (e.g. Abramowitz and Stegun, and the like) and they are a standard component of most mathematical software.
For your case, with a nonzero inner cutout, you're less lucky, because the situation is less common and it is therefore less widespread in both tabulations and mathematical software. The keyword to search for is the cross-product Bessel zeros: we do know a lot about them, but that knowledge doesn't always trickle down to applications.
When I work with those objects, I tend to do so in Mathematica, which does have a BesselJZero in-built function in the core language, but for which the rest of the BesselZeros package, including the cross-product zeros, didn't quite make it into the core system. The package is still available and it does work (though with some reliability problems on the cross-product zeros!) if that's the way you want to roll.

And, finally, that third constant: for the normalization constant $C_l$, which is fixed by requiring that
$$
\int_{R_1}^{R_2} |R(r)|^2 r\,\mathrm dr = 1
$$
or some similar requirement $-$ yeah, there's no chance of that being analytically integrable, much like the no-inner-cutout case. You just integrate it numerically when it comes down to it.
