[Date Prev][Date Next][Thread Prev][Thread Next][Date Index][Thread Index]

Re: a plea for more reliable mathematical routines



George White wrote:
> Add-ons are one thing, but first the basic math functions should be
> implemented properly.  It is interesting to read Cleve Moler's comments
> on http://math.nist.gov/javanumerics/
> 
> See Appendix 3 of the June 7th, 1999 "Recent Progress of the Java Grande
> Numerics Working Group".  Moler notes that although the Java specification
> calls for numeric functions that agree with Sun's "freely distributable
> math library" (FDLIBM), in fact both Sun and Microsoft use the Microsoft
> Visual C++ numerical library in their Java implementations for Win32.  A
> quote:
> 
>   "The exponential, sine, and cosine functions from Microsoft's math
>   library are so innaccurate for large arguments that the library is
>   unsuitable for general-purpose use."
> 
> Matlab uses FDLIBM to help insure the same results on all platforms.
> What does IDL use?

Moler's comments include the following:

---begin quote---
Let pi be the transcendental mathematical quantity usually denoted by a
Greek letter and let PI be Java.lang.Math.PI, the floating point number
closest to pi. What is the correct value for Java.lang.Math.sin(PI)? The
answer is not zero, because PI is not equal to pi. Here are two answers,
one obtained from the FSIN instruction on an Intel Pentium and the other
obtained from FDLIBM:

     1.224606353822377e-16
     1.224646799147353e-16

These two values agree to only five significant figures. If the FDLIBM
result is regarded as the correct answer, then the error in the FSIN
result is 1.64e+11 ulps, or 164 gigaulps. It can be argued that both
values are so small that either one is an acceptable result, but the
mere fact that more than one answer is possible violates both the
machine independence objective and the idealized model.
---begin quote---

I assume the first result is from a Pentium FSIN instruction, and the
second result is from FDLIBM. Here's what I get from IDL:

Windows NT4 and IDL 5.2:
IDL> print, sin(!dpi), format='(e22.15)'
1.224606353822377e-016

SGI Irix 6.5 and IDL 5.2:
IDL> print, sin(!dpi), format='(e22.15)'
 1.224646799147353e-16

I'd guess IDL is using the platform vendor supplied libm.

Cheers,
Liam.

-- 
Liam E. Gumley
Space Science and Engineering Center, UW-Madison
http://cimss.ssec.wisc.edu/~gumley