[prev in list] [next in list] [prev in thread] [next in thread]
List: netbsd-port-amd64
Subject: Re: Changing x87 precision to full 63bit as default
From: Joerg Sonnenberger <joerg () britannica ! bec ! de>
Date: 2013-11-13 11:25:20
Message-ID: 20131113112520.GA14677 () britannica ! bec ! de
[Download RAW message or body]
On Wed, Nov 13, 2013 at 07:40:31AM +0900, tsugutomo.enami@jp.sony.com wrote:
> Joerg Sonnenberger <joerg@britannica.bec.de> writes:
>
> > What exactly have you tried?
>
> This code:
>
> #include <stdio.h>
> #include <ieeefp.h>
>
> int
> main(int argc, char *argv[])
> {
> volatile double x, y, r;
>
> if (argc > 1) {
> fp_prec_t oldprec;
> oldprec = fpsetprec(FP_PD);
> printf("prec %d -> %d\n", oldprec, FP_PD);
> }
> x = 3002399751580332.0;
> y = 3002399751580331.0;
> r = x / y;
> printf("%.16f, %.16f, %.16f\n", x, y, r);
>
> return 0;
> }
>
> And result was
> % ./t
> 3002399751580332.0000000000000000, 3002399751580331.0000000000000000, 1.0000000000000004
> % ./t setprec
> prec 3 -> 2
> 3002399751580332.0000000000000000, 3002399751580331.0000000000000000, 1.0000000000000002
>
> Compiling with -ffloat-store and/or -fexcess-precision=standard was
> some result.
Ah, thanks. Took a moment to understand what is happening here. The
machine code produced is actually correct -- the problem is the division
itself. With the choosen values, the exact value of r is slightly
smaller than (1 + DBL_EPSILON / 2). In double precision mode, this gets
rounded down. In extended precision mode, this gets rounded *up* to (1+
DBL_EPSILON / 2) for the *internal* 80bit representation and then stored
to double, which results in another rounding up. So yes, I think for
this specific case wrapping the computations with fpsetprec if
FP_EVAL_METHOD == 2 and restoring the old value afterwards is the best
approach.
Joerg
[prev in list] [next in list] [prev in thread] [next in thread]
Configure |
About |
News |
Add a list |
Sponsored by KoreLogic