[788] in java-interest
Re: trig precision problems?
daemon@ATHENA.MIT.EDU (Marianne Mueller)
Thu Jul 20 02:57:04 1995
Date: Wed, 19 Jul 1995 07:57:15 -0700
From: Marianne.Mueller@eng.sun.com (Marianne Mueller)
To: rt@scndprsn.eng.sun.com, misfeldt@mda.ca
Cc: java-interest@java.Eng.Sun.COM
In-Reply-To: <9507182028.AA01565@ganymede.Eng.Sun.COM> (rt@scndprsn.Eng.Sun.COM)
Hi Trevor,
Like rt said ...
misfeldt> double y2 = y - 0.4794;
rt>Your problem is here ^^^^^^
misfeldt> System.out.println( "y2: " + y2 );
The value you get from (y - 0.4794) really has only passing
resemblance to sin(0.5).
Three things to keep in mind:
1) Every time a number is converted from base 10 to base 2, or vice
versa, it suffers at least one ULP (unit in the last place) rounding
error
2) Since SPARC machine instructions are IEEE arithmetic compliant, all
the low-level ops (add, subtract, multiply, divide) suffer at most one
ULP per operation. This is true for x86 boxes as well, and I believe
for modern Macs, but for fun, run the "paranoia" program from netlib
and see how well your box does on floating point
(http://www.netlib.org/paranoia/index.html). Note that sine is not
implemented in hardware on RISC machines, but instead, it's
implemented by a library function that executes probably dozens of
low-level floating point instructions. The Sun libm.a does contain
table-driven trig functions that deliver very accurate results for all
IEEE formats - less than one ULP in most cases, or, as accurate as you
might hope for.
3) Single precision floating point gets you about 6 or 7 decimal
digits, and double precision fp gets you about 15 or 16. But that
means: the first time you enter a base 10 number, then convert it to a
single precision base 2 number, you'll end up with 6 or 7 of the
original digits being faithfully represented at the machine level
(depending on where the original base 10 number falls on the real
axis.) Every time you operate on that number, you could lose an ULP.
Certainly, every time you convert that number back to base 10, you
could lose at least one ULP, maybe a couple. ULPs are measured as
base 2 units, so losing an ULP is not the same as losing a digit. But
analyzing the numerical accuracy of any but the most simple
mathematical progressions quickly gets really hairy.
If you're concerned about getting accuracy beyond 7 decimal digits,
you might be interested in reading "What Every Computer Scientist
Should Know About Floating Point Arithmetic" which was written a few
years ago by David Goldberg, who is now at Xerox PARC I think. The
Sun compiler group distributes this as a postscript file in
/opt/SUNWspro/READMEs/floating-point.ps if you have the unbundled
compilers. You can also get this by anonymous ftp to
sunsite.unc.edu,
ftp://sunsite.unc.edu//pub/sun-info/development-tools/multi-threaded/drop-site/floating-point.ps.Z.
That's a 75-page paper that does contain lots of equations. A more
user-friendly intro to numerical topics for programmers can be found
in the Sun manual "Numerical Computation Guide" which is distributed
as part of the compiler answerbook. The latest part number for that
manual is
801-7639-10 NUMERICAL COMPUTATION GUIDE REV A AUG 94
Marianne
p.s. I knew I shouldn't have gotten started on the fp thread but here
goes. If people are interested in a high quality public domain math
library, the floating point group at Sun did donate such a thing to
the public domain about a year ago. You can get this version of libm
from http://www.netlib.org/fdlibm/index.html
-
Note to Sun employees: this is an EXTERNAL mailing list!
Info: send 'help' to java-interest-request@java.sun.com