When printing out the value of a double, e.g. with %f, the rounding mode isn't taken into account. For instance, the program below: #include <stdio.h> #include <fenv.h> static volatile double x = 1.0 / 3.0; void out (const char *s, int r) { if (fesetround (r)) fprintf (stderr, "fesetround error\n"); else printf ("1/3 rounded %s: %1.6f\n", s, x); } int main (void) { out ("downward", FE_DOWNWARD); out (" upward", FE_UPWARD); return 0; } gives: 1/3 rounded downward: 0.333333 1/3 rounded upward: 0.333333 Since (1.0/3.0) isn't exactly representable in decimal, the rounded decimal values should be different (more precisely, they should be 0.333333 and 0.333334). The ISO C standard says in 7.19.6.1#13: "For e, E, f, F, g, and G conversions, if the number of significant decimal digits is at most DECIMAL_DIG, then the result should be correctly rounded.243) If the number of significant decimal digits is more than DECIMAL_DIG but the source value is exactly representable with DECIMAL_DIG digits, then the result should be an exact representation with trailing zeros. Otherwise, the source value is bounded by two adjacent decimal strings L < U, both having DECIMAL_DIG significant digits; the value of the resultant decimal string D should satisfy L <= D <= U, with the extra stipulation that the error should have a correct sign for the current rounding direction."
Checking the rounding mode from printf would be done using fegetround, which would create a dependency on libm from libc. Comments?
Checking the rounding mode can also be done with two well-chosen additions, for instance. So, no dependency on libm is needed.
You can skirt the libm dependency by using the internal glibc macros _FPU_GETCW and _FPU_SETCW. On architectures which have hardware rounding mode, like x86_64 and powerpc64 it'll access the rounding mode directly in the registers (cw register and FPSCR register respectively). On architectures without a hardware rounding mode it'll use the rounding mode that is already stored in the glibc thread-local-storage.
Created attachment 2706 [details] Quick patch to allow printf_fp to get rounding mode using _FPU_GETCW
I've attached a patch based on Ryan's suggestion to use _FPU_GETCW. It fixes the example given in the description. It may have some problems with the existing ties-to-even code. I've only tested it on POWER. Please check that the patch works on other architectures, and I'll work on the ties-to-even bit.
I haven't looked at the code itself, but the patch looks suspicious to me. Shouldn't L'0' be replaced by L'0' - 1 so that in directed rounding modes, the test "digit > rounddigit" is either always true or always false.
(In reply to comment #6) > Shouldn't L'0' be replaced by L'0' - 1 so that in directed rounding modes, the > test "digit > rounddigit" is either always true or always false. No. You never round when the digit is zero.
(In reply to comment #7) > No. You never round when the digit is zero. Hmm, I was assuming that the exact values were taken apart... There seems to be something missing in the code then. For instance, how does printf("%1.2f\n", 0.5001); round the value upward (0.51) in the FE_UPWARD rounding mode? I think you need some additional flag to say whether the converted value is exact or not, i.e. whether the remaining digits are all zeros.
Created attachment 2713 [details] Revised patch to allow printf_fp to get rounding mode using _FPU_GETCW Vincent, I see what you mean about taking into account whether the rounded digit is followed by zeroes. Rounding upward should never result in a value that is less than the original value, which wasn't caught when the rounded digit is zero followed by nonzero digits. I'm attaching a patch that modifies the existing ties-to-even code to work for FE_UPWARD and FE_DOWNWARD that checks for trailing zeroes when the rounded digit is zero.
the core problem is that you are trying to apply decimal rounding rules to IEEE754 Binary float. For example 0.5001 is really 0x1.000d1b71758e20p-1 (format %.14a which shows the mantisa in hex). Note that there is no exact binary representation for rounding value 0.005 implied by this request. To get what Vincent says he wants would require converting the IEEE754 double to decimal and perform the rounding in decimal (via quantize). It is is not clear that The ISO C standard says in 7.19.6.1#13 requires this: "Otherwise, the source value is bounded by two adjacent decimal strings L < U, both having DECIMAL_DIG significant digits; the value of the resultant decimal string D should satisfy L <= D <= U, with the extra stipulation that the error should have a correct sign for the current rounding direction." This implies that either L or U (.50 or .51) are allowed!
(In reply to comment #10) > For example 0.5001 is really 0x1.000d1b71758e20p-1 In my example, the exact value doesn't matter: the rounding to two digits will be the same for any good approximation to 0.5001 (denoted ~0.5001 below). > "Otherwise, the source value is bounded by two adjacent decimal strings L < U, > both having DECIMAL_DIG significant digits; the value of the resultant decimal > string D should satisfy L <= D <= U, with the extra stipulation that the error > should have a correct sign for the current rounding direction." > > This implies that either L or U (.50 or .51) are allowed! No: if .50 is output, the error does not have a correct sign in rounding upward, since 0.50 < ~0.5001. So, only .51 is allowed here. Now, this part of the ISO C standard does not apply anyway since 2 digits are required and 2 <= DECIMAL_DIG. You need to look at the beginning of the paragraph: "For e, E, f, F, g, and G conversions, if the number of significant decimal digits is at most DECIMAL_DIG, then the result should be correctly rounded." and the correctly-rounded value is .51 here.
Hi all, I have noticed the same problem with libc-2.3.6 (I'm using a precompiled libc6 debian(etch) package). In the example submitted by Vincent, maybe we should calculate x after choosing the rounding mode: since x is evaluated first, the default rounding mode (toward +oo) is used. Altering the rounding mode later won't change the internal binary representation for x (the type qualifier volatile here is without effect). But even doing this don't change the output of printf ! In fact, to get a correct (with respect to IEEE 754) result, I used an intermediate variable, say "b", set to 3.0, we get (see bottom for modified source code): (using %a) 1/3 rounded downward: 0x1.5555555555555p-2 1/3 rounded upward: 0x1.5555555555556p-2 (using %.17f) 1/3 rounded downward: 0.33333333333333331 1/3 rounded upward: 0.33333333333333337 A height precision (at least .17) is needed to see the difference using %f instead of %a (and this is not a normal behaviour, as one should see directly the difference, since we print two different numbers !). To summarize : - fesetround don't influence the internal representation of constants such as 1./3. (or 1./10.) - using %f seems to have some troubles ... Hope this helps. Regards, <<<<<< source modified >>>>>>>> #include <stdio.h> #include <fenv.h> static volatile double x; void out (const char *s, int r) { double b = 3.0; if (fesetround (r)) fprintf (stderr, "fesetround error\n"); else x = 1.0/b; printf ("1/3 rounded %s: %a\n", s, x); } int main (void) { out ("downward", FE_DOWNWARD); out (" upward", FE_UPWARD); return 0; }
The case of decimal output (%e/%E/%f/%F/%g/%G) is fixed for 2.17 by: commit 784761bee3828e4e21cbe15340c7a99e3fd1788b Author: Joseph Myers <joseph@codesourcery.com> Date: Fri Sep 14 20:18:49 2012 +0000 Make printf respect the rounding mode for decimal output (bug 5044). That is, decimal output should now be correctly rounding for all rounding modes, where previously it was for round-to-nearest only. Hex output (%a/%A) has the same issue and still needs to be fixed.
The remaining case, that of hex output (%a, %A), is fixed for 2.17 by: commit a9f8e53a5b14ba481999ded036b025554ea06362 Author: Joseph Myers <joseph@codesourcery.com> Date: Mon Sep 24 15:38:21 2012 +0000 Make printf respect the rounding mode for hex output (bug 5044).
*** Bug 260998 has been marked as a duplicate of this bug. *** Seen from the domain http://volichat.com Page where seen: http://volichat.com/adult-chat-rooms Marked for reference. Resolved as fixed @bugzilla.