soft-fp64/fadd: Rename aFrac and bFrac variables
[mesa.git] / src / compiler / glsl / float64.glsl
index 7f6d3a86e340b0b5de8711a6b8b5e8f798e97964..cdecfdb8dbbdacfb595a66b5b676db1832d8b2ff 100644 (file)
 #define FLOAT_ROUND_UP              3
 #define FLOAT_ROUNDING_MODE         FLOAT_ROUND_NEAREST_EVEN
 
+/* Relax propagation of NaN.  Binary operations with a NaN source will still
+ * produce a NaN result, but it won't follow strict IEEE rules.
+ */
+#define RELAXED_NAN_PROPAGATION
+
 /* Absolute value of a Float64 :
  * Clear the sign bit
  */
@@ -165,14 +170,14 @@ __extractFloat64Sign(uint64_t a)
    return unpackUint2x32(a).y & 0x80000000u;
 }
 
-/* Returns true if the 64-bit value formed by concatenating `a0' and `a1' is less
- * than the 64-bit value formed by concatenating `b0' and `b1'.  Otherwise,
- * returns false.
+/* Returns true if the signed 64-bit value formed by concatenating `a0' and
+ * `a1' is less than the signed 64-bit value formed by concatenating `b0' and
+ * `b1'.  Otherwise, returns false.
  */
 bool
-lt64(uint a0, uint a1, uint b0, uint b1)
+ilt64(uint a0, uint a1, uint b0, uint b1)
 {
-   return (a0 < b0) || ((a0 == b0) && (a1 < b1));
+   return (int(a0) < int(b0)) || ((a0 == b0) && (a1 < b1));
 }
 
 bool
@@ -180,12 +185,42 @@ __flt64_nonnan(uint64_t __a, uint64_t __b)
 {
    uvec2 a = unpackUint2x32(__a);
    uvec2 b = unpackUint2x32(__b);
-   uint aSign = __extractFloat64Sign(__a);
-   uint bSign = __extractFloat64Sign(__b);
-   if (aSign != bSign)
-      return (aSign != 0u) && ((((a.y | b.y)<<1) | a.x | b.x) != 0u);
 
-   return mix(lt64(a.y, a.x, b.y, b.x), lt64(b.y, b.x, a.y, a.x), aSign != 0u);
+   /* IEEE 754 floating point numbers are specifically designed so that, with
+    * two exceptions, values can be compared by bit-casting to signed integers
+    * with the same number of bits.
+    *
+    * From https://en.wikipedia.org/wiki/IEEE_754-1985#Comparing_floating-point_numbers:
+    *
+    *    When comparing as 2's-complement integers: If the sign bits differ,
+    *    the negative number precedes the positive number, so 2's complement
+    *    gives the correct result (except that negative zero and positive zero
+    *    should be considered equal). If both values are positive, the 2's
+    *    complement comparison again gives the correct result. Otherwise (two
+    *    negative numbers), the correct FP ordering is the opposite of the 2's
+    *    complement ordering.
+    *
+    * The logic implied by the above quotation is:
+    *
+    *    !both_are_zero(a, b) && (both_negative(a, b) ? a > b : a < b)
+    *
+    * This is equivalent to
+    *
+    *    fne(a, b) && (both_negative(a, b) ? a >= b : a < b)
+    *
+    *    fne(a, b) && (both_negative(a, b) ? !(a < b) : a < b)
+    *
+    *    fne(a, b) && ((both_negative(a, b) && !(a < b)) ||
+    *                  (!both_negative(a, b) && (a < b)))
+    *
+    * (A!|B)&(A|!B) is (A xor B) which is implemented here using !=.
+    *
+    *    fne(a, b) && (both_negative(a, b) != (a < b))
+    */
+   bool lt = ilt64(a.y, a.x, b.y, b.x);
+   bool both_negative = (a.y & b.y & 0x80000000u) != 0;
+
+   return !__feq64_nonnan(__a, __b) && (lt != both_negative);
 }
 
 /* Returns true if the double-precision floating-point value `a' is less than
@@ -195,10 +230,15 @@ __flt64_nonnan(uint64_t __a, uint64_t __b)
 bool
 __flt64(uint64_t a, uint64_t b)
 {
-   if (__is_nan(a) || __is_nan(b))
-      return false;
+   /* This weird layout matters.  Doing the "obvious" thing results in extra
+    * flow control being inserted to implement the short-circuit evaluation
+    * rules.  Flow control is bad!
+    */
+   bool x = !__is_nan(a);
+   bool y = !__is_nan(b);
+   bool z = __flt64_nonnan(a, b);
 
-   return __flt64_nonnan(a, b);
+   return (x && y && z);
 }
 
 /* Returns true if the double-precision floating-point value `a' is greater
@@ -209,19 +249,45 @@ __flt64(uint64_t a, uint64_t b)
 bool
 __fge64(uint64_t a, uint64_t b)
 {
-   if (__is_nan(a) || __is_nan(b))
-      return false;
+   /* This weird layout matters.  Doing the "obvious" thing results in extra
+    * flow control being inserted to implement the short-circuit evaluation
+    * rules.  Flow control is bad!
+    */
+   bool x = !__is_nan(a);
+   bool y = !__is_nan(b);
+   bool z = !__flt64_nonnan(a, b);
 
-   return !__flt64_nonnan(a, b);
+   return (x && y && z);
 }
 
 uint64_t
 __fsat64(uint64_t __a)
 {
-   if (__flt64(__a, 0ul))
+   uvec2 a = unpackUint2x32(__a);
+
+   /* fsat(NaN) should be zero. */
+   if (__is_nan(__a) || int(a.y) < 0)
       return 0ul;
 
-   if (__fge64(__a, 0x3FF0000000000000ul /* 1.0 */))
+   /* IEEE 754 floating point numbers are specifically designed so that, with
+    * two exceptions, values can be compared by bit-casting to signed integers
+    * with the same number of bits.
+    *
+    * From https://en.wikipedia.org/wiki/IEEE_754-1985#Comparing_floating-point_numbers:
+    *
+    *    When comparing as 2's-complement integers: If the sign bits differ,
+    *    the negative number precedes the positive number, so 2's complement
+    *    gives the correct result (except that negative zero and positive zero
+    *    should be considered equal). If both values are positive, the 2's
+    *    complement comparison again gives the correct result. Otherwise (two
+    *    negative numbers), the correct FP ordering is the opposite of the 2's
+    *    complement ordering.
+    *
+    * We know that both values are not negative, and we know that at least one
+    * value is not zero.  Therefore, we can just use the 2's complement
+    * comparison ordering.
+    */
+   if (ilt64(0x3FF00000, 0x00000000, a.y, a.x))
       return 0x3FF0000000000000ul;
 
    return __a;
@@ -578,6 +644,12 @@ __normalizeRoundAndPackFloat64(uint zSign,
 uint64_t
 __propagateFloat64NaN(uint64_t __a, uint64_t __b)
 {
+#if defined RELAXED_NAN_PROPAGATION
+   uvec2 a = unpackUint2x32(__a);
+   uvec2 b = unpackUint2x32(__b);
+
+   return packUint2x32(uvec2(a.x | b.x, a.y | b.y));
+#else
    bool aIsNaN = __is_nan(__a);
    bool bIsNaN = __is_nan(__b);
    uvec2 a = unpackUint2x32(__a);
@@ -586,8 +658,20 @@ __propagateFloat64NaN(uint64_t __a, uint64_t __b)
    b.y |= 0x00080000u;
 
    return packUint2x32(mix(b, mix(a, b, bvec2(bIsNaN, bIsNaN)), bvec2(aIsNaN, aIsNaN)));
+#endif
 }
 
+/* If a shader is in the soft-fp64 path, it almost certainly has register
+ * pressure problems.  Choose a method to exchange two values that does not
+ * require a temporary.
+ */
+#define EXCHANGE(a, b) \
+   do {                \
+       a ^= b;         \
+       b ^= a;         \
+       a ^= b;         \
+   } while (false)
+
 /* Returns the result of adding the double-precision floating-point values
  * `a' and `b'.  The operation is performed according to the IEEE Standard for
  * Floating-Point Arithmetic.
@@ -603,17 +687,16 @@ __fadd64(uint64_t a, uint64_t b)
    uint bFracHi = __extractFloat64FracHi(b);
    int aExp = __extractFloat64Exp(a);
    int bExp = __extractFloat64Exp(b);
-   uint zFrac0 = 0u;
-   uint zFrac1 = 0u;
    int expDiff = aExp - bExp;
    if (aSign == bSign) {
-      uint zFrac2 = 0u;
+      uint zFrac0;
+      uint zFrac1;
+      uint zFrac2;
       int zExp;
-      bool orig_exp_diff_is_zero = (expDiff == 0);
 
-      if (orig_exp_diff_is_zero) {
+      if (expDiff == 0) {
          if (aExp == 0x7FF) {
-            bool propagate = (aFracHi | aFracLo | bFracHi | bFracLo) != 0u;
+            bool propagate = ((aFracHi | bFracHi) | (aFracLo| bFracLo)) != 0u;
             return mix(a, __propagateFloat64NaN(a, b), propagate);
          }
          __add64(aFracHi, aFracLo, bFracHi, bFracLo, zFrac0, zFrac1);
@@ -624,29 +707,34 @@ __fadd64(uint64_t a, uint64_t b)
          zExp = aExp;
          __shift64ExtraRightJamming(
             zFrac0, zFrac1, zFrac2, 1, zFrac0, zFrac1, zFrac2);
-      } else if (0 < expDiff) {
-         if (aExp == 0x7FF) {
-            bool propagate = (aFracHi | aFracLo) != 0u;
-            return mix(a, __propagateFloat64NaN(a, b), propagate);
-         }
+      } else {
+         if (0 < expDiff) {
+            if (aExp == 0x7FF) {
+               bool propagate = (aFracHi | aFracLo) != 0u;
+               return mix(a, __propagateFloat64NaN(a, b), propagate);
+            }
 
-         expDiff = mix(expDiff, expDiff - 1, bExp == 0);
-         bFracHi = mix(bFracHi | 0x00100000u, bFracHi, bExp == 0);
-         __shift64ExtraRightJamming(
-            bFracHi, bFracLo, 0u, expDiff, bFracHi, bFracLo, zFrac2);
-         zExp = aExp;
-      } else if (expDiff < 0) {
-         if (bExp == 0x7FF) {
-            bool propagate = (bFracHi | bFracLo) != 0u;
-            return mix(__packFloat64(aSign, 0x7ff, 0u, 0u), __propagateFloat64NaN(a, b), propagate);
+            expDiff = mix(expDiff, expDiff - 1, bExp == 0);
+            bFracHi = mix(bFracHi | 0x00100000u, bFracHi, bExp == 0);
+            __shift64ExtraRightJamming(
+               bFracHi, bFracLo, 0u, expDiff, bFracHi, bFracLo, zFrac2);
+            zExp = aExp;
+         } else {
+            EXCHANGE(aFracHi, bFracHi);
+            EXCHANGE(aFracLo, bFracLo);
+            EXCHANGE(aExp, bExp);
+
+            if (aExp == 0x7FF) {
+               bool propagate = (aFracHi | aFracLo) != 0u;
+               return mix(__packFloat64(aSign, 0x7ff, 0u, 0u), __propagateFloat64NaN(a, b), propagate);
+            }
+            expDiff = mix(expDiff, expDiff + 1, bExp == 0);
+            bFracHi = mix(bFracHi | 0x00100000u, bFracHi, bExp == 0);
+            __shift64ExtraRightJamming(
+               bFracHi, bFracLo, 0u, - expDiff, bFracHi, bFracLo, zFrac2);
+            zExp = aExp;
          }
-         expDiff = mix(expDiff, expDiff + 1, aExp == 0);
-         aFracHi = mix(aFracHi | 0x00100000u, aFracHi, aExp == 0);
-         __shift64ExtraRightJamming(
-            aFracHi, aFracLo, 0u, - expDiff, aFracHi, aFracLo, zFrac2);
-         zExp = bExp;
-      }
-      if (!orig_exp_diff_is_zero) {
+
          aFracHi |= 0x00100000u;
          __add64(aFracHi, aFracLo, bFracHi, bFracLo, zFrac0, zFrac1);
          --zExp;
@@ -663,6 +751,9 @@ __fadd64(uint64_t a, uint64_t b)
       __shortShift64Left(aFracHi, aFracLo, 10, aFracHi, aFracLo);
       __shortShift64Left(bFracHi, bFracLo, 10, bFracHi, bFracLo);
       if (0 < expDiff) {
+         uint zFrac0;
+         uint zFrac1;
+
          if (aExp == 0x7FF) {
             bool propagate = (aFracHi | aFracLo) != 0u;
             return mix(a, __propagateFloat64NaN(a, b), propagate);
@@ -677,6 +768,9 @@ __fadd64(uint64_t a, uint64_t b)
          return __normalizeRoundAndPackFloat64(aSign, zExp - 10, zFrac0, zFrac1);
       }
       if (expDiff < 0) {
+         uint zFrac0;
+         uint zFrac1;
+
          if (bExp == 0x7FF) {
             bool propagate = (bFracHi | bFracLo) != 0u;
             return mix(__packFloat64(aSign ^ 0x80000000u, 0x7ff, 0u, 0u), __propagateFloat64NaN(a, b), propagate);
@@ -692,36 +786,35 @@ __fadd64(uint64_t a, uint64_t b)
          return __normalizeRoundAndPackFloat64(aSign, zExp - 10, zFrac0, zFrac1);
       }
       if (aExp == 0x7FF) {
-          bool propagate = (aFracHi | aFracLo | bFracHi | bFracLo) != 0u;
+         bool propagate = ((aFracHi | bFracHi) | (aFracLo | bFracLo)) != 0u;
          return mix(0xFFFFFFFFFFFFFFFFUL, __propagateFloat64NaN(a, b), propagate);
       }
       bExp = mix(bExp, 1, aExp == 0);
       aExp = mix(aExp, 1, aExp == 0);
-      bool zexp_normal = false;
-      bool blta = true;
+
+      uint zFrac0;
+      uint zFrac1;
+      uint sign_of_difference = 0;
       if (bFracHi < aFracHi) {
          __sub64(aFracHi, aFracLo, bFracHi, bFracLo, zFrac0, zFrac1);
-         zexp_normal = true;
       }
       else if (aFracHi < bFracHi) {
          __sub64(bFracHi, bFracLo, aFracHi, aFracLo, zFrac0, zFrac1);
-         blta = false;
-         zexp_normal = true;
+         sign_of_difference = 0x80000000;
       }
-      else if (bFracLo < aFracLo) {
+      else if (bFracLo <= aFracLo) {
+         /* It is possible that zFrac0 and zFrac1 may be zero after this. */
          __sub64(aFracHi, aFracLo, bFracHi, bFracLo, zFrac0, zFrac1);
-         zexp_normal = true;
       }
-      else if (aFracLo < bFracLo) {
+      else {
          __sub64(bFracHi, bFracLo, aFracHi, aFracLo, zFrac0, zFrac1);
-          blta = false;
-          zexp_normal = true;
+         sign_of_difference = 0x80000000;
       }
-      zExp = mix(bExp, aExp, blta);
-      aSign = mix(aSign ^ 0x80000000u, aSign, blta);
+      zExp = mix(bExp, aExp, sign_of_difference == 0u);
+      aSign ^= sign_of_difference;
       uint64_t retval_0 = __packFloat64(uint(FLOAT_ROUNDING_MODE == FLOAT_ROUND_DOWN) << 31, 0, 0u, 0u);
       uint64_t retval_1 = __normalizeRoundAndPackFloat64(aSign, zExp - 11, zFrac0, zFrac1);
-      return mix(retval_0, retval_1, zexp_normal);
+      return mix(retval_0, retval_1, zFrac0 != 0u || zFrac1 != 0u);
    }
 }
 
@@ -818,8 +911,13 @@ __fmul64(uint64_t a, uint64_t b)
       return __packFloat64(zSign, 0x7FF, 0u, 0u);
    }
    if (bExp == 0x7FF) {
+      /* a cannot be NaN, but is b NaN? */
       if ((bFracHi | bFracLo) != 0u)
+#if defined RELAXED_NAN_PROPAGATION
+         return b;
+#else
          return __propagateFloat64NaN(a, b);
+#endif
       if ((uint(aExp) | aFracHi | aFracLo) == 0u)
          return 0xFFFFFFFFFFFFFFFFUL;
       return __packFloat64(zSign, 0x7FF, 0u, 0u);
@@ -1636,7 +1734,19 @@ __ftrunc64(uint64_t __a)
 uint64_t
 __ffloor64(uint64_t a)
 {
-   bool is_positive = __fge64(a, 0ul);
+   /* The big assumtion is that when 'a' is NaN, __ftrunc(a) returns a.  Based
+    * on that assumption, NaN values that don't have the sign bit will safely
+    * return NaN (identity).  This is guarded by RELAXED_NAN_PROPAGATION
+    * because otherwise the NaN should have the "signal" bit set.  The
+    * __fadd64 will ensure that occurs.
+    */
+   bool is_positive =
+#if defined RELAXED_NAN_PROPAGATION
+      int(unpackUint2x32(a).y) >= 0
+#else
+      __fge64(a, 0ul)
+#endif
+      ;
    uint64_t tr = __ftrunc64(a);
 
    if (is_positive || __feq64(tr, a)) {
@@ -1694,21 +1804,29 @@ __fround64(uint64_t __a)
 uint64_t
 __fmin64(uint64_t a, uint64_t b)
 {
-   if (__is_nan(a)) return b;
-   if (__is_nan(b)) return a;
+   /* This weird layout matters.  Doing the "obvious" thing results in extra
+    * flow control being inserted to implement the short-circuit evaluation
+    * rules.  Flow control is bad!
+    */
+   bool b_nan = __is_nan(b);
+   bool a_lt_b = __flt64_nonnan(a, b);
+   bool a_nan = __is_nan(a);
 
-   if (__flt64_nonnan(a, b)) return a;
-   return b;
+   return (b_nan || a_lt_b) && !a_nan ? a : b;
 }
 
 uint64_t
 __fmax64(uint64_t a, uint64_t b)
 {
-   if (__is_nan(a)) return b;
-   if (__is_nan(b)) return a;
+   /* This weird layout matters.  Doing the "obvious" thing results in extra
+    * flow control being inserted to implement the short-circuit evaluation
+    * rules.  Flow control is bad!
+    */
+   bool b_nan = __is_nan(b);
+   bool a_lt_b = __flt64_nonnan(a, b);
+   bool a_nan = __is_nan(a);
 
-   if (__flt64_nonnan(a, b)) return b;
-   return a;
+   return (b_nan || a_lt_b) && !a_nan ? b : a;
 }
 
 uint64_t