This is the mail archive of the
`libc-alpha@sourceware.org`
mailing list for the glibc project.

Index Nav: | [Date Index] [Subject Index] [Author Index] [Thread Index] | |
---|---|---|

Message Nav: | [Date Prev] [Date Next] | [Thread Prev] [Thread Next] |

Other format: | [Raw text] |

*From*: Adhemerval Zanella <azanella at linux dot vnet dot ibm dot com>*To*: "GNU C. Library" <libc-alpha at sourceware dot org>*Date*: Mon, 09 Mar 2015 16:13:37 -0300*Subject*: RFC: powerpc: Incorrect results for pow when using FMA*Authentication-results*: sourceware.org; auth=none

Hi all, Checking on a internal issue report, the following testcase is showing precision errors only for powerpc builds: -- #include <math.h> #include <mpfr.h> int main () { union Di { double d; unsigned long long i; }; union Di x; x.i = 0x553B9BA800000000LL; mpfr_t xm; mpfr_init2 (xm, 512); mpfr_set_d (xm, x.d, MPFR_RNDD); mpfr_t squarem; mpfr_init2 (squarem, 512); mpfr_mul (squarem, xm, xm, MPFR_RNDD); double square = x.d * x.d; mpfr_printf ("mpfr: x*x = %.17Re\n", squarem); mpfr_printf (" : x*x = %.17e\n", square); mpfr_t cubem; mpfr_init2 (cubem, 512); mpfr_mul (cubem, squarem, xm, MPFR_RNDD); double cube = square * x.d; mpfr_printf ("mpfr: x*x*x = %.17Re\n", cubem); mpfr_printf (" : x*x*x = %.17e\n", cube); mpfr_t powm; mpfr_init2 (powm, 512); mpfr_t onehalfm; mpfr_init2 (onehalfm, 512); mpfr_set_d (onehalfm, 1.5, MPFR_RNDD); mpfr_pow (powm, squarem, onehalfm, MPFR_RNDD); double powc = pow (square, 1.5); mpfr_printf ("mpfr: pow (x*x, 1.5) = %.17Re\n", powm); mpfr_printf ("libc: pow (x*x, 1.5) = %.17e\n", powc); return 0; } -- With current master I am seeing: $ ./testrun.sh ./test-mpfr mpfr: x*x = 1.49357829132402765e+205 : x*x = 1.49357829132402765e+205 mpfr: x*x*x = 5.77220822056602683e+307 : x*x*x = 5.77220822056602633e+307 mpfr: pow (x*x, 1.5) = 5.77220822056602683e+307 libc: pow (x*x, 1.5) = 5.77220822056619098e+307 While I would expect the mpfr and libc pow result to have the same precision, as for x86_64. At first I though it was a compiler bug, bug looking into this deeper, this is probably a GLIBC build issue. The difference between the power1() -O1 assembler code and -O2 is that the -O2 code makes use of FMA instructions. This causes a slight difference in values of the second param to the __exp1() function call. The __exp1() return value in teh -O1 case returns a invalid value which forces us to call __slow_pow() leading to the "expected" value. In the -O2 case, the __exp1() value is considered valid and we just return it as is, leading to the slight difference in values. If I look at the comments at the top of the source file (sysdeps/ieee754/dbl-64/e_pow.c) where the power1() function lives, I see this comment: /* Assumption: Machine arithmetic operations are performed in */ /* round to nearest mode of IEEE 754 standard. */ It would seem to me that using FMAs would violate the assumption stated in that comment, since some of the internal FAM ops are not rounded before being used. I do notice that if I add the -ffp-contract=off, then we get the "expected" answer. Is the "correct" fix just to get this file compiled with -ffp-contract=off?

Index Nav: | [Date Index] [Subject Index] [Author Index] [Thread Index] | |
---|---|---|

Message Nav: | [Date Prev] [Date Next] | [Thread Prev] [Thread Next] |