Fun with Friedrich Bessel
in [Fortran]
|
From: Richard Maine on 27 Feb 2010 16:06 Richard Maine <nospam(a)see.signature> wrote: > steve <kargls(a)comcast.net> wrote: > > > Well, it was impossible to implement general purpose special > > functions in Fortran until the IEEE 754 intrinsic module was > > standardized. > > "Impossible" seems like a bit of an overstatement, particularly insomuch > as I'm sure counterexamples exist. Not that I'm going to bother to > search for suitable citations of such counterexamples. After saying I wasn't going to :-) See <http://jin.ece.uiuc.edu/routines/routines.html> which happens to be the second hit (after the wikipedia definition) when I Google "special function." I do not vouch for the accuracy, efficiency, or anything else of those particular implementations; I didn't bother to research that. If you want to start quibbling about such things, that would be a different topic from the above overly broad claim, which is unqualified with any such conditions. They are certainly sufficient as an existance proof by example. -- Richard Maine | Good judgment comes from experience; email: last name at domain . net | experience comes from bad judgment. domain: summertriangle | -- Mark Twain
From: Tobias Burnus on 27 Feb 2010 16:37 glen herrmannsfeldt wrote: > There have been some bessel function routines in many C libraries > over the years. As far as I know, not part of the standard. I think Bessel functions are part of C since the C99 standard - I am too lazy to check, but "man 3p jn" shows that they part part of POSIX standard. In case of Fortran, the Bessel functions (BESSEL_J0 etc.) are part of the Fortran 2008, which has almost reached the FDIS stage (final draft international standard) - the FDIS ballot is scheduled for May. > It is interesting that they have been associated with Unix/C > longer than with Fortran... My feeling is that the new functions only ended up in the Fortran standard after being standardized in the C and thus part of the systems C library (e.g. glibc). On the other hand, despite that C99 is now about 10 years old, there are many systems which do not have these functions. GNU Fortran provides a c99 library, which offers replacement functions for those which are not part of the system - thus you should be able to use all of Fortran 90/95/2003/2008 (as far as implemented in gfortran) even those which you cannot use in C due to the incomplete libc/libm. For instance: Complex argument to inverse trigonometric functions, which were added in Fortran 2008 and are supported in gfortran 4.5 [also as initialization expressions]. Admittedly, the used algorithm in gfortran's run-time library is not perfectly tuned - especially for IEEE numbers, but usually you only get the expected ULP error. Frankly, I have not checked what gfortran's C99 library does for jn as it uses on my system the (buggy) system's GLIBC function. * * * Talking about math functions. Does anyone know what happened to the proposal to support the extra math functions to Fortran, which were specified in C's TR 24747 [1,2]? I remember seeing something about it on the J3 mailing list, though I have not heard anything lately. Does anyone know a (C) library / C compiler which implements those already? (The TR provides Laguerre/Legendre/Hermite polynominals, beta/zeta function, elliptical integrals, cylindrical/spherical Bessel/Neumann functions.) [1] http://www.iso.org/iso/catalogue_detail.htm?csnumber=38857 [2] http://www.open-std.org/jtc1/sc22/wg21/docs/papers/2009/n2850.pdf Tobias
From: Richard Maine on 27 Feb 2010 16:49 Tobias Burnus <burnus(a)net-b.de> wrote: > My feeling is that the new functions only ended up in the Fortran > standard after being standardized in the C and thus part of the systems > C library (e.g. glibc). That might have been part of the rationale; I couldn't say for sure. But I don't think it holds water if one is viewing it as "this is minimal work because the system will already have it." In particular, I haven't checked the details, but I'd expect to find that glibc didn't have implementations for all of the precisions required by the draft F2008 standard. Does it have both single and double precision versions? How about quad? I don't recall that f2008 requires quad support, but I bet that if a compilere does support quad, it would be required to have quad versions of those functions as well; anyway, that's the way most of the intrinsics are. -- Richard Maine | Good judgment comes from experience; email: last name at domain . net | experience comes from bad judgment. domain: summertriangle | -- Mark Twain
From: glen herrmannsfeldt on 27 Feb 2010 17:04 Richard Maine <nospam(a)see.signature> wrote: (snip) > That might have been part of the rationale; I couldn't say for sure. But > I don't think it holds water if one is viewing it as "this is minimal > work because the system will already have it." In particular, I haven't > checked the details, but I'd expect to find that glibc didn't have > implementations for all of the precisions required by the draft F2008 > standard. Does it have both single and double precision versions? How > about quad? I don't recall that f2008 requires quad support, but I bet > that if a compilere does support quad, it would be required to have quad > versions of those functions as well; anyway, that's the way most of the > intrinsics are. C99 now has single precision versions of all the math functions, unlike previous C versions. There is also a (long double) type which may or may not be longer than double. It likely is if the hardware implements something longer. C99 also has complex types and, I believe, the appropriate math functions. -- glen
From: glen herrmannsfeldt on 27 Feb 2010 17:08
Richard Maine <nospam(a)see.signature> wrote: (snip) > After saying I wasn't going to :-) See > <http://jin.ece.uiuc.edu/routines/routines.html> > which happens to be the second hit (after the wikipedia definition) when > I Google "special function." I do not vouch for the accuracy, > efficiency, or anything else of those particular implementations; I > didn't bother to research that. If you want to start quibbling about > such things, that would be a different topic from the above overly broad > claim, which is unqualified with any such conditions. They are certainly > sufficient as an existance proof by example. I believe those are the ones from the book I mentioned previously. The authors are Zhang and Jin, and from UIUC. The book has Fortran 77 versions. The book comes with a floppy disk with the usual "install on one computer" copyright restriction. The question of copyright restrictions on code in books is always complicated. Since copyright applies to the expression of the idea, not the idea itself, there is always the suggestion that if you write a new version based on the code, then the previous copyright doesn't apply. (That is, modify it sufficiently for some definition of sufficiently.) -- glen |