Description Usage Arguments Value Author(s) References See Also Examples
Computes the modified Bessel J function, using one of its basic
definitions as an infinite series, e.g. A. & S., p.360, (9.1.10). The
implementation is pure R, working for numeric
,
complex
, but also e.g., for objects of class
"mpfr"
from package Rmpfr.
1 2 |
x |
numeric or complex vector, or of another |
nu |
non-negative numeric (scalar). |
nterm |
integer indicating the number of terms to be used.
Should be in the order of |
log |
logical indicating if the logarithm log J.() is required. |
Ceps |
a relative error tolerance for checking if |
a “numeric” (or complex or "mpfr"
)
vector of the same class and length as x
.
Martin Maechler
Abramowitz, M., and Stegun, I. A. (1955, etc). Handbook of mathematical functions (NBS AMS series 55, U.S. Dept. of Commerce). http://people.math.sfu.ca/~cbm/aands/page_360.htm
This package BesselJ
, base besselJ
, etc
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 | stopifnot(all.equal(besselJs(1:10, 1), # our R code
besselJ (1:10, 1)))# internal C code w/ different algorithm
## Large 'nu' ...
x <- (0:20)/4
(bx <- besselJ(x, nu=200))# base R's -- gives (mostly wrong) warnings
if(require("Rmpfr")) { ## Use high precision, notably large exponent range, numbers:
Bx <- besselJs(mpfr(x, 64), nu=200)
all.equal(Bx, bx, tol = 1e-15)# TRUE -- warnings were mostly wrong; specifically:
cbind(bx, Bx)
signif(asNumeric(1 - (bx/Bx)[19:21]), 4) # only [19] had lost accuracy
## With*out* mpfr numbers -- using log -- is accurate (here)
lbx <- besselJs( x, nu=200, log=TRUE)
lBx <- besselJs(mpfr(x, 64), nu=200, log=TRUE)
cbind(x, lbx, lBx)
stopifnot(all.equal(asNumeric(log(Bx)), lbx, tol=1e-15),
all.equal(lBx, lbx, tol=4e-16))
} # Rmpfr
|
Warning messages:
1: In FUN(X[[i]], ...) : 'nterm=800' may be too small for x=4
2: In FUN(X[[i]], ...) : 'nterm=800' may be too small for x=5
3: In FUN(X[[i]], ...) : 'nterm=800' may be too small for x=6
4: In FUN(X[[i]], ...) : 'nterm=800' may be too small for x=7
[1] 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00
[6] 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00
[11] 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00 0.000000e+00
[16] 0.000000e+00 0.000000e+00 0.000000e+00 3.378127e-305 1.673577e-300
[21] 4.760010e-296
There were 19 warnings (use warnings() to see them)
Loading required package: Rmpfr
Loading required package: gmp
Attaching package: 'gmp'
The following objects are masked from 'package:base':
%*%, apply, crossprod, matrix, tcrossprod
C code of R package 'Rmpfr': GMP using 64 bits per limb
Attaching package: 'Rmpfr'
The following object is masked from 'package:gmp':
outer
The following objects are masked from 'package:stats':
dbinom, dnorm, dpois, pnorm
The following objects are masked from 'package:base':
cbind, pmax, pmin, rbind
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.