[flang][runtime] Fixed identity value for REAL(16) == __float128.

[vzakhari]

Summary

std::numeric_limits<__float128>::max/lowest return 0.0, so recreate
value of FLT128_MAX ourselves to avoid using quadmath.h's FLT128_MAX
that is currently causes warnings with GCC -Wpedantic.

[klausler]

I consider it a bug in the Fortran standard that the identity values for real minval/maxval are not signed infinities on IEEE-754 targets.

Perhaps you could compute the largest finite number with scalbn or ldexp applied to nexttoward((quad precision)1.0, 0.0) on hosts whose numeric_limits header is broken.

[vzakhari]

Thanks!

I will use ldexpf128(nextafterf128(1.0, 0.0), __FLT128_MAX_EXP__). These functions are defined by glibc headers for GCC compiler. I could have used GCC's nextafterq/ldexpq, but I think it is better to use glibc names so that when they add support for clang it will hopefully work as-is. In any case, I will have to check in CMake file that these entry points are available and only then enable the REAL(16) specialization. Runtime will not work correctly, when built with clang, and the unittest will fail, so I will probably have to make the same checks in CMake files for the unittests.

This will allow getting rid of these weird -Wpedantic workaround.

[klausler]

You can also just initialize the value with a union with a 128-bit integer, using 0x7ffeffffffffffff for the most significant 64 bits and 0xffffffffffffffff for the least significant 64 bits. This would be even better since that's a known permanent constant.

[vzakhari]

In D134496#3810602, @klausler wrote:

You can also just initialize the value with a union with a 128-bit integer, using 0x7ffeffffffffffff for the most significant 64 bits and 0xffffffffffffffff for the least significant 64 bits. This would be even better since that's a known permanent constant.

This is a possibility, but I am not sure if it is safe to rely on particular __float128 representation. Is it guaranteed to comply with IEEE binary128 format in all implementations? I know we only set HAS_FLOAT128 to true for clang and gcc and I guess it should be a safe assumption, but I am not quite comfortable with the hard-coded constant.

Also, can there be different endianness for integer and floating point numbers on a target?

[klausler]

In D134496#3810661, @vzakhari wrote:

In D134496#3810602, @klausler wrote:

You can also just initialize the value with a union with a 128-bit integer, using 0x7ffeffffffffffff for the most significant 64 bits and 0xffffffffffffffff for the least significant 64 bits. This would be even better since that's a known permanent constant.

This is a possibility, but I am not sure if it is safe to rely on particular __float128 representation. Is it guaranteed to comply with IEEE binary128 format in all implementations? I know we only set HAS_FLOAT128 to true for clang and gcc and I guess it should be a safe assumption, but I am not quite comfortable with the hard-coded constant.

Also, can there be different endianness for integer and floating point numbers on a target?

Yes, there is only one IEEE-754 representation for each size of numbers; implementations may vary their encodings only for signaling vs quiet NaNs. Yes, the PDP-11 and successors had a byte order for floating-point numbers that was neither little-endian nor big-endian.

If you use an expression to compute a known constant value, please at least check the result to ensure that it equals the known constant value.

[vzakhari]

In D134496#3811756, @klausler wrote:

In D134496#3810661, @vzakhari wrote:

In D134496#3810602, @klausler wrote:

You can also just initialize the value with a union with a 128-bit integer, using 0x7ffeffffffffffff for the most significant 64 bits and 0xffffffffffffffff for the least significant 64 bits. This would be even better since that's a known permanent constant.

This is a possibility, but I am not sure if it is safe to rely on particular __float128 representation. Is it guaranteed to comply with IEEE binary128 format in all implementations? I know we only set HAS_FLOAT128 to true for clang and gcc and I guess it should be a safe assumption, but I am not quite comfortable with the hard-coded constant.

Also, can there be different endianness for integer and floating point numbers on a target?

Yes, there is only one IEEE-754 representation for each size of numbers; implementations may vary their encodings only for signaling vs quiet NaNs. Yes, the PDP-11 and successors had a byte order for floating-point numbers that was neither little-endian nor big-endian.

If you use an expression to compute a known constant value, please at least check the result to ensure that it equals the known constant value.

Yes, it produces expected result on x86/arm/power with gcc.

[klausler]

In D134496#3811945, @vzakhari wrote:

In D134496#3811756, @klausler wrote:

In D134496#3810661, @vzakhari wrote:

In D134496#3810602, @klausler wrote:

You can also just initialize the value with a union with a 128-bit integer, using 0x7ffeffffffffffff for the most significant 64 bits and 0xffffffffffffffff for the least significant 64 bits. This would be even better since that's a known permanent constant.

This is a possibility, but I am not sure if it is safe to rely on particular __float128 representation. Is it guaranteed to comply with IEEE binary128 format in all implementations? I know we only set HAS_FLOAT128 to true for clang and gcc and I guess it should be a safe assumption, but I am not quite comfortable with the hard-coded constant.

Also, can there be different endianness for integer and floating point numbers on a target?

Yes, there is only one IEEE-754 representation for each size of numbers; implementations may vary their encodings only for signaling vs quiet NaNs. Yes, the PDP-11 and successors had a byte order for floating-point numbers that was neither little-endian nor big-endian.

If you use an expression to compute a known constant value, please at least check the result to ensure that it equals the known constant value.

Yes, it produces expected result on x86/arm/power with gcc.

I mean that you should check it in code with an assert.

[vzakhari]

But then I will have to use FLT128_MAX and quadmath.h, and I am back at the weird -Wpedantic issue :) can you please clarify how you suggest doing this check?

[klausler]

In D134496#3812138, @vzakhari wrote:

But then I will have to use FLT128_MAX and quadmath.h, and I am back at the weird -Wpedantic issue :) can you please clarify how you suggest doing this check?

Use hex.

[klausler]

In D134496#3812151, @klausler wrote:

In D134496#3812138, @vzakhari wrote:

But then I will have to use FLT128_MAX and quadmath.h, and I am back at the weird -Wpedantic issue :) can you please clarify how you suggest doing this check?

Use hex.

Or verify that the value is finite and that adding one ULP to it makes it infinite (or equivalently raises an overflow).

[vzakhari]

Ok, then I will just use hex with some sanity verification.

flang/runtime/CMakeLists.txt

check_cxx_source_compiles(
    }
   "
   HAVE_DECL_STRERROR_S)
 
 if (NOT (HAVE_STRERROR OR HAVE_STRERROR_R OR HAVE_DECL_STRERROR_S))
   message(FATAL_ERROR "None of strerror, strerror_r, strerror_s found.")
 endif()
 
+include(CheckCXXCompilerFlag)
+check_cxx_compiler_flag(-fext-numeric-literals SUPPORTS_FEXT_NUMERIC_LITERALS)
+check_cxx_compiler_flag(-Wno-error=pedantic SUPPORTS_W_NO_ERROR_PEDANTIC)
+if (SUPPORTS_FEXT_NUMERIC_LITERALS AND SUPPORTS_W_NO_ERROR_PEDANTIC)
+  set(CMAKE_REQUIRED_FLAGS
+    "${CMAKE_REQUIRED_FLAGS} -Wno-error=pedantic -fext-numeric-literals")
+  check_cxx_source_compiles(
+    "#include <quadmath.h>
+     int main() {
+       __float128 x = FLT128_MAX;
+       return 0;
+     }
+    "
+    HAVE_QUADMATH_H)
+
+  # Even with -fext-numeric-literals, GCC -Wpedantic warns about
+  # non-standard suffix on floating constant expanded from FLT128_MAX.
+  # Neither __extension__ nor 'GCC diagnostic ignored "-Wpedantic"' work,
+  # so usage of FLT128_MAX will cause warnings with -Wpedantic.
+  # -Werror will turn this warning into an error, so we have to
+  # override it with -Wno-error=pedantic both here and for the files
+  # that use FLT128_MAX (e.g. extrema.cpp).
+  # Also see: https://gcc.gnu.org/bugzilla/show_bug.cgi?id=87274
+  if (HAVE_QUADMATH_H)
+    set(files_with_quadmath extrema.cpp)
+    set_source_files_properties(${files_with_quadmath}
+      PROPERTIES COMPILE_OPTIONS "-Wno-error=pedantic;-fext-numeric-literals")
+    set_source_files_properties(${files_with_quadmath}
+      PROPERTIES COMPILE_DEFINITIONS "HAVE_QUADMATH_H")
+  endif()
+endif()
+
 configure_file(config.h.cmake config.h)
 # include_directories is used here instead of target_include_directories
 # because add_flang_library creates multiple objects (STATIC/SHARED, OBJECT)
 # with different names
 include_directories(AFTER ${CMAKE_CURRENT_BINARY_DIR})
 
 add_subdirectory(FortranMain)
 

flang/runtime/extrema.cpp

 #include "flang/Runtime/character.h"
 #include "flang/Runtime/float128.h"
 #include "flang/Runtime/reduction.h"
 #include <algorithm>
 #include <cfloat>
 #include <cinttypes>
 #include <cmath>
 #include <optional>
+#if HAVE_QUADMATH_H
+#include <quadmath.h>
+#endif // HAVE_QUADMATH_H
 
 namespace Fortran::runtime {
 
 // MAXLOC & MINLOC
 
 template <typename T, bool IS_MAX, bool BACK> struct NumericCompare {
   using Type = T;
   explicit NumericCompare(std::size_t /*elemLen; ignored*/) {}
void RTNAME(MinlocDim)(Descriptor &result, const Descriptor &x, int kind,
     int dim, const char *source, int line, const Descriptor *mask, bool back) {
   TypedPartialMaxOrMinLoc<false>(
       "MINLOC", result, x, kind, dim, source, line, mask, back);
 }
 } // extern "C"
 
 // MAXVAL and MINVAL
 
-template <TypeCategory CAT, int KIND, bool IS_MAXVAL> struct MaxOrMinIdentity {
+template <TypeCategory CAT, int KIND, bool IS_MAXVAL, typename Enable = void>
+struct MaxOrMinIdentity {
   using Type = CppTypeFor<CAT, KIND>;
   static constexpr Type Value() {
     return IS_MAXVAL ? std::numeric_limits<Type>::lowest()
                      : std::numeric_limits<Type>::max();
   }
 };
 
 // std::numeric_limits<> may not know int128_t
 template <bool IS_MAXVAL>
 struct MaxOrMinIdentity<TypeCategory::Integer, 16, IS_MAXVAL> {
   using Type = CppTypeFor<TypeCategory::Integer, 16>;
   static constexpr Type Value() {
     return IS_MAXVAL ? Type{1} << 127 : ~Type{0} >> 1;
   }
 };
 
+#if HAVE_QUADMATH_H
+// std::numeric_limits<> may not support __float128,
+// so use macros defined in quadmath.h (if available).
+//
+// Note that HAVE_QUADMATH_H implies that __float128 is available,
+// but it does not necessarily mean that __float128 is used for REAL(16).
+// This specialization must engage only when
+// CppTypeFor<TypeCategory::Real, 16> is __float128.
+template <bool IS_MAXVAL>
+struct MaxOrMinIdentity<TypeCategory::Real, 16, IS_MAXVAL,
+    typename std::enable_if_t<
+        std::is_same_v<CppTypeFor<TypeCategory::Real, 16>, __float128>>> {
+  using Type = CppTypeFor<TypeCategory::Real, 16>;
+  static constexpr Type Value() { return IS_MAXVAL ? -FLT128_MAX : FLT128_MAX; }
+};
+#endif // HAVE_QUADMATH_H
+
 template <TypeCategory CAT, int KIND, bool IS_MAXVAL>
 class NumericExtremumAccumulator {
 public:
   using Type = CppTypeFor<CAT, KIND>;
   explicit NumericExtremumAccumulator(const Descriptor &array)
       : array_{array} {}
   void Reinitialize() {
     extremum_ = MaxOrMinIdentity<CAT, KIND, IS_MAXVAL>::Value();

flang/unittests/Runtime/Reduction.cpp

auto logicalVector2{MakeArray<TypeCategory::Logical, 1>(
       std::vector<int>{4}, std::vector<bool>{true, true, false, false})};
   EXPECT_TRUE(RTNAME(DotProductLogical)(
       *logicalVector2, *logicalVector2, __FILE__, __LINE__));
   EXPECT_FALSE(RTNAME(DotProductLogical)(
       *logicalVector1, *logicalVector2, __FILE__, __LINE__));
   EXPECT_FALSE(RTNAME(DotProductLogical)(
       *logicalVector2, *logicalVector1, __FILE__, __LINE__));
 }
+
+#if LDBL_MANT_DIG == 113 || HAS_FLOAT128
+TEST(Reductions, ExtremaReal16) {
+  // The identity value for Min/Maxval for REAL(16) was mistakenly
+  // set to 0.0.
+  using ElemType = CppTypeFor<TypeCategory::Real, 16>;
+  std::vector<int> shape{3};
+  //   1.0  2.0  3.0
+  std::vector<ElemType> rawMinData{1.0, 2.0, 3.0};
+  auto minArray{MakeArray<TypeCategory::Real, 16>(shape, rawMinData)};
+  EXPECT_EQ(RTNAME(MinvalReal16)(*minArray, __FILE__, __LINE__), 1.0);
+  //   -1.0  -2.0  -3.0
+  std::vector<ElemType> rawMaxData{-1.0, -2.0, -3.0};
+  auto maxArray{MakeArray<TypeCategory::Real, 16>(shape, rawMaxData)};
+  EXPECT_EQ(RTNAME(MaxvalReal16)(*maxArray, __FILE__, __LINE__), -1.0);
+}
+#endif // LDBL_MANT_DIG == 113 || HAS_FLOAT128