Incorrect rounding in floating-point operations with gcc/x87

[mdickinson]

BPO	2937
Nosy	@loewis, @mdickinson

^{Note: these values reflect the state of the issue at the time it was migrated and might not reflect the current state.}

Show more details

GitHub fields:

assignee = None
closed_at = <Date 2008-12-05.21:32:58.705>
created_at = <Date 2008-05-21.20:36:56.457>
labels = ['interpreter-core', 'type-bug', 'invalid']
title = 'Incorrect rounding in floating-point operations with gcc/x87'
updated_at = <Date 2008-12-05.21:32:58.704>
user = 'https://github.com/mdickinson'

bugs.python.org fields:

activity = <Date 2008-12-05.21:32:58.704>
actor = 'mark.dickinson'
assignee = 'none'
closed = True
closed_date = <Date 2008-12-05.21:32:58.705>
closer = 'mark.dickinson'
components = ['Interpreter Core']
creation = <Date 2008-05-21.20:36:56.457>
creator = 'mark.dickinson'
dependencies = []
files = []
hgrepos = []
issue_num = 2937
keywords = []
message_count = 4.0
messages = ['67174', '67180', '67194', '77073']
nosy_count = 2.0
nosy_names = ['loewis', 'mark.dickinson']
pr_nums = []
priority = 'normal'
resolution = 'not a bug'
stage = None
status = 'closed'
superseder = None
type = 'behavior'
url = 'https://bugs.python.org/issue2937'
versions = ['Python 2.6', 'Python 3.0']

[mdickinson]

On some older Intel 32-bit hardware, under Linux, floating-point
operations don't always give correctly rounded results. Here's an
example involving addition, on SuSE Linux 10.2/Xeon.

Python 2.6a3+ (trunk:63521, May 21 2008, 15:40:39) 
[GCC 4.1.2 20061115 (prerelease) (SUSE Linux)] on linux2
Type "help", "copyright", "credits" or "license" for more information.
>>> 1e16 + 2.999
10000000000000002.0
>>> 1e16 + 2.9999
10000000000000004.0

The second result should really be 1e16+2., not 1e16+4. This appears to
be related to this GCC issue:

http://gcc.gnu.org/bugzilla/show_bug.cgi?id=323

Various fixes are possible. One possible fix is to add the -ffloat-
store flag to the gcc options. Another is to use the information in
fpu_control.h, if available, to set the precision control. Yet another
is to sprinkle some 'volatile' modifiers throughout floatobject.c.

It's not clear to me that this *should* be fixed, but I think the
problem should at least be documented. Hence this bug report.

[loewis]

Python (the language) makes no guarantee itself on what the precise
semantics of floating-point operations is. This is documented in

http://docs.python.org/ref/types.html

"These represent machine-level double precision floating point numbers.
You are at the mercy of the underlying machine architecture (and C or
Java implementation) for the accepted range and handling of overflow."

If you want to, one could add ", precision" in the sentence; I think it
is fine as it stands.

[mdickinson]

Okay; so this is definitely not a Python bug---it's a well-known
and well-documented problem with IA32 floating-point. And I accept
that it's really not Python's responsibility to document this, either.

Nevertheless, it was a surprise to me when my (supposedly IEEE 754
compliant) Pentium 4 box produced this. I probably shouldn't have
been surprised. I'm aware of issues with 80-bit extended precision when
programming in C, but naively expected that Python would be largely
immune from these, since it's always going to force intermediate results
from (80-bit) floating-point registers into (64-bit) memory slots.

There's an excellent recent article by David Monniaux, "The pitfalls of
verifying floating-point computations.", that's available online at

http://hal.archives-ouvertes.fr/hal-00128124

that explains exactly what's going on here (it's a case of double-
rounding, as described in section 3.1.2 of that paper).

Do you think a documentation patch that added this reference, along with
the oft-quoted "What Every Computer Scientist Should Know About
Floating-Point Arithmetic" by David Goldberg, to Appendix B of the
tutorial would be acceptable?

One other thing that's worth mentioning: on Pentium 4 and later, the
gcc flags "-mfpmath=sse -msse2" appear to fix the problem, by forcing
gcc to use the SSE floating-point unit instead of the x87-derived one.

In any case, I guess this report should be closed as 'invalid', but I
hope that at least others who encounter this problem manage to find this
bug report.

[mdickinson]

Closing as invalid.