I have the following system of coupled differential equations (the Bogoliubov-de Gennes equations for a certain geometry):
$$ \begin{cases} -\frac{1}{2}u''(r) + \big[\frac{1}{2r^2} + f_1(r)\big]u(r) + f_2(r)v(r)=\omega u(r)\\ -\frac{1}{2}v''(r) + \big[\frac{1}{2r^2} + f_1(r)\big]v(r) + f_2(r)u(r)=-\omega v(r), \end{cases} $$ where $f_i(r)$ are known functions which go to zero for large $r$ and are finite for small $r$. This is in fact an eigenvalue problem, because I need to find the value of $\omega$ for which this system has (normalizable) solution.
Now I would like to solve this system numerically: hence find the value of $\omega$ such that the system has a solution and then also the functions $u(r),v(r)$. I am familiar with numerical recipes to solve coupled differential equations, but now that this is an eigenvalue problem as well, I am at a loss. Could anyone give me some pointers on how to tackle such a problem? Some references where this is described would be appreciated as well!
