[flang-commits] [flang] [flang] Subnormal arguments to and results from SPACING (PR #108861)

[Peter Klausler via flang-commits]

https://github.com/klausler created https://github.com/llvm/llvm-project/pull/108861

The standards aren't clear about how IEEE-754 subnormal values interact with the intrinsic function SPACING.  Four compilers interpret the standard such that SPACING(x) will return a value never less than TINY(x); one compiler returns TINY(x) for ABS(x) <= TINY(x) but can return SPACING(x) < TINY(x) for some ABS(x) > TINY(x); one other compiler works similarly, but also oddly returns SPACING(x) < TINY(x) for ABS(x) >= TINY(x)/2.

Follow the most common precedent.

>From 9373f897123f2bfa6ac141b9233624696a123db2 Mon Sep 17 00:00:00 2001
From: Peter Klausler <pklausler at nvidia.com>
Date: Mon, 16 Sep 2024 10:14:23 -0700
Subject: [PATCH] [flang] Subnormal arguments to and results from SPACING

The standards aren't clear about how IEEE-754 subnormal values
interact with the intrinsic function SPACING.  Four compilers
interpret the standard such that SPACING(x) will return a value
never less than TINY(x); one compiler returns TINY(x) for
ABS(x) <= TINY(x) but can return SPACING(x) < TINY(x) for some
ABS(x) > TINY(x); one other compiler works similarly, but also
oddly returns SPACING(x) < TINY(x) for ABS(x) >= TINY(x)/2.

Follow the most common precedent.
---
 flang/lib/Evaluate/real.cpp          |  6 ++++--
 flang/runtime/numeric-templates.h    | 27 ++++++++++++++-------------
 flang/test/Evaluate/fold-spacing.f90 |  4 +++-
 3 files changed, 21 insertions(+), 16 deletions(-)

diff --git a/flang/lib/Evaluate/real.cpp b/flang/lib/Evaluate/real.cpp
index 642490646dcaf3..7ea58297c12e8d 100644
--- a/flang/lib/Evaluate/real.cpp
+++ b/flang/lib/Evaluate/real.cpp
@@ -770,11 +770,13 @@ template <typename W, int P> Real<W, P> Real<W, P>::SPACING() const {
   } else if (IsInfinite()) {
     return NotANumber();
   } else if (IsZero() || IsSubnormal()) {
-    return TINY(); // mandated by standard
+    return TINY(); // standard & 100% portable
   } else {
     Real result;
     result.Normalize(false, Exponent(), Fraction::MASKR(1));
-    return result.IsZero() ? TINY() : result;
+    // Can the result be less than TINY()?  No, with four commonly
+    // used compilers; yes, with two less commonly used ones.
+    return result.IsZero() || result.IsSubnormal() ? TINY() : result;
   }
 }
 
diff --git a/flang/runtime/numeric-templates.h b/flang/runtime/numeric-templates.h
index 1b43498a6bfd12..8cc63859ab7d4f 100644
--- a/flang/runtime/numeric-templates.h
+++ b/flang/runtime/numeric-templates.h
@@ -88,12 +88,16 @@ struct MaxOrMinIdentity<TypeCategory::Real, 16, IS_MAXVAL,
 
 // Minimum finite representable value.
 // For floating-point types, returns minimum positive normalized value.
-template <typename T> struct MinValue {
+template <int PREC, typename T> struct MinValue {
   static RT_API_ATTRS T get() { return std::numeric_limits<T>::min(); }
 };
+template <typename T> struct MinValue<11, T> {
+  // TINY(0._2)
+  static constexpr RT_API_ATTRS T get() { return 0.00006103515625E-04; }
+};
 
 #if HAS_FLOAT128
-template <> struct MinValue<CppTypeFor<TypeCategory::Real, 16>> {
+template <> struct MinValue<113, CppTypeFor<TypeCategory::Real, 16>> {
   using Type = CppTypeFor<TypeCategory::Real, 16>;
   static RT_API_ATTRS Type get() {
     // Create a buffer to store binary representation of __float128 constant.
@@ -167,8 +171,8 @@ template <> struct MAXTy<CppTypeFor<TypeCategory::Real, 16>> {
 };
 #endif
 
-template <typename T> struct MINTy {
-  static constexpr RT_API_ATTRS T compute() { return MinValue<T>::get(); }
+template <int PREC, typename T> struct MINTy {
+  static constexpr RT_API_ATTRS T compute() { return MinValue<PREC, T>::get(); }
 };
 
 template <typename T> struct QNANTy {
@@ -339,23 +343,20 @@ template <int PREC, typename T> inline RT_API_ATTRS T RRSpacing(T x) {
 
 // SPACING (16.9.180)
 template <int PREC, typename T> inline RT_API_ATTRS T Spacing(T x) {
+  T tiny{MINTy<PREC, T>::compute()};
   if (ISNANTy<T>::compute(x)) {
     return x; // NaN -> same NaN
   } else if (ISINFTy<T>::compute(x)) {
     return QNANTy<T>::compute(); // +/-Inf -> NaN
   } else if (x == 0) { // 0 -> TINY(x)
-    // The standard-mandated behavior seems broken, since TINY() can't be
-    // subnormal.
-    if constexpr (PREC == 11) { // REAL(2)
-      return 0.00006103515625E-04; // TINY(0._2)
-    } else {
-      // N.B. TINY(0._3) == TINY(0._4) so this works even if no std::bfloat16_t.
-      return MINTy<T>::compute();
-    }
+    return tiny;
   } else {
     T result{LDEXPTy<T>::compute(
         static_cast<T>(1.0), ILOGBTy<T>::compute(x) + 1 - PREC)}; // 2**(e-p)
-    return result == 0 ? /*TINY(x)*/ MINTy<T>::compute() : result;
+    // All compilers return TINY(x) for |x| <= TINY(x), but differ over whether
+    // SPACING(x) can be < TINY(x) for |x| > TINY(x).  The most common precedent
+    // is to never return a value < TINY(x).
+    return result <= tiny ? tiny : result;
   }
 }
 
diff --git a/flang/test/Evaluate/fold-spacing.f90 b/flang/test/Evaluate/fold-spacing.f90
index 1d7b58081d70fe..cfae6eab8b9b8b 100644
--- a/flang/test/Evaluate/fold-spacing.f90
+++ b/flang/test/Evaluate/fold-spacing.f90
@@ -5,7 +5,9 @@ module m
   logical, parameter :: test_2 = spacing(-3.0) == scale(1.0, -22)
   logical, parameter :: test_3 = spacing(3.0d0) == scale(1.0, -51)
   logical, parameter :: test_4 = spacing(0.) == tiny(0.)
-  logical, parameter :: test_5 = spacing(tiny(0.)) == 1.e-45
+  logical, parameter :: test_5a = spacing(tiny(0.)) == tiny(0.)
+  logical, parameter :: test_5b = spacing(tiny(0.)/2) == tiny(0.)
+  logical, parameter :: test_5c = spacing(tiny(0.)*2) == tiny(0.)
   logical, parameter :: test_6 = spacing(8388608.) == 1.
   logical, parameter :: test_7 = spacing(spacing(tiny(.0))) == tiny(0.)
   logical, parameter :: test_11 = rrspacing(3.0) == scale(0.75, 24)