[libc] Add hardware implementations of fma and fmaf for x86_64 and aarch64.

[sivachandra]

Summary

The current generic implementation of the fmaf function has been moved
to the FPUtil directory. This allows one use the fma operation from
implementations of other math functions like the trignometric functions
without depending on/requiring the fma/fmaf/fmal function targets. If
this pattern ends being convenient, we will switch all generic math
implementations to this pattern.

[sivachandra]

I have not yet fully tested this patch on x86_64. Will update soon if it requires any changes.

[sivachandra]

Add fma to the x86_64 entrypoint list.

[lntue]

libc/utils/FPUtil/FMA.h

15

Do we also need to check for _ _FMA__ flag?

[tschuett]

Stupid question: why assembler instead of builtins:
https://github.com/gcc-mirror/gcc/blob/master/gcc/testsuite/gcc.target/ia64/builtin-fma-1.c
I guess there is a tradeoff between control over the process and helping the register allocator.

[sivachandra]

In D100811#2700415, @tschuett wrote:

Stupid question: why assembler instead of builtins:
https://github.com/gcc-mirror/gcc/blob/master/gcc/testsuite/gcc.target/ia64/builtin-fma-1.c
I guess there is a tradeoff between control over the process and helping the register allocator.

Builtins can potentially silently call into the libc. I preferred inline assembly to avoid such a possibility. I think using target specific compiler intrinsics is a better option but they don't have convenient ways to convert from intrinsic types like __m128 to and from float/double.

libc/utils/FPUtil/FMA.h

15

I agree that not checking on __FMA__ is not ideal. But, I preferred not to condition on any other thing because it complicates the build system. For example, if __FMA__ was not defined, then we wouldn't get the hardware fma and fmaf but only get the integer fmaf implementation. But, other places will expect that fma and fmaf both are available. We can of course teach the build system to check for __FMA__ and propagate that info suitably to the rest of the libc build. That is a change we can take up separately if really required. I would prefer to see more examples of that pattern before building the appropriate machinery in to the build system.

libc/config/linux/aarch64/entrypoints.txt

set(TARGET_LIBM_ENTRYPOINTS
     libc.src.math.fabsf
     libc.src.math.fabsl
     libc.src.math.fdim
     libc.src.math.fdimf
     libc.src.math.fdiml
     libc.src.math.floor
     libc.src.math.floorf
     libc.src.math.floorl
+    libc.src.math.fma
     libc.src.math.fmaf
     libc.src.math.fmax
     libc.src.math.fmaxf
     libc.src.math.fmaxl
     libc.src.math.fmin
     libc.src.math.fminf
     libc.src.math.fminl
     libc.src.math.frexp

libc/config/linux/x86_64/entrypoints.txt

set(TARGET_LIBM_ENTRYPOINTS
     libc.src.math.fabsf
     libc.src.math.fabsl
     libc.src.math.fdim
     libc.src.math.fdimf
     libc.src.math.fdiml
     libc.src.math.floor
     libc.src.math.floorf
     libc.src.math.floorl
+    libc.src.math.fma
     libc.src.math.fmaf
     libc.src.math.fmin
     libc.src.math.fminf
     libc.src.math.fminl
     libc.src.math.fmax
     libc.src.math.fmaxf
     libc.src.math.fmaxl
     libc.src.math.frexp

libc/src/math/CMakeLists.txt

add_entrypoint_object(
         .generic.${name}
     )
     return()
   endif()
 
   message(FATAL_ERROR "No machine specific or generic implementation found for ${name}.")
 endfunction()
 
+add_entrypoint_object(
+  fmaf
+  SRCS
+    fmaf.cpp
+  HDRS
+    fmaf.h
+  DEPENDS
+    libc.utils.FPUtil.fputil
+  COMPILE_OPTIONS
+    -O2
+)
+
+add_entrypoint_object(
+  fma
+  SRCS
+    fma.cpp
+  HDRS
+    fma.h
+  DEPENDS
+    libc.utils.FPUtil.fputil
+  COMPILE_OPTIONS
+    -O2
+)
+
 add_math_entrypoint_object(ceil)
 add_math_entrypoint_object(ceilf)
 add_math_entrypoint_object(ceill)
 
 add_math_entrypoint_object(copysign)
 add_math_entrypoint_object(copysignf)
 add_math_entrypoint_object(copysignl)
 

 add_math_entrypoint_object(fdim)
 add_math_entrypoint_object(fdimf)
 add_math_entrypoint_object(fdiml)
 
 add_math_entrypoint_object(floor)
 add_math_entrypoint_object(floorf)
 add_math_entrypoint_object(floorl)
 
-add_math_entrypoint_object(fmaf)
-
 add_math_entrypoint_object(fmax)
 add_math_entrypoint_object(fmaxf)
 add_math_entrypoint_object(fmaxl)
 
 add_math_entrypoint_object(fmin)
 add_math_entrypoint_object(fminf)
 add_math_entrypoint_object(fminl)
 

libc/src/math/fma.h

+//===-- Implementation header for fma ---------------------------*- C++ -*-===//
+//
+// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
+// See https://llvm.org/LICENSE.txt for license information.
+// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
+//
+//===----------------------------------------------------------------------===//
+
+#ifndef LLVM_LIBC_SRC_MATH_FMA_H
+#define LLVM_LIBC_SRC_MATH_FMA_H
+
+namespace __llvm_libc {
+
+double fma(double x, double y, double z);
+
+} // namespace __llvm_libc
+
+#endif // LLVM_LIBC_SRC_MATH_FMA_H

libc/src/math/fma.cpp

+//===-- Implementation of fma function ------------------------------------===//
+//
+// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
+// See https://llvm.org/LICENSE.txt for license information.
+// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
+//
+//===----------------------------------------------------------------------===//
+
+#include "src/math/fma.h"
+#include "src/__support/common.h"
+
+#include "utils/FPUtil/FMA.h"
+
+namespace __llvm_libc {
+
+LLVM_LIBC_FUNCTION(double, fma, (double x, double y, double z)) {
+  return fputil::fma(x, y, z);
+}
+
+} // namespace __llvm_libc

libc/src/math/fmaf.cpp

+//===-- Implementation of fmaf function -----------------------------------===//
+//
+// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
+// See https://llvm.org/LICENSE.txt for license information.
+// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
+//
+//===----------------------------------------------------------------------===//
+
+#include "src/math/fmaf.h"
+#include "src/__support/common.h"
+
+#include "utils/FPUtil/FMA.h"
+
+namespace __llvm_libc {
+
+LLVM_LIBC_FUNCTION(float, fmaf, (float x, float y, float z)) {
+  return fputil::fma(x, y, z);
+}
+
+} // namespace __llvm_libc

libc/src/math/generic/CMakeLists.txt

add_entrypoint_object(
   SRCS
     sinf.cpp
   HDRS
     ../sinf.h
   DEPENDS
     .sincosf_utils
     libc.include.math
     libc.src.errno.__errno_location
+    libc.utils.FPUtil.fputil
+  COMPILE_OPTIONS
+    -O3
 )
 
 add_entrypoint_object(
   sincosf
   SRCS
     sincosf.cpp
   HDRS
     ../sincosf.h
SRCS
     nextafterl.cpp
   HDRS
     ../nextafterl.h
   DEPENDS
     libc.utils.FPUtil.fputil
   COMPILE_OPTIONS
     -O2
 )
-
-add_entrypoint_object(
-  fmaf
-  SRCS
-    fmaf.cpp
-  HDRS
-    ../fmaf.h
-  DEPENDS
-    libc.utils.FPUtil.fputil
-  COMPILE_OPTIONS
-    -O2
-)

libc/test/src/math/CMakeLists.txt

add_fp_unittest(
   SRCS
     fmaf_test.cpp
   DEPENDS
     libc.include.math
     libc.src.math.fmaf
     libc.utils.FPUtil.fputil
 )
 
+add_fp_unittest(
+  fma_test
+  NEED_MPFR
+  SUITE
+    libc_math_unittests
+  SRCS
+    fma_test.cpp
+  DEPENDS
+    libc.include.math
+    libc.src.math.fma
+    libc.utils.FPUtil.fputil
+)
+
 add_subdirectory(generic)
 add_subdirectory(exhaustive)

libc/test/src/math/fma_test.cpp

+//===-- Unittests for fma ------------------------------------------------===//
+//
+// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
+// See https://llvm.org/LICENSE.txt for license information.
+// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
+//
+//===----------------------------------------------------------------------===//
+
+#include "FmaTest.h"
+
+#include "src/math/fma.h"
+
+using LlvmLibcFmaTest = FmaTestTemplate<double>;
+
+TEST_F(LlvmLibcFmaTest, SpecialNumbers) {
+  testSpecialNumbers(&__llvm_libc::fma);
+}
+
+TEST_F(LlvmLibcFmaTest, SubnormalRange) {
+  testSubnormalRange(&__llvm_libc::fma);
+}
+
+TEST_F(LlvmLibcFmaTest, NormalRange) { testNormalRange(&__llvm_libc::fma); }

libc/utils/FPUtil/FMA.h

+//===-- Common header for FMA implementations -------------------*- C++ -*-===//
+//
+// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
+// See https://llvm.org/LICENSE.txt for license information.
+// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
+//
+//===----------------------------------------------------------------------===//
+
+#ifndef LLVM_LIBC_UTILS_FPUTIL_FMA_H
+#define LLVM_LIBC_UTILS_FPUTIL_FMA_H
+
+#include "utils/CPP/TypeTraits.h"
+
+#ifdef __x86_64__
+#include "x86_64/FMA.h"
+#elif defined(__aarch64__)
+#include "aarch64/FMA.h"
+#else
+#include "generic/FMA.h"
+
+namespace __llvm_libc {
+namespace fputil {
+
+// We have a generic implementation available only for single precision fma os
+// we restrict it to float values for now.
+template <typename T>
+static inline cpp::EnableIfType<cpp::IsSame<T, float>::Value, T> fma(T x, T y,
+                                                                     T z) {
+  return generic::fma(x, y, z);
+}
+
+} // namespace fputil
+} // namespace __llvm_libc
+
+#endif
+
+#endif // LLVM_LIBC_UTILS_FPUTIL_FMA_H

libc/utils/FPUtil/aarch64/FMA.h

+//===-- Aarch64 implementations of the fma function -------------*- C++ -*-===//
+//
+// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
+// See https://llvm.org/LICENSE.txt for license information.
+// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
+//
+//===----------------------------------------------------------------------===//
+
+#ifndef LLVM_LIBC_UTILS_FPUTIL_AARCH64_FMA_H
+#define LLVM_LIBC_UTILS_FPUTIL_AARCH64_FMA_H
+
+#include "utils/CPP/TypeTraits.h"
+
+namespace __llvm_libc {
+namespace fputil {
+
+template <typename T>
+cpp::EnableIfType<cpp::IsSame<T, float>::Value, T> fma(T x, T y, T z) {
+  float result;
+  __asm__ __volatile__("fmadd %s0, %s1, %s2, %s3\n\t"
+                       : "=w"(result)
+                       : "w"(x), "w"(y), "w"(z));
+  return result;
+}
+
+template <typename T>
+cpp::EnableIfType<cpp::IsSame<T, double>::Value, T> fma(T x, T y, T z) {
+  double result;
+  __asm__ __volatile__("fmadd %d0, %d1, %d2, %d3\n\t"
+                       : "=w"(result)
+                       : "w"(x), "w"(y), "w"(z));
+  return result;
+}
+
+} // namespace fputil
+} // namespace __llvm_libc
+
+#endif // LLVM_LIBC_UTILS_FPUTIL_AARCH64_FMA_H

libc/utils/FPUtil/generic/FMA.h

-//===-- Implementation of fmaf function -----------------------------------===//
+//===-- Common header for FMA implementations -------------------*- C++ -*-===//
 //
 // Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
 // See https://llvm.org/LICENSE.txt for license information.
 // SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
 //
 //===----------------------------------------------------------------------===//
 
-#include "src/math/fmaf.h"
-#include "src/__support/common.h"
+#ifndef LLVM_LIBC_UTILS_FPUTIL_GENERIC_FMA_H
+#define LLVM_LIBC_UTILS_FPUTIL_GENERIC_FMA_H
 
-#include "utils/FPUtil/FEnv.h"
-#include "utils/FPUtil/FPBits.h"
+#include "utils/CPP/TypeTraits.h"
 
 namespace __llvm_libc {
+namespace fputil {
+namespace generic {
 
-LLVM_LIBC_FUNCTION(float, fmaf, (float x, float y, float z)) {
+template <typename T>
+static inline cpp::EnableIfType<cpp::IsSame<T, float>::Value, T> fma(T x, T y,
+                                                                     T z) {
   // Product is exact.
   double prod = static_cast<double>(x) * static_cast<double>(y);
   double z_d = static_cast<double>(z);
   double sum = prod + z_d;
   fputil::FPBits<double> bit_prod(prod), bitz(z_d), bit_sum(sum);
 
   if (!(bit_sum.isInfOrNaN() || bit_sum.isZero())) {
     // Since the sum is computed in double precision, rounding might happen
if (!t.isZero() && ((bit_sum.encoding.mantissa & 0xfff'ffffULL) == 0)) {
         --bit_sum.encoding.mantissa;
       }
     }
   }
 
   return static_cast<float>(static_cast<double>(bit_sum));
 }
 
+} // namespace generic
+} // namespace fputil
 } // namespace __llvm_libc
+
+#endif // Generic fma implementations
+
+#endif // LLVM_LIBC_UTILS_FPUTIL_GENERIC_FMA_H

libc/utils/FPUtil/generic/README.md

+This directory contains machine independent implementations of floating point
+operations. The implementations are nested in the namespace
+`__llvm_libc::fputil::generic`. This is to facilitate calling these generic
+implementations from machine dependent implementations. Consider the example of
+the fuse-multiply-add operation (FMA). The C standard library requires three
+different flavors, `fma` which operates double precsion numbers, `fmaf` which
+operates on single precision numbers, and `fmal` which operates on `lond double`
+numbers. On Aarch64, there are hardware instructions which implement the single
+and double precision flavors but not the `long double` flavor. For such targets,
+we want to be able to call the generic `long double` implementation from the
+`long double` flavor. By putting the generic implementations in a separate
+nested namespace, we will be to call them as follows:
+
+```
+namespace __llvm_libc {
+namespace fputil {
+
+long double fmal(long double x, long double y, long double z) {
+  return generic::fmal(x, y, z);
+}
+
+} // namespace fputil
+} // namespace __llvm_libc
+```
+
+Note that actual code might not be as straightforward as above (for example,
+we might want to prevent implicit type promotions by using some template
+facilities). But, the general idea is very similar.

libc/utils/FPUtil/x86_64/FMA.h

+//===-- x86_64 implementations of the fma function --------------*- C++ -*-===//
+//
+// Part of the LLVM Project, under the Apache License v2.0 with LLVM Exceptions.
+// See https://llvm.org/LICENSE.txt for license information.
+// SPDX-License-Identifier: Apache-2.0 WITH LLVM-exception
+//
+//===----------------------------------------------------------------------===//
+
+#ifndef LLVM_LIBC_UTILS_FPUTIL_X86_64_FMA_H
+#define LLVM_LIBC_UTILS_FPUTIL_X86_64_FMA_H
+
+#include "utils/CPP/TypeTraits.h"
+
+namespace __llvm_libc {
+namespace fputil {
+
+template <typename T>
+static inline cpp::EnableIfType<cpp::IsSame<T, float>::Value, T> fma(T x, T y,
+                                                                     T z) {
+  float result = x;
+  __asm__ __volatile__("vfmadd213ss %x2, %x1, %x0"
+                       : "+x"(result)
+                       : "x"(y), "x"(z));
+  return result;
+}
+
+template <typename T>
+static inline cpp::EnableIfType<cpp::IsSame<T, double>::Value, T> fma(T x, T y,
+                                                                      T z) {
+  double result = x;
+  __asm__ __volatile__("vfmadd213sd %x2, %x1, %x0"
+                       : "+x"(result)
+                       : "x"(y), "x"(z));
+  return result;
+}
+
+} // namespace fputil
+} // namespace __llvm_libc
+
+#endif // LLVM_LIBC_UTILS_FPUTIL_X86_64_FMA_H