5288
Analyzing the Modified Bessel Function of the First Kind with IDL
This Help Artcile discusses "Analyzing the Modified 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:
Analyzing the Modified Bessel Function of the First Kind. In the following example, the recurrence relationship
where I(x) is the modified 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 modified 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 analyzingBESELI
; Derive x values.
x = (DINDGEN(1000) + 1.)/200.
; Initialize display window.
WINDOW, 0, TITLE = 'Modified Bessel Functions'
; Display the first 8 orders of the modified Bessel
; function of the first kind.
PLOT, x, BESELI(x, 0), /XSTYLE, /YSTYLE, $
XTITLE = 'x', YTITLE = 'f(x)', $
TITLE = 'Modified Bessel Functions of the First Kind'
OPLOT, x, BESELI(x, 1), LINESTYLE = 1
OPLOT, x, BESELI(x, 2), LINESTYLE = 2
OPLOT, x, BESELI(x, 3), LINESTYLE = 3
OPLOT, x, BESELI(x, 4), LINESTYLE = 4
OPLOT, x, BESELI(x, 5), LINESTYLE = 5
OPLOT, x, BESELI(x, 6), LINESTYLE = 0
OPLOT, x, BESELI(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*(BESELI(x, (n - 1)) - $
BESELI(x, (n + 1))) - 2.*FLOAT(n)*BESELI(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