Bessel
Bessel Functions
Bessel Functions of integer and fractional order, of first and second kind, \(J_{\nu}\) and \(Y_{\nu}\), and Modified Bessel functions (of first and third kind), \(I_{\nu}\) and \(K_{\nu}\).
- Keywords
- math
Usage
besselI(x, nu, expon.scaled = FALSE)
besselK(x, nu, expon.scaled = FALSE)
besselJ(x, nu)
besselY(x, nu)Arguments
- x
numeric, \(\ge 0\).
- nu
numeric; The order (maybe fractional and negative) of the corresponding Bessel function.
- expon.scaled
logical; if
TRUE, the results are exponentially scaled in order to avoid overflow (\(I_{\nu}\)) or underflow (\(K_{\nu}\)), respectively.
Details
If expon.scaled = TRUE, \(e^{-x} I_{\nu}(x)\),
or \(e^{x} K_{\nu}(x)\) are returned.
For \(\nu < 0\), formulae 9.1.2 and 9.6.2 from Abramowitz &
Stegun are applied (which is probably suboptimal), except for
besselK which is symmetric in nu.
The current algorithms will give warnings about accuracy loss for
large arguments. In some cases, these warnings are exaggerated, and
the precision is perfect. For large nu, say in the order of
millions, the current algorithms are rarely useful.
Value
Numeric vector with the (scaled, if expon.scaled = TRUE)
values of the corresponding Bessel function.
The length of the result is the maximum of the lengths of the parameters. All parameters are recycled to that length.
References
Abramowitz, M. and Stegun, I. A. (1972). Handbook of Mathematical Functions. Dover, New York; Chapter 9: Bessel Functions of Integer Order.
In order of “Source” citation above:
Sockne, David J. (1973). Bessel Functions of Real Argument and Integer Order. Journal of Research of the National Bureau of Standards, 77B, 125--132.
Cody, William J. (1983). Algorithm 597: Sequence of modified Bessel functions of the first kind. ACM Transactions on Mathematical Software, 9(2), 242--245. 10.1145/357456.357462.
Campbell, J.B. (1980). On Temme's algorithm for the modified Bessel function of the third kind. ACM Transactions on Mathematical Software, 6(4), 581--586. 10.1145/355921.355928.
Campbell, J.B. (1979). Bessel functions J_nu(x) and Y_nu(x) of float order and float argument. Computer Physics Communications, 18, 133--142. 10.1016/0010-4655(79)90030-4.
Temme, Nico M. (1976). On the numerical evaluation of the ordinary Bessel function of the second kind. Journal of Computational Physics, 21, 343--350. 10.1016/0021-9991(76)90032-2.
See Also
Other special mathematical functions, such as
gamma, \(\Gamma(x)\), and beta,
\(B(x)\).
Examples
library(base)
# NOT RUN {
require(graphics)
nus <- c(0:5, 10, 20)
x <- seq(0, 4, length.out = 501)
plot(x, x, ylim = c(0, 6), ylab = "", type = "n",
main = "Bessel Functions I_nu(x)")
for(nu in nus) lines(x, besselI(x, nu = nu), col = nu + 2)
legend(0, 6, legend = paste("nu=", nus), col = nus + 2, lwd = 1)
x <- seq(0, 40, length.out = 801); yl <- c(-.5, 1)
plot(x, x, ylim = yl, ylab = "", type = "n",
main = "Bessel Functions J_nu(x)")
abline(h=0, v=0, lty=3)
for(nu in nus) lines(x, besselJ(x, nu = nu), col = nu + 2)
legend("topright", legend = paste("nu=", nus), col = nus + 2, lwd = 1, bty="n")
## Negative nu's --------------------------------------------------
xx <- 2:7
nu <- seq(-10, 9, length.out = 2001)
## --- I() --- --- --- ---
matplot(nu, t(outer(xx, nu, besselI)), type = "l", ylim = c(-50, 200),
main = expression(paste("Bessel ", I[nu](x), " for fixed ", x,
", as ", f(nu))),
xlab = expression(nu))
abline(v = 0, col = "light gray", lty = 3)
legend(5, 200, legend = paste("x=", xx), col=seq(xx), lty=1:5)
## --- J() --- --- --- ---
bJ <- t(outer(xx, nu, besselJ))
matplot(nu, bJ, type = "l", ylim = c(-500, 200),
xlab = quote(nu), ylab = quote(J[nu](x)),
main = expression(paste("Bessel ", J[nu](x), " for fixed ", x)))
abline(v = 0, col = "light gray", lty = 3)
legend("topright", legend = paste("x=", xx), col=seq(xx), lty=1:5)
## ZOOM into right part:
matplot(nu[nu > -2], bJ[nu > -2,], type = "l",
xlab = quote(nu), ylab = quote(J[nu](x)),
main = expression(paste("Bessel ", J[nu](x), " for fixed ", x)))
abline(h=0, v = 0, col = "gray60", lty = 3)
legend("topright", legend = paste("x=", xx), col=seq(xx), lty=1:5)
##--------------- x --> 0 -----------------------------
x0 <- 2^seq(-16, 5, length.out=256)
plot(range(x0), c(1e-40, 1), log = "xy", xlab = "x", ylab = "", type = "n",
main = "Bessel Functions J_nu(x) near 0\n log - log scale") ; axis(2, at=1)
for(nu in sort(c(nus, nus+0.5)))
lines(x0, besselJ(x0, nu = nu), col = nu + 2, lty= 1+ (nu%%1 > 0))
legend("right", legend = paste("nu=", paste(nus, nus+0.5, sep=", ")),
col = nus + 2, lwd = 1, bty="n")
x0 <- 2^seq(-10, 8, length.out=256)
plot(range(x0), 10^c(-100, 80), log = "xy", xlab = "x", ylab = "", type = "n",
main = "Bessel Functions K_nu(x) near 0\n log - log scale") ; axis(2, at=1)
for(nu in sort(c(nus, nus+0.5)))
lines(x0, besselK(x0, nu = nu), col = nu + 2, lty= 1+ (nu%%1 > 0))
legend("topright", legend = paste("nu=", paste(nus, nus + 0.5, sep = ", ")),
col = nus + 2, lwd = 1, bty="n")
x <- x[x > 0]
plot(x, x, ylim = c(1e-18, 1e11), log = "y", ylab = "", type = "n",
main = "Bessel Functions K_nu(x)"); axis(2, at=1)
for(nu in nus) lines(x, besselK(x, nu = nu), col = nu + 2)
legend(0, 1e-5, legend=paste("nu=", nus), col = nus + 2, lwd = 1)
yl <- c(-1.6, .6)
plot(x, x, ylim = yl, ylab = "", type = "n",
main = "Bessel Functions Y_nu(x)")
for(nu in nus){
xx <- x[x > .6*nu]
lines(xx, besselY(xx, nu=nu), col = nu+2)
}
legend(25, -.5, legend = paste("nu=", nus), col = nus+2, lwd = 1)
## negative nu in bessel_Y -- was bogus for a long time
curve(besselY(x, -0.1), 0, 10, ylim = c(-3,1), ylab = "")
for(nu in c(seq(-0.2, -2, by = -0.1)))
curve(besselY(x, nu), add = TRUE)
title(expression(besselY(x, nu) * " " *
{nu == list(-0.1, -0.2, ..., -2)}))
# }