I am trying to calculate the logarithm of a modified Bessel function of second type in MATLAB, i.e. something like that:
log(besselk(nu, Z))
where e.g.
nu = 750;
Z = 1;
I have a problem because the value of log(besselk(nu, Z)) goes to infinity, because besselk(nu, Z) is infinity. However, log(besselk(nu, Z)) should be small indeed.
I am trying to write something like
f = double(sym('ln(besselk(double(nu), double(Z)))'));
However, I get the following error:
Error using mupadmex Error in MuPAD command: DOUBLE cannot convert the input expression into a double array. If the input expression contains a symbolic variable, use the VPA function instead.
Error in sym/double (line 514) Xstr = mupadmex('symobj::double', S.s, 0)`;
How can I avoid this error?
As njuffa pointed out, DLMF gives asymptotic expansions of K_nu(z) for large nu. From 10.41.2 we find for real positive arguments z:
besselk(nu,z) ~ sqrt(pi/(2nu)) (e z/(2nu))^-nu
which gives after some simplification
log( besselk(nu,z) ) ~ 1/2*log(pi) + (nu-1/2)*log(2nu) - nu(1 + log(z))
So it is O(nu log(nu)). No surprise the direct calculation fails for nu > 750.
I dont know how accurate this approximation is. Perhaps you can compare it for the values where besselk is smaller than the numerical infinity, to see if it fits your purpose?
EDIT: I just tried for nu=750 and z=1: The above approximation gives 4.7318e+03, while with the result of horchler we get log(1.02*10^2055) = 2055*log(10) + log(1.02) = 4.7318e+03. So it is correct to at least 5 significant digits, for nu >= 750 and z=1! If this is good enough for you this will be much faster than symbolic math.