Analyzing the Bessel Function of the First Kind with IDL

This Help Article discusses "Analyzing the Bessel Function of the First Kind."

Determining Bessel Function Accuracy. A recurrence relationship between Bessel functions of differing order can be used to determine how accurately IDL is computing Bessel functions.

For Related information, see the following Help Artciles:
• “#3184: Analyzing the Bessel Function of the Second Kind”
• “#3185: Analyzing the Modified Bessel Function of the First Kind”
• “#3186: Analyzing the Modified Bessel Function of the Second Kind”

Analyzing the Bessel Function of the First Kind. In the following example, the recurrence relationship:



where J(x) is the Bessel function of the first kind of order n - 1, n, or n + 1 is used. Subtracting the Bessel function of order n from both sides results in the following equation.



Evaluating the left side of this equation will reveal how accurately IDL computes the Bessel function of the first kind.

The resulting plots are for n equal to 1 through 6. All of these plots show that this Bessel function is calculated within machine tolerance.

Code example: 

PRO analyzingBESELJ

    ; Derive x values.
    x = (DINDGEN(1000) + 1.)/100.
 
    ; Initialize display window.
    WINDOW, 0, TITLE = 'Bessel Functions'
 
    ; Display the first 8 orders of the Bessel function of
    ; the first kind.
    PLOT, x, BESELJ(x, 0), /XSTYLE, /YSTYLE, $
       XTITLE = 'x', YTITLE = 'f(x)', $
       TITLE = 'Bessel Functions of the First Kind'
    OPLOT, x, BESELJ(x, 1), LINESTYLE = 1
    OPLOT, x, BESELJ(x, 2), LINESTYLE = 2
    OPLOT, x, BESELJ(x, 3), LINESTYLE = 3
    OPLOT, x, BESELJ(x, 4), LINESTYLE = 4
    OPLOT, x, BESELJ(x, 5), LINESTYLE = 5
    OPLOT, x, BESELJ(x, 6), LINESTYLE = 0
    OPLOT, x, BESELJ(x, 7), LINESTYLE = 1
 
    ; Initialize display window for recurrence relations.
    WINDOW, 1, XSIZE = 896, YSIZE = 512, $
       TITLE = 'Testing the Recurrence Relations'
    !P.MULTI = [0, 2, 3, 0, 0]
 
    ; Initialize title variable.
    nString = ['0', '1', '2', '3', '4', '5', '6', '7']
 
    ; Display recurrence relationships for order 1 to 6.
    ; NOTE:  the results of these relationships should be
    ; very close to zero.
    FOR n = 1, 6 DO BEGIN
       equation = x*(BESELJ(x, (n - 1)) + $
          BESELJ(x, (n + 1))) - 2.*FLOAT(n)*BESELJ(x, n)
       PLOT, x, equation, /XSTYLE, /YSTYLE, CHARSIZE = 1.5, $
          TITLE = 'n = ' + nString[n] + ':  Orders of ' + $
          nString[n - 1] + ', ' + nString[n] + ', and ' + $
          nString[n + 1]
       PRINT, 'n = ' + nString[n] + ':  '
       PRINT, 'minimum = ', MIN(equation)
       PRINT, 'maximum = ', MAX(equation)
    ENDFOR
 
    ; Return display window back to its default setting, one
    ; display per window.
    !P.MULTI = 0

 
END

_____________________________________________________
Reviewed by BC on 09/05/2014