Bessel Function of the First Kind  Jν(x)
Graphical Error Analysis



D. Baruth, x87@iging.com



Bessel functions' analytical and numerical characteristics are well documented.[1,2]  When getting into the particulars of the numerical evaluation of different algorithms, however, error analyses are sometimes cumbersome or only partially available.  Using the capabilities of computer graphics, it is possible to visualize errors and estimate their magnitude "empirically".  To demonstrate the power of this method, graphical error analysis of the procedure  _JaX [3]  that evaluates the Bessel Function of the First Kind, Jν(x), is performed.



Accuracy and Relative Error

Depending on the real parameters pair (ν,x), different algorithms are used to achieve an optimal mix of accuracy and function evaluation time.  This is achieved by utilizing the computer's 80-bit FPU registers that can handle binary exponents of up to 2±14 and a mantissa of 64 bits.  Values exceeding these extremes are set to infinity and zero, respectively.

The relative error of the function  f ≡ Jν(x)  is defined as  ε = |Δf/f|,  where  Δf  is the absolute error to be visually determined.  The analysis at the convergence-point  f = J100(95)  is captured below:


Convergence-point at J100(95)


From this graph one finds that  ε ≤ 4x10-16.  The time required to evaluate the function,[4] left and right of x = 95, is approximately 25 µs and 3 µs respectively (1 µs = 0.000001 seconds).


Function Evaluation Map 

The procedure  _JaX  evaluates  Jν(x)  for a wide range of arguments and indices.  The Function Evaluation Map below highlights several results:


|x|
 ν 
t  [µs]
(≤)  ε
0
−∞  −  +∞
0.3  −  0.4
0
1
0  −  1000
1.0  −  5.1
10-19
1
−1000.5  −  −0.5
2.0  −  4.5
10-19  −  10-18
10
0  −  1000
4.3  −  5.3
10-19  −  10-18
10
−0.5
3.0
2x10-16
10
−1000.5  −  −10.5
3.8  −  63.9
10-18
100
0  −  10
1.1
10-19  −  10-18
100
100
3.1
10-17
100
1000
1.0
10-18
100
−10.5  −  −0.5
0.8  −  1.2
10-19
100
−100.5
3.2
10-17
100
−1000.5
8.7
10-18
1000
0  −  100
0.9  −  2.6
10-19  −  10-18
1000
1000
8.2
10-16
1000
−1000.5  −  −0.5
0.8  −  8.3
10-18



  1. Handbook of Mathematical Functions, Abramowitz & Stegun, Dover Publications, fifth printing, Chapters 9,10.
  2. Table of Integrals Series and Products, Gradshteyn & Rizhik, Academic Press, 1965, pp951.
  3. The proprietary procedure  _JaX  can be called, through appropriate interfaces, by FORTRAN, C and other X86/87 compilers.
  4. Calculations were performed on an AMD Athlon XP 1700+, 1.47 GHz, microprocessor.  To assure optimum performance, computer programs were written in Assembly (MASM 5.1).



Copyright Dan Baruth © 2006.  All rights reserved.