I'm trying to calculate a function that includes Bessel functions and am running into precision issues, generally due to a small difference between large numbers.
The function is
$$F=\frac{J_0(K\phi)Y_1(\phi)-J_1(\phi)Y_0(K\phi)}{J_0(\phi)Y_1(\phi)-J_1(\phi)Y_0(\phi)} $$
where $J$ and $Y$ are the first and second Bessel functions. The overall function is on the order of unity, but both numerator and denominator tend to zero for the values I'm interested in. I've tried various rearrangements (including using a normalized Bessel function) but to no avail. I feel there should be some simplification that can be made to help with precision.
As an example, for k=1.05 and $\phi$= -40 000 + 20i, with double precision I get INF because the denominator is zero, whereas when using quadruple precision I get a reasonable results of -0.5533 + 1.0666i.
Any ideas?
edit: Thanks Paul Enta for the Wronskian link (denominator = $-2/(\pi\phi)$). I think that solves half of the problem. However, the numerator is still problematic. If I calculate it as shown above I get 0 (quad precision gives -8.8151e-06 + 1.6972e-05i). If I rewrite it in terms of ratios: $$J_0(k\phi)Y_1(\phi)\left[1-\frac{J_1(\phi)}{Y_1(\phi)}\frac{Y_0(k\phi)}{J_0(k\phi)}\right]$$ then I also get poor results as the ratios are close to unity.
Evaluating the numerator for small values of $\phi$ and $k\phi$ should not pose numerical problems. For large values, we can use the asymptotic expansions of $J_\nu$ and $Y_\nu$. After simplification of the trigonometric expressions and by keeping two terms, one obtains of the following asymptotic expansion \begin{equation} J_0(K\phi)Y_1(\phi)-J_1(\phi)Y_0(K\phi)\simeq -\frac{2}{\pi\phi\sqrt{K}}\left[\cos\left( \phi(K-1) \right)+(3K+1)\frac{\sin\left( \phi(K-1) \right)}{8\phi K}\right] \end{equation} Then, for large arguments, taking into account the simplification of the denominator which is a Wronskian, \begin{equation} F\simeq\frac{1}{\sqrt{K}}\left[\cos\left( \phi(K-1) \right)+(3K+1)\frac{\sin\left( \phi(K-1) \right)}{8\phi K}\right] \end{equation}