Bug 43358 - fma(x,y,z) does not compute x*y+z as a single ternary operation
fma(x,y,z) does not compute x*y+z as a single ternary operation
Status: CLOSED ERRATA
Product: Fedora
Classification: Fedora
Component: glibc (Show other bugs)
rawhide
All Linux
low Severity medium
: ---
: ---
Assigned To: Andreas Schwab
: FutureFeature
Depends On:
Blocks:
  Show dependency treegraph
 
Reported: 2001-06-03 16:11 EDT by Trevin Beattie
Modified: 2016-11-24 10:46 EST (History)
6 users (show)

See Also:
Fixed In Version: glibc-2.12.90-17
Doc Type: Enhancement
Doc Text:
Story Points: ---
Clone Of:
Environment:
Last Closed: 2010-10-19 18:22:42 EDT
Type: ---
Regression: ---
Mount Type: ---
Documentation: ---
CRM:
Verified Versions:
Category: ---
oVirt Team: ---
RHEL 7.3 requirements from Atomic Host:
Cloudforms Team: ---


Attachments (Terms of Use)
Fused multiply and add - long double version for i386, written for MASM (22.25 KB, text/plain)
2001-06-07 20:55 EDT, Trevin Beattie
no flags Details
Test cases for fmal() (i386-specific). If anyone has ideas for additional tests, let me know. The above assembly file passes all these tests. (6.77 KB, text/plain)
2001-06-07 20:58 EDT, Trevin Beattie
no flags Details

  None (edit)
Description Trevin Beattie 2001-06-03 16:11:59 EDT
From Bugzilla Helper:
User-Agent: Mozilla/4.77 [en] (X11; U; Linux 2.4.2-2 i586)

Description of problem:
fma() is currently implemented as two binary operations, which may result
in a loss of precision as x*y is rounded before adding z.


How reproducible:
Always

Steps to Reproduce:
Compile and run the following program:

#include <stdio.h>
#define __USE_ISOC99 1
#include <math.h>

main ()
{
  const double x = 4660.3377777777768642408773303032;
  const double y = 8173.7955555555554383317939937115;
  const double z = -38092648.215387649834156036376953125;
  double a;

  a = fma (x, y, z);
  printf ("%f * %f + %f = %.9g\n", x, y, z, a);
  printf ("True result rounds to %.9g\n",
	  -3.526534563229407686121644747209e-9);
}


Actual Results:  4660.337778 * 8173.795556 + -38092648.215388 =
-3.52520146e-09
True result rounds to -3.52653456e-09


Expected Results:  4660.337778 * 8173.795556 + -38092648.215388 =
-3.52653456e-09
True result rounds to -3.52653456e-09


Additional info:

The fma() function is new in ISO C99.  According to WG14/N869 7.12.13.1:

2	"The fma functions compute the sum z plus the product x times y, rounded
as one ternary operation: they computes the sum z plus the product x times
y (as if) to infinite precision and round once to the result format,
according to the rounding mode characterized by the value of FLT_ROUNDS."

A correct implementation requires the intermediate product x*y to have at
least double the number of significant bits of the operands.  This is
easily achieved in the float version of fmaf by computing the intermediate
result using double precision, but with double and long double, you need
greater precision than most hardware can handle.
Comment 1 Jakub Jelinek 2001-06-05 04:06:48 EDT
On ia64 it works correctly, the rest of architectures need better fma
implementations. I think the current one can well stay in sysdeps/generic,
plus there should be fmaf which does the computation on double (thus is
FP_FAST_FMA) - 53 is more than twice as large as 24, so it will be correct.
For fma and fmal somebody (perhaps me, but not this week) will need to write
a 0.5 ulp precise implementation.
Comment 2 Trevin Beattie 2001-06-07 20:53:26 EDT
I've finished coding and testing an i386-specific implementation that passes all
of the tests I gave it.  It's written in assembly, and uses the Micros**t
AsSeMbler syntax (which is what my employer's custom O/S is written in), but I
think it shouldn't be hard to port to other architectures or generic C code. 
Feel free to study the code and borrow any ideas you like.
Comment 3 Trevin Beattie 2001-06-07 20:55:13 EDT
Created attachment 20593 [details]
Fused multiply and add - long double version for i386, written for MASM
Comment 4 Trevin Beattie 2001-06-07 20:58:57 EDT
Created attachment 20594 [details]
Test cases for fmal() (i386-specific).  If anyone has ideas for additional tests, let me know.  The above assembly file passes all these tests.
Comment 5 Jakub Jelinek 2001-06-08 04:28:02 EDT
We have already a framework for this in glibc (soft-fp), I was just trying
to think of an implementation which would work without the need to decompose
the numbers, do the math in and compose them again.
But you're right, it is probably best to code the generic fma and fmal
using soft-fp (after checking for exceptional values) and arches can then
optimize if they want.
Comment 6 Bill Nottingham 2006-08-07 14:31:17 EDT
Red Hat Linux is no longer supported by Red Hat, Inc. If you are still
running Red Hat Linux, you are strongly advised to upgrade to a
current Fedora Core release or Red Hat Enterprise Linux or comparable.
Some information on which option may be right for you is available at
http://www.redhat.com/rhel/migrate/redhatlinux/.

Red Hat apologizes that these issues have not been resolved yet. We do
want to make sure that no important bugs slip through the cracks.
Please check if this issue is still present in a current Fedora Core
release. If so, please change the product and version to match, and
check the box indicating that the requested information has been
provided. Note that any bug still open against Red Hat Linux on will be
closed as 'CANTFIX' on September 30, 2006. Thanks again for your help.
Comment 7 Trevin Beattie 2006-08-08 12:49:15 EDT
I've verified that the bug is still present in Fedora Core 6 RC2 (glibc-2.4.90-17).
Comment 8 Matthew Miller 2007-04-06 12:11:47 EDT
Fedora Core 5 and Fedora Core 6 are, as we're sure you've noticed, no longer
test releases. We're cleaning up the bug database and making sure important bug
reports filed against these test releases don't get lost. It would be helpful if
you could test this issue with a released version of Fedora or with the latest
development / test release. Thanks for your help and for your patience.

[This is a bulk message for all open FC5/FC6 test release bugs. I'm adding
myself to the CC list for each bug, so I'll see any comments you make after this
and do my best to make sure every issue gets proper attention.]
Comment 9 Trevin Beattie 2007-04-06 22:23:23 EDT
Using the same program, I get different results this time.  It's worse:

[trevin@hydra src]$ gcc -lm -o bug-43358{,.c}
[trevin@hydra src]$ ./bug-43358 
4660.337778 * 8173.795556 + -38092648.215388 = 0
True result rounds to -3.52653456e-09

It may make a difference that my current computer is an AMD64 system, rather
than the Athlon I had before (it died) or the K6-2 I had when I first reported
the bug.  So I compiled and ran the program again on a virtual machine running
the 32-bit version of FC6:

[trevin@fc6-i386 src]$ gcc -lm -o bug-43358{,.c}
[trevin@fc6-i386 src]$ ./bug-43358 
4660.337778 * 8173.795556 + -38092648.215388 = -3.52520146e-09
True result rounds to -3.52653456e-09

Well, what do you know?  No change from the original report.  So it looks like
this bug exists in two forms -- one for the i386, and another for the x86_64.
Comment 10 Matthew Miller 2007-04-06 23:55:00 EDT
This shouldn't be marked as "verified", since that means a *fix* has been
verified. Putting to "assigned".
Comment 11 petrosyan 2008-02-13 19:57:10 EST
I can confirm this bug on Fedora 8, glibc-2.7-2, x86_64 architecture.

$ gcc -lm -o bug-43358{,.c}
$ ./bug-43358
4660.337778 * 8173.795556 + -38092648.215388 = 0
True result rounds to -3.52653456e-09
Comment 12 petrosyan 2008-08-15 01:28:44 EDT
this bug is still present in Fedora 9. glibc-2.8-8

$ gcc -lm -o bug-43358{,.c}
$ ./bug-43358 
4660.337778 * 8173.795556 + -38092648.215388 = 0
True result rounds to -3.52653456e-09
Comment 13 petrosyan 2008-10-01 13:51:11 EDT
this bug is still present in glibc-2.8.90-12.x86_64
Comment 14 Andreas Schwab 2009-06-25 10:40:42 EDT
The compiler now expands fma inline when optimizing.
Comment 15 Fedora Update System 2010-10-14 08:03:47 EDT
glibc-2.12.90-16 has been submitted as an update for Fedora 14.
https://admin.fedoraproject.org/updates/glibc-2.12.90-16
Comment 16 Fedora Update System 2010-10-14 09:51:38 EDT
glibc-2.12.90-16 has been pushed to the Fedora 14 testing repository.  If problems still persist, please make note of it in this bug report.
 If you want to test the update, you can install it with 
 su -c 'yum --enablerepo=updates-testing update glibc'.  You can provide feedback for this update here: https://admin.fedoraproject.org/updates/glibc-2.12.90-16
Comment 17 Kamil Páral 2010-10-19 13:44:37 EDT
Confirmed fix with glibc-2.12.90-16.fc14.x86_64.
Comment 18 Fedora Update System 2010-10-19 18:22:20 EDT
glibc-2.12.90-17 has been pushed to the Fedora 14 stable repository.  If problems still persist, please make note of it in this bug report.

Note You need to log in before you can comment on or make changes to this bug.