Bessel                 package:base                 R Documentation

_B_e_s_s_e_l _F_u_n_c_t_i_o_n_s

_D_e_s_c_r_i_p_t_i_o_n:

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

_U_s_a_g_e:

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

_A_r_g_u_m_e_n_t_s:

       x: numeric, >= 0.

      nu: numeric; The _order_ (maybe fractional!) 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.

_D_e_t_a_i_l_s:

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

_V_a_l_u_e:

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

_A_u_t_h_o_r(_s):

     Original Fortran code: W. J. Cody, Argonne National Laboratory 
      Translation to C and adaption to R: Martin Maechler
     maechler@stat.math.ethz.ch.

_S_o_u_r_c_e:

     The C code is a translation of Fortran routines from <URL:
     http://www.netlib.org/specfun/r[ijky]besl>.

_R_e_f_e_r_e_n_c_e_s:

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

_S_e_e _A_l_s_o:

     Other special mathematical functions, such as 'gamma', Gamma(x),
     and 'beta', B(x).

_E_x_a_m_p_l_e_s:

     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+.5)))
         lines(x0, besselJ(x0, nu=nu), col = nu+2)
     legend(3, 1e50, legend = 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, legend = 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, 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)}))