Bessel {base} | R Documentation |

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)*.

besselI(x, nu, expon.scaled = FALSE) besselK(x, nu, expon.scaled = FALSE) besselJ(x, nu) besselY(x, nu)

`x` |
numeric, |

`nu` |
numeric; The |

`expon.scaled` |
logical; if |

If `expon.scaled = TRUE`

, *exp(-x) I(x;nu)*,
or *exp(x) K(x;nu)* 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.

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.

Original Fortran code:
W. J. Cody, Argonne National Laboratory

Translation to C and adaptation to **R**:
Martin Maechler maechler@stat.math.ethz.ch.

The C code is a translation of Fortran routines from http://www.netlib.org/specfun/ribesl, ../rjbesl, etc. The four source code files for bessel[IJKY] each contain a paragraph “Acknowledgement” and “References”, a short summary of which is

- besselI
based on (code) by David J. Sookne, see Sookne (1973)... Modifications... An earlier version was published in Cody (1983).

- besselJ
as

`besselI`

- besselK
based on (code) by J. B. Campbell (1980)... Modifications...

- besselY
draws heavily on Temme's Algol program for

*Y*... and on Campbell's programs for*Y_ν(x)*.... ... heavily modified.

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.
doi: 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.
doi: 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.
doi: 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.
doi: 10.1016/0021-9991(76)90032-2.

Other special mathematical functions, such as
`gamma`

, *Γ(x)*, and `beta`

,
*B(x)*.

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(-.8, .8) plot(x, x, ylim = yl, ylab = "", type = "n", main = "Bessel Functions J_nu(x)") for(nu in nus) lines(x, besselJ(x, nu = nu), col = nu + 2) legend(32, -.18, legend = paste("nu=", nus), col = nus + 2, lwd = 1) ## Negative nu's : xx <- 2:7 nu <- seq(-10, 9, length.out = 2001) op <- par(lab = c(16, 5, 7)) 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=seq(xx)) par(op) x0 <- 2^(-20:10) plot(x0, x0^-8, log = "xy", ylab = "", type = "n", main = "Bessel Functions J_nu(x) near 0\n log - log scale") for(nu in sort(c(nus, nus+0.5))) lines(x0, besselJ(x0, nu = nu), col = nu + 2) legend(3, 1e50, legend = paste("nu=", paste(nus, nus+0.5, sep=",")), col = nus + 2, lwd = 1) plot(x0, x0^-8, log = "xy", ylab = "", type = "n", main = "Bessel Functions K_nu(x) near 0\n log - log scale") for(nu in sort(c(nus, nus+0.5))) lines(x0, besselK(x0, nu = nu), col = nu + 2) legend(3, 1e50, legend = paste("nu=", paste(nus, nus + 0.5, sep = ",")), col = nus + 2, lwd = 1) x <- x[x > 0] plot(x, x, ylim = c(1e-18, 1e11), log = "y", ylab = "", type = "n", main = "Bessel Functions K_nu(x)") 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)}))

[Package *base* version 3.5.0 Index]