[flang] More real work. All ops, rounding modes, and flags good except for division.

Original-commit: flang-compiler/f18@ea697295db84b64b1c13144e3b5eaef9a879e319
Reviewed-on: https://github.com/flang-compiler/f18/pull/101
Tree-same-pre-rewrite: false

diff --git a/flang/lib/evaluate/real.h b/flang/lib/evaluate/real.h
index bf69a9f..4ac465d 100644 (file)
--- a/flang/lib/evaluate/real.h
+++ b/flang/lib/evaluate/real.h
@@ -267,8 +267,11 @@ public:
         return y.Add(*this, rounding);
       }
       if (order == Ordering::Equal) {
-        // x + (-x) -> +0.0, never -0.0
-        return {};
+        // x + (-x) -> +0.0 unless rounding is directed downwards
+        if (rounding == Rounding::Down) {
+          result.value.Normalize(true, 0, Fraction{}, rounding);  // -0.0
+        }
+        return result;
       }
     }
     // Our exponent is greater than y's, or the exponents match and y is not
@@ -296,8 +299,8 @@ public:
       fraction = fraction.SHIFTR(1).IBSET(fraction.bits - 1);
       ++exponent;
     }
-    result.flags |=
-        result.value.Normalize(isNegative, exponent, fraction, &roundingBits);
+    result.flags |= result.value.Normalize(
+        isNegative, exponent, fraction, rounding, &roundingBits);
     result.flags |= result.value.Round(rounding, roundingBits);
     return result;
   }
@@ -312,7 +315,9 @@ public:
     ValueWithRealFlags<Real> result;
     if (IsNotANumber() || y.IsNotANumber()) {
       result.value.word_ = NaNWord();  // NaN * x -> NaN
-      result.flags.set(RealFlag::InvalidArgument);
+      if (IsSignalingNaN() || y.IsSignalingNaN()) {
+        result.flags.set(RealFlag::InvalidArgument);
+      }
     } else {
       bool isNegative{IsNegative() != y.IsNegative()};
       if (IsInfinite() || y.IsInfinite()) {
@@ -338,7 +343,10 @@ public:
           if (rshift >= product.upper.bits + product.lower.bits) {
             sticky = !product.lower.IsZero() || !product.upper.IsZero();
           } else if (rshift >= product.lower.bits) {
-            sticky = !product.lower.IsZero();
+            sticky = !product.lower.IsZero() ||
+                !product.upper
+                     .IAND(product.upper.MASKR(rshift - product.lower.bits))
+                     .IsZero();
           } else {
             sticky = !product.lower.IAND(product.lower.MASKR(rshift)).IsZero();
           }
@@ -359,10 +367,14 @@ public:
         exponent -= lshift;
         product.upper = product.upper.DSHIFTL(product.lower, lshift);
         product.lower = product.lower.SHIFTL(lshift);
-        RoundingBits roundingBits{product.lower, product.upper.bits};
+        RoundingBits roundingBits{product.lower, product.lower.bits};
         result.flags |= result.value.Normalize(
-            isNegative, exponent, product.upper, &roundingBits);
+            isNegative, exponent, product.upper, rounding, &roundingBits);
         result.flags |= result.value.Round(rounding, roundingBits);
+        if (result.flags.test(RealFlag::Inexact) &&
+            result.value.Exponent() == 0) {
+          result.flags.set(RealFlag::Underflow);
+        }
       }
     }
     return result;
@@ -514,9 +526,19 @@ private:
   }
 
   constexpr RealFlags Normalize(bool negative, std::uint64_t exponent,
-      const Fraction &fraction, RoundingBits *roundingBits = nullptr) {
+      const Fraction &fraction, Rounding rounding = Rounding::TiesToEven,
+      RoundingBits *roundingBits = nullptr) {
     if (exponent >= maxExponent) {
-      word_ = Word{maxExponent}.SHIFTL(significandBits);  // Inf
+      if (rounding == Rounding::TiesToEven ||
+          rounding == Rounding::TiesAwayFromZero ||
+          (rounding == Rounding::Up && !negative) ||
+          (rounding == Rounding::Down && negative)) {
+        word_ = Word{maxExponent}.SHIFTL(significandBits);  // Inf
+      } else {
+        // directed rounding: round to largest finite value rather than infinity
+        // (x86 does this, not sure whether it's standard or not)
+        word_ = Word{word_.MASKR(word_.bits - 1)}.IBCLR(significandBits);
+      }
       if (negative) {
         word_ = word_.IBSET(bits - 1);
       }
@@ -608,7 +630,7 @@ extern template class Real<Integer<32>, 24>;
 using RealKind8 = Real<Integer<64>, 53>;
 extern template class Real<Integer<64>, 53>;
 
-using RealKind10 = Real<Integer<80>, 64, false>;  // 80387
+using RealKind10 = Real<Integer<80>, 64, false>;  // 80387 extended precision
 extern template class Real<Integer<80>, 64, false>;
 
 using RealKind16 = Real<Integer<128>, 112>;

diff --git a/flang/test/evaluate/fp-testing.cc b/flang/test/evaluate/fp-testing.cc
index a4cbc9e..5dc718f 100644 (file)
--- a/flang/test/evaluate/fp-testing.cc
+++ b/flang/test/evaluate/fp-testing.cc
@@ -63,7 +63,7 @@ void ScopedHostFloatingPointEnvironment::ClearFlags() const {
   feclearexcept(FE_ALL_EXCEPT);
 }
 
-RealFlags ScopedHostFloatingPointEnvironment::CurrentFlags() const {
+RealFlags ScopedHostFloatingPointEnvironment::CurrentFlags() {
   int exceptions = fetestexcept(FE_ALL_EXCEPT);
   RealFlags flags;
   if (exceptions & FE_INVALID) {
@@ -83,3 +83,16 @@ RealFlags ScopedHostFloatingPointEnvironment::CurrentFlags() const {
   }
   return flags;
 }
+
+void ScopedHostFloatingPointEnvironment::SetRounding(Rounding rounding) {
+  switch (rounding) {
+  case Rounding::TiesToEven: fesetround(FE_TONEAREST); break;
+  case Rounding::ToZero: fesetround(FE_TOWARDZERO); break;
+  case Rounding::Up: fesetround(FE_UPWARD); break;
+  case Rounding::Down: fesetround(FE_DOWNWARD); break;
+  case Rounding::TiesAwayFromZero:
+    std::fprintf(stderr, "SetRounding: TiesAwayFromZero not available");
+    std::abort();
+    break;
+  }
+}

diff --git a/flang/test/evaluate/fp-testing.h b/flang/test/evaluate/fp-testing.h
index 8ee53ad..afa03ed 100644 (file)
--- a/flang/test/evaluate/fp-testing.h
+++ b/flang/test/evaluate/fp-testing.h
@@ -19,6 +19,7 @@
 #include <fenv.h>
 
 using Fortran::evaluate::RealFlags;
+using Fortran::evaluate::Rounding;
 
 class ScopedHostFloatingPointEnvironment {
 public:
@@ -26,7 +27,8 @@ public:
       bool flushDenormalResultsToZero = false);
   ~ScopedHostFloatingPointEnvironment();
   void ClearFlags() const;
-  RealFlags CurrentFlags() const;
+  static RealFlags CurrentFlags();
+  static void SetRounding(Rounding rounding);
 
 private:
   fenv_t originalFenv_;

diff --git a/flang/test/evaluate/real.cc b/flang/test/evaluate/real.cc
index b93e793..2507c15 100644 (file)
--- a/flang/test/evaluate/real.cc
+++ b/flang/test/evaluate/real.cc
@@ -20,7 +20,7 @@
 
 using namespace Fortran::evaluate;
 
-template<typename R> void tests() {
+template<typename R> void basicTests(int rm, Rounding rounding) {
   char desc[64];
   using Word = typename R::Word;
   std::snprintf(
@@ -99,51 +99,52 @@ template<typename R> void tests() {
   TEST(inf.Compare(negInf) == Relation::Greater)(desc);
   TEST(negInf.Compare(negInf) == Relation::Equal)(desc);
   for (std::uint64_t j{0}; j < 63; ++j) {
+    char ldesc[128];
     std::uint64_t x{1};
     x <<= j;
+    std::snprintf(ldesc, sizeof ldesc, "%s j=%d x=0x%llx rm=%d", desc,
+        static_cast<int>(j), static_cast<unsigned long long>(x), rm);
     Integer<64> ix{x};
-    TEST(!ix.IsNegative())("%s,%d,0x%llx", desc, j, x);
-    MATCH(x, ix.ToUInt64())("%s,%d,0x%llx", desc, j, x);
-    vr = R::ConvertSigned(ix);
-    TEST(!vr.value.IsNegative())("%s,%d,0x%llx", desc, j, x);
-    TEST(!vr.value.IsNotANumber())("%s,%d,0x%llx", desc, j, x);
-    TEST(!vr.value.IsZero())("%s,%d,0x%llx", desc, j, x);
+    TEST(!ix.IsNegative())(ldesc);
+    MATCH(x, ix.ToUInt64())(ldesc);
+    vr = R::ConvertSigned(ix, rounding);
+    TEST(!vr.value.IsNegative())(ldesc);
+    TEST(!vr.value.IsNotANumber())(ldesc);
+    TEST(!vr.value.IsZero())(ldesc);
     auto ivf = vr.value.template ToInteger<Integer<64>>();
     if (j > (maxExponent / 2)) {
-      TEST(vr.flags.test(RealFlag::Overflow))(desc);
-      TEST(vr.value.IsInfinite())("%s,%d,0x%llx", desc, j, x);
-      TEST(ivf.flags.test(RealFlag::Overflow))("%s,%d,0x%llx", desc, j, x);
-      MATCH(0x7fffffffffffffff, ivf.value.ToUInt64())
-      ("%s,%d,0x%llx", desc, j, x);
+      TEST(vr.flags.test(RealFlag::Overflow))(ldesc);
+      TEST(vr.value.IsInfinite())(ldesc);
+      TEST(ivf.flags.test(RealFlag::Overflow))(ldesc);
+      MATCH(0x7fffffffffffffff, ivf.value.ToUInt64())(ldesc);
     } else {
-      TEST(vr.flags.empty())(desc);
-      TEST(!vr.value.IsInfinite())("%s,%d,0x%llx", desc, j, x);
-      TEST(ivf.flags.empty())("%s,%d,0x%llx", desc, j, x);
-      MATCH(x, ivf.value.ToUInt64())("%s,%d,0x%llx", desc, j, x);
+      TEST(vr.flags.empty())(ldesc);
+      TEST(!vr.value.IsInfinite())(ldesc);
+      TEST(ivf.flags.empty())(ldesc);
+      MATCH(x, ivf.value.ToUInt64())(ldesc);
     }
     ix = ix.Negate().value;
-    TEST(ix.IsNegative())("%s,%d,0x%llx", desc, j, x);
+    TEST(ix.IsNegative())(ldesc);
     x = -x;
     std::int64_t nx = x;
-    MATCH(x, ix.ToUInt64())("%s,%d,0x%llx", desc, j, x);
-    MATCH(nx, ix.ToInt64())("%s,%d,0x%llx", desc, j, x);
+    MATCH(x, ix.ToUInt64())(ldesc);
+    MATCH(nx, ix.ToInt64())(ldesc);
     vr = R::ConvertSigned(ix);
-    TEST(vr.value.IsNegative())("%s,%d,0x%llx", desc, j, x);
-    TEST(!vr.value.IsNotANumber())("%s,%d,0x%llx", desc, j, x);
-    TEST(!vr.value.IsZero())("%s,%d,0x%llx", desc, j, x);
+    TEST(vr.value.IsNegative())(ldesc);
+    TEST(!vr.value.IsNotANumber())(ldesc);
+    TEST(!vr.value.IsZero())(ldesc);
     ivf = vr.value.template ToInteger<Integer<64>>();
     if (j > (maxExponent / 2)) {
-      TEST(vr.flags.test(RealFlag::Overflow))(desc);
-      TEST(vr.value.IsInfinite())("%s,%d,0x%llx", desc, j, x);
-      TEST(ivf.flags.test(RealFlag::Overflow))("%s,%d,0x%llx", desc, j, x);
-      MATCH(0x8000000000000000, ivf.value.ToUInt64())
-      ("%s,%d,0x%llx", desc, j, x);
+      TEST(vr.flags.test(RealFlag::Overflow))(ldesc);
+      TEST(vr.value.IsInfinite())(ldesc);
+      TEST(ivf.flags.test(RealFlag::Overflow))(ldesc);
+      MATCH(0x8000000000000000, ivf.value.ToUInt64())(ldesc);
     } else {
-      TEST(vr.flags.empty())(desc);
-      TEST(!vr.value.IsInfinite())("%s,%d,0x%llx", desc, j, x);
-      TEST(ivf.flags.empty())("%s,%d,0x%llx", desc, j, x);
-      MATCH(x, ivf.value.ToUInt64())("%s,%d,0x%llx", desc, j, x);
-      MATCH(nx, ivf.value.ToInt64())("%s,%d,0x%llx", desc, j, x);
+      TEST(vr.flags.empty())(ldesc);
+      TEST(!vr.value.IsInfinite())(ldesc);
+      TEST(ivf.flags.empty())(ldesc);
+      MATCH(x, ivf.value.ToUInt64())(ldesc);
+      MATCH(nx, ivf.value.ToInt64())(ldesc);
     }
   }
 }
@@ -183,12 +184,41 @@ std::uint32_t FlagsToBits(const RealFlags &flags) {
   return bits;
 }
 
-void subset32bit() {
+void inttest(std::int64_t x, int pass, Rounding rounding) {
+  union {
+    std::uint32_t u32;
+    float f;
+  } u;
+  ScopedHostFloatingPointEnvironment fpenv;
+  Integer<64> ix{x};
+  ValueWithRealFlags<RealKind4> real;
+  real = real.value.ConvertSigned(ix, rounding);
+  fpenv.ClearFlags();
+  float fcheck = x;  // TODO unsigned too
+  auto actualFlags{FlagsToBits(fpenv.CurrentFlags())};
+  u.f = fcheck;
+  std::uint32_t rcheck{NormalizeNaN(u.u32)};
+  std::uint32_t check = real.value.RawBits().ToUInt64();
+  MATCH(rcheck, check)("%d 0x%llx", pass, x);
+  MATCH(actualFlags, FlagsToBits(real.flags))("%d 0x%llx", pass, x);
+}
+
+void subset32bit(int pass, Rounding rounding) {
+  for (int j{0}; j < 63; ++j) {
+    std::int64_t x{1};
+    x <<= j;
+    inttest(x, pass, rounding);
+    inttest(-x, pass, rounding);
+  }
+  inttest(0, pass, rounding);
+  inttest(static_cast<std::int64_t>(0x8000000000000000), pass, rounding);
+
   union {
     std::uint32_t u32;
     float f;
   } u;
   ScopedHostFloatingPointEnvironment fpenv;
+
   for (std::uint32_t j{0}; j < 8192; ++j) {
     std::uint32_t rj{MakeReal(j)};
     u.u32 = rj;
@@ -200,51 +230,71 @@ void subset32bit() {
       float fk{u.f};
       RealKind4 y{Integer<32>{std::uint64_t{rk}}};
       {
-        ValueWithRealFlags<RealKind4> sum{x.Add(y)};
+        ValueWithRealFlags<RealKind4> sum{x.Add(y, rounding)};
         fpenv.ClearFlags();
         float fcheck{fj + fk};
         auto actualFlags{FlagsToBits(fpenv.CurrentFlags())};
         u.f = fcheck;
         std::uint32_t rcheck{NormalizeNaN(u.u32)};
         std::uint32_t check = sum.value.RawBits().ToUInt64();
-        MATCH(rcheck, check)("0x%x + 0x%x", rj, rk);
-        MATCH(actualFlags, FlagsToBits(sum.flags))("0x%x + 0x%x", rj, rk);
+        MATCH(rcheck, check)("%d 0x%x + 0x%x", pass, rj, rk);
+        MATCH(actualFlags, FlagsToBits(sum.flags))
+        ("%d 0x%x + 0x%x", pass, rj, rk);
       }
       {
-        ValueWithRealFlags<RealKind4> diff{x.Subtract(y)};
+        ValueWithRealFlags<RealKind4> diff{x.Subtract(y, rounding)};
+        fpenv.ClearFlags();
         float fcheck{fj - fk};
+        auto actualFlags{FlagsToBits(fpenv.CurrentFlags())};
         u.f = fcheck;
         std::uint32_t rcheck{NormalizeNaN(u.u32)};
         std::uint32_t check = diff.value.RawBits().ToUInt64();
-        MATCH(rcheck, check)("0x%x - 0x%x", rj, rk);
+        MATCH(rcheck, check)("%d 0x%x - 0x%x", pass, rj, rk);
+        MATCH(actualFlags, FlagsToBits(diff.flags))
+        ("%d 0x%x - 0x%x", pass, rj, rk);
       }
       {
-        ValueWithRealFlags<RealKind4> prod{x.Multiply(y)};
+        ValueWithRealFlags<RealKind4> prod{x.Multiply(y, rounding)};
+        fpenv.ClearFlags();
         float fcheck{fj * fk};
+        auto actualFlags{FlagsToBits(fpenv.CurrentFlags())};
         u.f = fcheck;
         std::uint32_t rcheck{NormalizeNaN(u.u32)};
         std::uint32_t check = prod.value.RawBits().ToUInt64();
-        MATCH(rcheck, check)("0x%x * 0x%x", rj, rk);
+        MATCH(rcheck, check)("%d 0x%x * 0x%x", pass, rj, rk);
+        MATCH(actualFlags, FlagsToBits(prod.flags))
+        ("%d 0x%x * 0x%x -> 0x%x", pass, rj, rk, rcheck);
       }
 #if 0
-      { ValueWithRealFlags<RealKind4> quot{x.Divide(y)};
-        float fcheck{fj * fk};
+      { ValueWithRealFlags<RealKind4> quot{x.Divide(y, rounding)};
+        fpenv.ClearFlags();
+        float fcheck{fj / fk};
+        auto actualFlags{FlagsToBits(fpenv.CurrentFlags())};
         u.f = fcheck;
         std::uint32_t rcheck{NormalizeNaN(u.u32)};
         std::uint32_t check = quot.value.RawBits().ToUInt64();
-        MATCH(rcheck, check)("0x%x / 0x%x", rj, rk);
+        MATCH(rcheck, check)("%d 0x%x / 0x%x", pass, rj, rk);
+        MATCH(actualFlags, FlagsToBits(quot.flags))("%d 0x%x / 0x%x", pass, rj, rk);
       }
 #endif
     }
   }
 }
 
+void roundTest(int rm, Rounding rounding) {
+  basicTests<RealKind2>(rm, rounding);
+  basicTests<RealKind4>(rm, rounding);
+  basicTests<RealKind8>(rm, rounding);
+  basicTests<RealKind10>(rm, rounding);
+  basicTests<RealKind16>(rm, rounding);
+  ScopedHostFloatingPointEnvironment::SetRounding(rounding);
+  subset32bit(rm, rounding);
+}
+
 int main() {
-  tests<RealKind2>();
-  tests<RealKind4>();
-  tests<RealKind8>();
-  tests<RealKind10>();
-  tests<RealKind16>();
-  subset32bit();  // TODO rounding modes
-  return testing::Complete();
+  roundTest(0, Rounding::TiesToEven);
+  roundTest(1, Rounding::ToZero);
+  roundTest(2, Rounding::Up);
+  roundTest(3, Rounding::Down);
+  // TODO: how to test Rounding::TiesAwayFromZero on x86?
 }