This is the mail archive of the
`gcc@gcc.gnu.org`
mailing list for the GCC 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*: Alan Modra <amodra at bigpond dot net dot au>*To*: gcc at gcc dot gnu dot org*Cc*: David Edelsohn <dje at watson dot ibm dot com>*Date*: Wed, 18 Sep 2002 00:36:03 +0930*Subject*: libgcc2.c floatdisf is broken

floatdisf, 64 bit int to float (32 bit IEEE, ie. 23 bit fraction) If implemented as DI -> DF DF -> SF then there is a possible double rounding problem, when the first conversion is not exact, and the rounding mode is "round to nearest". Rounding modes of "toward zero", "toward +inf" or "toward -inf" do not suffer this problem. The following was written with rs6000.md floatdisf2 pattern in mind, but the same analysis applies to libgcc2.c __floatdisf function where the DI -> DF rounding occurs when the low part of the DImode value is added to the high part. Examples: 2**60 is representable as a SF number, the next higher SF value being 2**60 + 2**37. Considering integers in this range, according to the "round to nearest" rules, values in the range [2**60, 2**60 + 2**36] should be converted to a result of 2**60, while integers in the range [2**60 + 2**36 + 1, 2**60 + 2**37] should convert to 2**60 + 2**37. Note that since powerpc hardware uses "round to nearest even", the midpoint, 2**60 + 2**36, is converted to the lower value in this case as the lower value has a zero LSB. If the lower value had the LSB set, we'd convert the midpoint up. Let's see what happens with 2**60 + 2**36 + 1. The nearest DF value (52 bit fraction) is 2**60 + 2**36 (the next highest value is 2**60 + 2**36 + 2**8), so the result of the first conversion is a number that is further rounded down in the DF->SF conversion, and we end up with an incorrect result of 2**60. It's fairly easy to see that integers in the range [2**60 + 2**36 + 1, 2**60 + 2**36 + 2**7] will similarly produce the wrong result. The next higher SF range is [2**60 + 2**37, 2**60 + 2**38]. Here, the midpoint value of 2**60 + 3*2**36 converts up to 2**60 + 2**38. What happens with 2**60 + 3*2**36 - 1? Well, the DI->DF conversion results in 2**60 + 3*2**36, which goes to 2**60 + 2**38 on the DF->SF conversion. Wrong again. This will happen for integers in the range [2**60 + 3*2**36 - 2**7, 2**60 + 3*2**36 - 1]. Details: As the examples indicate, the problem occurs when the first conversion results in a number that falls on the midpoint between two SF values, and the midpoint is rounded away from the original value. ie. The original integer is a positive number of the form (A) 0..01XX..X010000000000000000000000000000Y..Y where XX..X is a string of 22 binary digits, and Y..Y is a string of 1 or more digits with value 0..01 to 10..0 or (B) 0..01XX..X101111111111111111111111111111Y..Y where XX..X is a string of 22 binary digits, and Y..Y is a string of 1 or more digits with value 10..0 to 11..1 or a negative number with a bit pattern matching the complement of the above patterns. Values matching (A) will be incorrectly rounded down, while values matching (B) will be incorrectly rounded up. Conclusion: The attempted fix in libgcc2.c:__floatdisf for the double rounding problem won't work. Case (B) cannot be fixed by setting a bit. I think the following patch will cure case (B) as well as (A), but haven't properly tested for compliance with IEEE rounding requirements yet. ie. I'm not at this point requesting permission to apply the patch. It's more a case of inviting comment. Also, can someone suggest a test program? * libgcc2.c (__floatdisf): Mask off bits below REP_BIT. * config/rs6000/rs6000.md (floatdisf2): Rename to floatdisf2_int1. (floatdisf2): New define_expand. (floatdisf2_int2): Likewise. Index: gcc/libgcc2.c =================================================================== RCS file: /cvs/gcc/gcc/gcc/libgcc2.c,v retrieving revision 1.149 diff -u -p -r1.149 libgcc2.c --- gcc/libgcc2.c 8 Sep 2002 12:47:26 -0000 1.149 +++ gcc/libgcc2.c 17 Sep 2002 13:42:09 -0000 @@ -1091,7 +1091,10 @@ __floatdisf (DWtype u) && u < ((DWtype) 1 << DF_SIZE))) { if ((UDWtype) u & (REP_BIT - 1)) - u |= REP_BIT; + { + u &= ~ (REP_BIT - 1); + u |= REP_BIT; + } } } f = (Wtype) (u >> WORD_SIZE); Index: gcc/config/rs6000/rs6000.md =================================================================== RCS file: /cvs/gcc/gcc/gcc/config/rs6000/rs6000.md,v retrieving revision 1.206 diff -u -p -r1.206 rs6000.md --- gcc/config/rs6000/rs6000.md 13 Sep 2002 01:40:44 -0000 1.206 +++ gcc/config/rs6000/rs6000.md 17 Sep 2002 14:55:58 -0000 @@ -5890,13 +5890,30 @@ "fctidz %0,%1" [(set_attr "type" "fp")]) -;; This only is safe if rounding mode set appropriately. -(define_insn_and_split "floatdisf2" +(define_expand "floatdisf2" + [(set (match_operand:SF 0 "gpc_reg_operand" "") + (float:SF (match_operand:DI 1 "gpc_reg_operand" "")))] + "TARGET_POWERPC64 && TARGET_HARD_FLOAT && TARGET_FPRS" + " +{ + if (!flag_unsafe_math_optimizations) + { + rtx label = gen_label_rtx (); + emit_insn (gen_floatdisf2_int2 (operands[1], label)); + emit_label (label); + } + emit_insn (gen_floatdisf2_int1 (operands[0], operands[1])); + DONE; +}") + +;; This is not IEEE compliant if rounding mode is "round to nearest". +;; If the DI->DF conversion is inexact, then it's possible to suffer +;; from double rounding. +(define_insn_and_split "floatdisf2_int1" [(set (match_operand:SF 0 "gpc_reg_operand" "=f") (float:SF (match_operand:DI 1 "gpc_reg_operand" "*f"))) (clobber (match_scratch:DF 2 "=f"))] - "TARGET_POWERPC64 && TARGET_HARD_FLOAT && TARGET_FPRS - && flag_unsafe_math_optimizations" + "TARGET_POWERPC64 && TARGET_HARD_FLOAT && TARGET_FPRS" "#" "&& reload_completed" [(set (match_dup 2) @@ -5904,6 +5921,31 @@ (set (match_dup 0) (float_truncate:SF (match_dup 2)))] "") + +;; Twiddles bits to avoid double rounding. +;; See libgcc2.c __floatdisf, except that this is somewhat optimized. +(define_expand "floatdisf2_int2" + [(set (match_dup 2) (and:DI (match_operand:DI 0 "" "") (const_int 2047))) + (set (match_dup 4) (compare:CC (match_dup 2) (const_int 0))) + (set (match_dup 3) (ashiftrt:DI (match_dup 0) (const_int 53))) + (set (match_dup 3) (plus:DI (match_dup 3) (const_int 1))) + (set (pc) (if_then_else (eq (match_dup 4) (const_int 0)) + (label_ref (match_operand:DI 1 "" "")) + (pc))) + (set (match_dup 5) (compare:CCUNS (match_dup 3) (const_int 2))) + (set (pc) (if_then_else (ltu (match_dup 5) (const_int 0)) + (label_ref (match_dup 1)) + (pc))) + (set (match_dup 0) (xor:DI (match_dup 0) (match_dup 2))) + (set (match_dup 0) (ior:DI (match_dup 0) (const_int 2048)))] + "TARGET_POWERPC64 && TARGET_HARD_FLOAT && TARGET_FPRS" + " +{ + operands[2] = gen_reg_rtx (DImode); + operands[3] = gen_reg_rtx (DImode); + operands[4] = gen_reg_rtx (CCmode); + operands[5] = gen_reg_rtx (CCUNSmode); +}") ;; Define the DImode operations that can be done in a small number ;; of instructions. The & constraints are to prevent the register -- Alan Modra IBM OzLabs - Linux Technology Centre

**Follow-Ups**:**Re: libgcc2.c floatdisf is broken***From:*Richard Henderson

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

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