[Top][All Lists]

[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]
## [Help-gsl] Hankel transform and Bessel functions zeros

**From**: |
Michael Carley |

**Subject**: |
[Help-gsl] Hankel transform and Bessel functions zeros |

**Date**: |
Fri, 26 Aug 2011 11:28:14 +0100 (BST) |

**User-agent**: |
Alpine 2.00 (LNX 1167 2008-08-23) |

`I had a query about the Discrete Hankel Transform a while ago which turned
``out to be an error on my side rather than in GSL but investigating it has
``turned up an issue.
`

`There is a loss of accuracy in computing the zeros of Bessel functions in
``gsl_sf_bessel_zero_Jnu. Running the code at the end of this email to
``estimate the error in the eighth order Bessel function at the computed
``zeros gives:
`
i log10(fabs(Jn))
1 -14.69
2 -14.75
3 -14.81
4 -15.42
5 -8.16
6 -14.67
7 -14.89
8 -15.04
9 -15.61
10 -15.54
11 -14.20
12 -14.36
13 -14.60
14 -14.67
15 -14.70

`and you can see that the Bessel function at the fifth zero is about 10^-8,
``rather than 10^-15. There are similar results for different orders of
``Bessel function. Looking at the code for computing the zeros, it seems
``that this is the point where the method for estimating the zeros changes.
``There is a similar (though not so large) variation in the error at the
``higher roots.
`

`If I use the method of `A Fast Algorithm for the Calculation of the Roots
``of Special Functions':
`
http://dx.doi.org/10.1137/06067016X

`I get much better results with a roughly constant error and if I use this
``method to compute the Bessel function zeros for the Hankel transform, that
``is also improved. I had already implemented the root finding function in
``my Gaussian Quadrature library:
`
http://www.paraffinalia.co.uk/Software/gqr.shtml

`It seems to me that it might be useful to include an implementation of the
``root finding algorithm in GSL since it could be used for a wide range of
``special functions, not just Bessel functions, and would also improve the
``accuracy of the Hankel transform, and I would be happy to recode my
``implementation to GSL standards or pass it on to someone else who could do
``so.
`
/*test code for zeros of Bessel functions*/
#include <math.h>
#include <stdio.h>
#include <gsl/gsl_sf_bessel.h>
int main(int argc, char **argv)
{
double x0, Jn ;
int n = 8, i ;
for ( i = 1 ; i < 16 ; i ++ ) {
x0 = gsl_sf_bessel_zero_Jnu(n, i) ;
Jn = gsl_sf_bessel_Jn(n, x0) ;
fprintf(stdout, "%d %1.2f\n", i, log10(fabs(Jn))) ;
}
return 0 ;
}
gives
--
`To tell the truth, let us be honest at least, it is some considerable
time since I last knew what I was talking about.'
No MS attachments: http://www.gnu.org/philosophy/no-word-attachments.html
Home page: http://people.bath.ac.uk/ensmjc/

[Prev in Thread] |
**Current Thread** |
[Next in Thread] |

**[Help-gsl] Hankel transform and Bessel functions zeros**,
*Michael Carley* **<=**