Re: gnuplot numerical errors on Intel MacBook

[Joe K. <jko...@co...>]

on 9/24/06 5:23 PM, Hans-Bernhard Br=F6ker at br...@ph...
wrote:

> CM wrote:
>=20
>> configure:8573: checking whether signgam is declared
>> configure:8607: gcc -c -g -O2   -I/usr/X11R6/include conftest.c >&5
>> configure:8613: $? =3D 0
>> configure:8620: test -z "$ac_c_werror_flag" || test ! -s conftest.err
>> configure:8623: $? =3D 0
>> configure:8630: test -s conftest.o
>> configure:8633: $? =3D 0
>> configure:8645: result: yes
>>=20
>> Is this what you expected?
>=20
> Not really.  What this set of observations means is that the signgam
> aspect of the gamma() implementation on Intel Macs is rather
> interestingly broken.  That's about as far as gnuplot's end of the
> problem can be traced.  For the rest you'll have to bring this to the
> attention of Mac experts, or maybe Apple themselves.  It might be some
> unwanted side effect of someone's try to make gamma() thread-safe.
>=20
My total knowledge of this subject consists of owning a copy of Abramowitz
and Stegun, so please excuse me it the following is a red herring.

On his daily blog on the g95 website, g95.sourceforge.net, Andy Vaught wrot=
e
=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D
September 6
Bill McKie pointed out that x86/OSX wasn't trapping overflows during real
exponentiation. It turns out that OSX uses the MMX unit to do exponentiatio=
n
instead of the x87. So I added code to set the bits for MMX traps. Not only
does the MMX exponentiation use about six times the code that the x87 uses,
but it turns out that OSX doesn't handle traps from the MMX unit correctly
either, reporting only a generic floating point exception. To add insult to
injury, it takes about 3/4 of second (!) to process a floating point
exception on OSX, be it x86 or powerpc.

September 13
The origin of this problem was detecting overflow during exponentiation on
OSX/386, which uses the SSE unit to do the math. I ran a speed test
comparing exponentiation using the SSE unit with x87, and the x87 is 25%
faster, and delivers a more accurate result because the intermediate
quantities are carried to 64 bits in the mantissa. So the solution is not t=
o
bend over backwards for the SSE unit, but rather use the x87 for
exponentiation where it is available.
=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D=3D
So the Intel Mac does math somewhat differently than its Windows cousins.
This may or may not be relevant here, but I just thought I'd point it out.

Joe