R: Bessel Functions

Bessel {base}

R Documentation

Bessel Functions

Description

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}.

gammaCody is the (\Gamma) function as from the Specfun package and originally used in the Bessel code.

Usage

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

Arguments

x

numeric, \ge 0.

nu

numeric; \ge 0 unless in besselK which is symmetric in nu. The order 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

The underlying C code stems from Netlib (http://www.netlib.org/specfun/r[ijky]besl).

If expon.scaled = TRUE, e^{-x} I_{\nu}(x), or e^{x} K_{\nu}(x) are returned.

gammaCody may be somewhat faster but less precise and/or robust than R's standard gamma. It is here for experimental purpose mainly, and may be defunct very soon.

Value

Numeric vector of the same length of x with the (scaled, if expon.scale=TRUE) values of the corresponding Bessel function.

Author(s)

Original Fortran code: W. J. Cody, Argonne National Laboratory
Translation to C and adaption to R: Martin Maechler maechler@stat.math.ethz.ch. nus <- c(0:5,10,20)

x <- seq(0,4, len= 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, leg=paste("nu=",nus), col = nus+2, lwd=1)

x <- seq(0,40,len=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, leg=paste("nu=",nus), col = nus+2, lwd=1)

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+.5))) lines(x0,besselJ(x0,nu=nu), col = nu+2) legend(3,1e50, leg=paste("nu=", paste(nus,nus+.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+.5))) lines(x0,besselK(x0,nu=nu), col = nu+2) legend(3,1e50, leg=paste("nu=", paste(nus,nus+.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, leg=paste("nu=",nus), col = nus+2, lwd=1)

## Check the Scaling : for(nu in nus) print(all(abs(1- besselK(x,nu)*exp( x) / besselK(x,nu,expo=TRUE)) < 2e-15)) for(nu in nus) print(all(abs(1- besselI(x,nu)*exp(-x) / besselI(x,nu,expo=TRUE)) < 1e-15))

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, leg=paste("nu=",nus), col = nus+2, lwd=1)

math

References

Abramowitz, M. and Stegun, I. A. (1972) Handbook of Mathematical Functions. Dover, New York; Chapter 9: Bessel Functions of Integer Order.