Tim Peters <t...@python.org> added the comment:

When Sun was developing fdlibm, I was (among other things) working on a 
proprietary libm for Kendall Square Research.  I corresponded with fdlibm's 
primary author (KC Ng) often at the time.  There's no way he would have left 
errors this egregious slip by, or neglect careful testing near singularities.

But writing portable code at this level is very delicate, and fdlibm failed in 
many ways when I first compiled it for the KSR machine.  That was the first 
64-bit box fdlibm was tried on, & it exposed many unwarranted assumptions in 
the code.

Browsing the OpenBSD port, they bumped into more, seemingly driven by gcc 
violating assumptions about how tricky casting would work - there are many 
places in the code that, e.g., want to view a single double as a pair of int32s 
(or vice versa), at the bit level.  An aggressive compiler can easily be 
provoked into reordering operations in fatally damaging ways, if it doesn't 
realize that, under the maze of cast and union tricks, successive operations 
are actually operating on the _same_ bits.

So, if possible, it would be great to try compiling OpenBSD's libm with every 
conceivable optimization disabled, to see if that changes results.

Why I suspect that may be promising:  as Stefan pointed out, for arguments near 
pi/2 fdlibm's tan() subtracts pi/2 and evaluates the result as tan(x) = 
-1/tan(x - pi/2).  But these arguments are _so_ close to pi/2 that |x - pi/2| 
is "tiny".  But for tiny arguments z, tan(z) == z to machine precision (tan(z) 
= sin(z) / cos(z) ~= z/1 = z for tiny |z|).  So it's effectively using tan(x) = 
-1/(x - pi/2) in this range.

In the original failing, test, behold:

>>> (1.0 / 1978937966095219.0).hex() # correct
'0x1.234c4c6628b81p-51'

>>> (1.0 / 1978945885716843.0).hex() # OpenBSD
'0x1.234c000000000p-51'

It's very hard to believe all those trailing zero bits in the OpenBSD 
reciprocal are a coincidence.  It's extremely hard to believe that fdlibm 
_intended_ them to be all zeroes.  But it's pretty easy to believe an 
aggressive optimizer managed to lose the intent of the code.

----------

_______________________________________
Python tracker <rep...@bugs.python.org>
<https://bugs.python.org/issue31630>
_______________________________________
_______________________________________________
Python-bugs-list mailing list
Unsubscribe: 
https://mail.python.org/mailman/options/python-bugs-list/archive%40mail-archive.com

Reply via email to