[libc-commits] [libc] [libc][math][c23] fmul correcly rounded to all rounding modes (PR #91537)

Job Henandez Lara via libc-commits libc-commits at lists.llvm.org
Sun Jun 2 17:17:12 PDT 2024


https://github.com/Jobhdez updated https://github.com/llvm/llvm-project/pull/91537

>From 9f16ac2f1d4b31aa0af88300dd70738fa8b95b44 Mon Sep 17 00:00:00 2001
From: Job Hernandez <hj93 at protonmail.com>
Date: Wed, 8 May 2024 14:27:15 -0700
Subject: [PATCH 1/8] start fmul template

---
 libc/config/linux/aarch64/entrypoints.txt   |  3 ++
 libc/config/linux/arm/entrypoints.txt       |  3 ++
 libc/config/linux/riscv/entrypoints.txt     |  3 ++
 libc/config/linux/x86_64/entrypoints.txt    |  3 ++
 libc/config/windows/entrypoints.txt         |  3 ++
 libc/spec/stdc.td                           |  5 +++
 libc/src/__support/FPUtil/BasicOperations.h | 12 +++++++
 libc/src/math/CMakeLists.txt                |  4 +++
 libc/src/math/dmull.h                       | 18 ++++++++++
 libc/src/math/fmul.h                        | 18 ++++++++++
 libc/src/math/fmull.h                       | 18 ++++++++++
 libc/src/math/generic/CMakeLists.txt        | 36 +++++++++++++++++++
 libc/src/math/generic/dmull.cpp             | 19 ++++++++++
 libc/src/math/generic/fmul.cpp              | 19 ++++++++++
 libc/src/math/generic/fmull.cpp             | 19 ++++++++++
 libc/test/src/math/smoke/CMakeLists.txt     | 40 ++++++++++++++++++++-
 libc/test/src/math/smoke/FMulTest.h         | 32 +++++++++++++++++
 libc/test/src/math/smoke/dmull_test.cpp     | 13 +++++++
 libc/test/src/math/smoke/fmul_test.cpp      | 13 +++++++
 libc/test/src/math/smoke/fmull_test.cpp     | 13 +++++++
 20 files changed, 293 insertions(+), 1 deletion(-)
 create mode 100644 libc/src/math/dmull.h
 create mode 100644 libc/src/math/fmul.h
 create mode 100644 libc/src/math/fmull.h
 create mode 100644 libc/src/math/generic/dmull.cpp
 create mode 100644 libc/src/math/generic/fmul.cpp
 create mode 100644 libc/src/math/generic/fmull.cpp
 create mode 100644 libc/test/src/math/smoke/FMulTest.h
 create mode 100644 libc/test/src/math/smoke/dmull_test.cpp
 create mode 100644 libc/test/src/math/smoke/fmul_test.cpp
 create mode 100644 libc/test/src/math/smoke/fmull_test.cpp

diff --git a/libc/config/linux/aarch64/entrypoints.txt b/libc/config/linux/aarch64/entrypoints.txt
index ad50b6f59cdcc..c92d238f14f7d 100644
--- a/libc/config/linux/aarch64/entrypoints.txt
+++ b/libc/config/linux/aarch64/entrypoints.txt
@@ -394,6 +394,9 @@ set(TARGET_LIBM_ENTRYPOINTS
     libc.src.math.fminimum_mag_num
     libc.src.math.fminimum_mag_numf
     libc.src.math.fminimum_mag_numl
+    libc.src.math.fmul
+    libc.src.math.fmull
+    libc.src.math.dmull
     libc.src.math.fmod
     libc.src.math.fmodf
     libc.src.math.fmodl
diff --git a/libc/config/linux/arm/entrypoints.txt b/libc/config/linux/arm/entrypoints.txt
index 335981ff7dc7c..c3bb9bf666764 100644
--- a/libc/config/linux/arm/entrypoints.txt
+++ b/libc/config/linux/arm/entrypoints.txt
@@ -261,6 +261,9 @@ set(TARGET_LIBM_ENTRYPOINTS
     libc.src.math.fminimum_mag_num
     libc.src.math.fminimum_mag_numf
     libc.src.math.fminimum_mag_numl
+    libc.src.math.fmul
+    libc.src.math.fmull
+    libc.src.math.dmull
     libc.src.math.fmod
     libc.src.math.fmodf
     libc.src.math.frexp
diff --git a/libc/config/linux/riscv/entrypoints.txt b/libc/config/linux/riscv/entrypoints.txt
index 479af40b5b26b..8b4744445ae01 100644
--- a/libc/config/linux/riscv/entrypoints.txt
+++ b/libc/config/linux/riscv/entrypoints.txt
@@ -402,6 +402,9 @@ set(TARGET_LIBM_ENTRYPOINTS
     libc.src.math.fminimum_mag_num
     libc.src.math.fminimum_mag_numf
     libc.src.math.fminimum_mag_numl
+    libc.src.math.fmul
+    libc.src.math.fmull
+    libc.src.math.dmull
     libc.src.math.fmod
     libc.src.math.fmodf
     libc.src.math.fmodl
diff --git a/libc/config/linux/x86_64/entrypoints.txt b/libc/config/linux/x86_64/entrypoints.txt
index 5e3ddd34fb4dc..9623f0b967c40 100644
--- a/libc/config/linux/x86_64/entrypoints.txt
+++ b/libc/config/linux/x86_64/entrypoints.txt
@@ -421,6 +421,9 @@ set(TARGET_LIBM_ENTRYPOINTS
     libc.src.math.fminimum_mag_num
     libc.src.math.fminimum_mag_numf
     libc.src.math.fminimum_mag_numl
+    libc.src.math.fmul
+    libc.src.math.fmull
+    libc.src.math.dmull
     libc.src.math.fmod
     libc.src.math.fmodf
     libc.src.math.fmodl
diff --git a/libc/config/windows/entrypoints.txt b/libc/config/windows/entrypoints.txt
index 71216530c4041..199ce9bc91cd9 100644
--- a/libc/config/windows/entrypoints.txt
+++ b/libc/config/windows/entrypoints.txt
@@ -180,6 +180,9 @@ set(TARGET_LIBM_ENTRYPOINTS
     libc.src.math.fminimum_mag_num
     libc.src.math.fminimum_mag_numf
     libc.src.math.fminimum_mag_numl
+    libc.src.math.fmul
+    libc.src.math.fmull
+    libc.src.math.dmull
     libc.src.math.fmod
     libc.src.math.fmodf
     libc.src.math.fmodl
diff --git a/libc/spec/stdc.td b/libc/spec/stdc.td
index eb67c9b0b009b..0130ace283a55 100644
--- a/libc/spec/stdc.td
+++ b/libc/spec/stdc.td
@@ -457,6 +457,11 @@ def StdC : StandardSpec<"stdc"> {
           FunctionSpec<"fminimum_mag_numl", RetValSpec<LongDoubleType>, [ArgSpec<LongDoubleType>, ArgSpec<LongDoubleType>]>,
           GuardedFunctionSpec<"fminimum_mag_numf128", RetValSpec<Float128Type>, [ArgSpec<Float128Type>, ArgSpec<Float128Type>], "LIBC_TYPES_HAS_FLOAT128">,
 
+          FunctionSpec<"fmul", RetValSpec<FloatType>, [ArgSpec<DoubleType>, ArgSpec<DoubleType>]>,
+	  FunctionSpec<"fmull", RetValSpec<FloatType>, [ArgSpec<LongDoubleType>, ArgSpec<LongDoubleType>]>,
+	  FunctionSpec<"dmull", RetValSpec<DoubleType>, [ArgSpec<LongDoubleType>, ArgSpec<LongDoubleType>]>,
+	    
+
           FunctionSpec<"fma", RetValSpec<DoubleType>, [ArgSpec<DoubleType>, ArgSpec<DoubleType>, ArgSpec<DoubleType>]>,
           FunctionSpec<"fmaf", RetValSpec<FloatType>, [ArgSpec<FloatType>, ArgSpec<FloatType>, ArgSpec<FloatType>]>,
 
diff --git a/libc/src/__support/FPUtil/BasicOperations.h b/libc/src/__support/FPUtil/BasicOperations.h
index e5ac101fedc0e..e75b4b9391bb8 100644
--- a/libc/src/__support/FPUtil/BasicOperations.h
+++ b/libc/src/__support/FPUtil/BasicOperations.h
@@ -176,6 +176,18 @@ LIBC_INLINE T fdim(T x, T y) {
 
   return (x > y ? x - y : 0);
 }
+  
+LIBC_INLINE float fmul(double x, double y) {
+  return static_cast<float>(x * y);
+}
+
+LIBC_INLINE float fmull(long double x, long double y) {
+  return static_cast<float>(x * y);
+}
+
+LIBC_INLINE double dmull(long double x, long double y) {
+  return static_cast<double>(x * y);
+}
 
 template <typename T, cpp::enable_if_t<cpp::is_floating_point_v<T>, int> = 0>
 LIBC_INLINE int canonicalize(T &cx, const T &x) {
diff --git a/libc/src/math/CMakeLists.txt b/libc/src/math/CMakeLists.txt
index c34c58575441d..dd26109d9db4d 100644
--- a/libc/src/math/CMakeLists.txt
+++ b/libc/src/math/CMakeLists.txt
@@ -164,6 +164,10 @@ add_math_entrypoint_object(fminimum_mag_numf)
 add_math_entrypoint_object(fminimum_mag_numl)
 add_math_entrypoint_object(fminimum_mag_numf128)
 
+add_math_entrypoint_object(fmul)
+add_math_entrypoint_object(fmull)
+add_math_entrypoint_object(dmull)
+
 add_math_entrypoint_object(fmod)
 add_math_entrypoint_object(fmodf)
 add_math_entrypoint_object(fmodl)
diff --git a/libc/src/math/dmull.h b/libc/src/math/dmull.h
new file mode 100644
index 0000000000000..556eb1559d535
--- /dev/null
+++ b/libc/src/math/dmull.h
@@ -0,0 +1,18 @@
+//===-- Implementation header for dmull -------------------------*- 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_DMULL_H
+#define LLVM_LIBC_SRC_MATH_DMULL_H
+
+namespace LIBC_NAMESPACE {
+
+double dmull(long double x, long double y);
+
+} // namespace LIBC_NAMESPACE
+
+#endif // LLVM_LIBC_SRC_MATH_DMULL_H
diff --git a/libc/src/math/fmul.h b/libc/src/math/fmul.h
new file mode 100644
index 0000000000000..fbc1069db733e
--- /dev/null
+++ b/libc/src/math/fmul.h
@@ -0,0 +1,18 @@
+//===-- Implementation header for fmul --------------------------*- 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_FMUL_H
+#define LLVM_LIBC_SRC_MATH_FMUL_H
+
+namespace LIBC_NAMESPACE {
+
+float fmul(double x, double y);
+
+} // namespace LIBC_NAMESPACE
+
+#endif // LLVM_LIBC_SRC_MATH_FMUL_H
diff --git a/libc/src/math/fmull.h b/libc/src/math/fmull.h
new file mode 100644
index 0000000000000..aebf52a9b0c0a
--- /dev/null
+++ b/libc/src/math/fmull.h
@@ -0,0 +1,18 @@
+//===-- Implementation header for fmull -------------------------*- 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_FMULL_H
+#define LLVM_LIBC_SRC_MATH_FMULL_H
+
+namespace LIBC_NAMESPACE {
+
+float fmull(long double x, long double y);
+
+} // namespace LIBC_NAMESPACE
+
+#endif // LLVM_LIBC_SRC_MATH_FMULL_H
diff --git a/libc/src/math/generic/CMakeLists.txt b/libc/src/math/generic/CMakeLists.txt
index daaf505008ca1..3554ec088c49f 100644
--- a/libc/src/math/generic/CMakeLists.txt
+++ b/libc/src/math/generic/CMakeLists.txt
@@ -2042,6 +2042,42 @@ add_entrypoint_object(
     -O3
 )
 
+add_entrypoint_object(
+  fmul
+  SRCS
+    fmul.cpp
+  HDRS
+    ../fmul.h
+  DEPENDS
+    libc.src.__support.FPUtil.basic_operations
+  COMPILE_OPTIONS
+    -O2
+)
+
+add_entrypoint_object(
+  fmull
+  SRCS
+    fmull.cpp
+  HDRS
+    ../fmull.h
+  DEPENDS
+    libc.src.__support.FPUtil.basic_operations
+  COMPILE_OPTIONS
+    -O2
+)
+
+add_entrypoint_object(
+  dmull
+  SRCS
+    dmull.cpp
+  HDRS
+    ../dmull.h
+  DEPENDS
+    libc.src.__support.FPUtil.basic_operations
+  COMPILE_OPTIONS
+    -O2
+)
+
 add_entrypoint_object(
   sqrt
   SRCS
diff --git a/libc/src/math/generic/dmull.cpp b/libc/src/math/generic/dmull.cpp
new file mode 100644
index 0000000000000..2c70841dc2d28
--- /dev/null
+++ b/libc/src/math/generic/dmull.cpp
@@ -0,0 +1,19 @@
+//===-- Implementation of dmull 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/dmull.h"
+#include "src/__support/FPUtil/BasicOperations.h"
+#include "src/__support/common.h"
+
+namespace LIBC_NAMESPACE {
+
+LLVM_LIBC_FUNCTION(double, dmull, (long double x, long double y)) {
+  return fputil::dmull(x, y);
+}
+
+} // namespace LIBC_NAMESPACE
diff --git a/libc/src/math/generic/fmul.cpp b/libc/src/math/generic/fmul.cpp
new file mode 100644
index 0000000000000..77452013951d7
--- /dev/null
+++ b/libc/src/math/generic/fmul.cpp
@@ -0,0 +1,19 @@
+//===-- Implementation of fmul 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/fmul.h"
+#include "src/__support/FPUtil/BasicOperations.h"
+#include "src/__support/common.h"
+
+namespace LIBC_NAMESPACE {
+
+LLVM_LIBC_FUNCTION(float, fmul, (double x, double y)) {
+  return fputil::fmul(x, y);
+}
+
+} // namespace LIBC_NAMESPACE
diff --git a/libc/src/math/generic/fmull.cpp b/libc/src/math/generic/fmull.cpp
new file mode 100644
index 0000000000000..f0d0e0d8da22a
--- /dev/null
+++ b/libc/src/math/generic/fmull.cpp
@@ -0,0 +1,19 @@
+//===-- Implementation of fmull 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/fmull.h"
+#include "src/__support/FPUtil/BasicOperations.h"
+#include "src/__support/common.h"
+
+namespace LIBC_NAMESPACE {
+
+LLVM_LIBC_FUNCTION(float, fmull, (long double x, long double y)) {
+  return fputil::fmull(x, y);
+}
+
+} // namespace LIBC_NAMESPACE
diff --git a/libc/test/src/math/smoke/CMakeLists.txt b/libc/test/src/math/smoke/CMakeLists.txt
index 112b2985829ca..30e6a3cb98184 100644
--- a/libc/test/src/math/smoke/CMakeLists.txt
+++ b/libc/test/src/math/smoke/CMakeLists.txt
@@ -1911,7 +1911,6 @@ add_fp_unittest(
     libc.src.__support.FPUtil.fp_bits
 )
 
-
 add_fp_unittest(
   fminimum_mag_numf_test
   SUITE
@@ -1964,6 +1963,45 @@ add_fp_unittest(
     libc.src.__support.FPUtil.fp_bits
 )
 
+add_fp_unittest(
+  fmul_test
+  SUITE
+    libc-math-smoke-tests
+  SRCS
+    fmul_test.cpp
+  HDRS
+    FMulTest.h
+  DEPENDS
+    libc.src.math.fmul
+    libc.src.__support.FPUtil.fp_bits
+)
+
+add_fp_unittest(
+  fmull_test
+  SUITE
+    libc-math-smoke-tests
+  SRCS
+    fmull_test.cpp
+  HDRS
+    FMulTest.h
+  DEPENDS
+    libc.src.math.fmull
+    libc.src.__support.FPUtil.fp_bits
+)
+
+add_fp_unittest(
+  dmull_test
+  SUITE
+    libc-math-smoke-tests
+  SRCS
+    dmull_test.cpp
+  HDRS
+    FMulTest.h
+  DEPENDS
+    libc.src.math.dmull
+    libc.src.__support.FPUtil.fp_bits
+)
+
 add_fp_unittest(
   sqrtf_test
   SUITE
diff --git a/libc/test/src/math/smoke/FMulTest.h b/libc/test/src/math/smoke/FMulTest.h
new file mode 100644
index 0000000000000..d217dc369ef81
--- /dev/null
+++ b/libc/test/src/math/smoke/FMulTest.h
@@ -0,0 +1,32 @@
+//===-- Utility class to test fmul[f|l] ---------------------*- 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_TEST_SRC_MATH_SMOKE_FMULTEST_H
+#define LLVM_LIBC_TEST_SRC_MATH_SMOKE_FMULTEST_H
+
+#include "test/UnitTest/FEnvSafeTest.h"
+#include "test/UnitTest/FPMatcher.h"
+#include "test/UnitTest/Test.h"
+
+template <typename T, typename R>
+class FmulTest : public LIBC_NAMESPACE::testing::FEnvSafeTest {
+
+public:
+  typedef T (*FMulFunc)(R, R);
+
+  void testMul(FMulFunc func) {
+    EXPECT_FP_EQ(T(1.0), func(1.0, 1.0));
+  }
+  
+};
+
+#define LIST_FMUL_TESTS(T, R, func)				       \
+  using LlvmLibcFmulTest = FmulTest<T,R>;                                \
+  TEST_F(LlvmLibcFmulTest, Mul) { testMul(&func); }                
+
+#endif // LLVM_LIBC_TEST_SRC_MATH_SMOKE_FMULTEST_H
diff --git a/libc/test/src/math/smoke/dmull_test.cpp b/libc/test/src/math/smoke/dmull_test.cpp
new file mode 100644
index 0000000000000..f094b2578158f
--- /dev/null
+++ b/libc/test/src/math/smoke/dmull_test.cpp
@@ -0,0 +1,13 @@
+//===-- Unittests for dmull------------------------------------------------===//
+//
+// 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 "FMulTest.h"
+
+#include "src/math/dmull.h"
+
+LIST_FMUL_TESTS(double, long double, LIBC_NAMESPACE::dmull)
diff --git a/libc/test/src/math/smoke/fmul_test.cpp b/libc/test/src/math/smoke/fmul_test.cpp
new file mode 100644
index 0000000000000..0eb664f7411ee
--- /dev/null
+++ b/libc/test/src/math/smoke/fmul_test.cpp
@@ -0,0 +1,13 @@
+//===-- Unittests for fmul-------------------------------------------------===//
+//
+// 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 "FMulTest.h"
+
+#include "src/math/fmul.h"
+
+LIST_FMUL_TESTS(float, double, LIBC_NAMESPACE::fmul)
diff --git a/libc/test/src/math/smoke/fmull_test.cpp b/libc/test/src/math/smoke/fmull_test.cpp
new file mode 100644
index 0000000000000..d1151e298c0b6
--- /dev/null
+++ b/libc/test/src/math/smoke/fmull_test.cpp
@@ -0,0 +1,13 @@
+//===-- Unittests for fmull------------------------------------------------===//
+//
+// 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 "FMulTest.h"
+
+#include "src/math/fmull.h"
+
+LIST_FMUL_TESTS(float, long double, LIBC_NAMESPACE::fmull)

>From 92d66d8a93d33599ff8ba286f6860838899f494b Mon Sep 17 00:00:00 2001
From: Job Hernandez <hj93 at protonmail.com>
Date: Wed, 8 May 2024 14:27:55 -0700
Subject: [PATCH 2/8] format code

---
 libc/src/__support/FPUtil/BasicOperations.h |  6 ++----
 libc/test/src/math/smoke/FMulTest.h         | 11 ++++-------
 2 files changed, 6 insertions(+), 11 deletions(-)

diff --git a/libc/src/__support/FPUtil/BasicOperations.h b/libc/src/__support/FPUtil/BasicOperations.h
index e75b4b9391bb8..e3581910e9799 100644
--- a/libc/src/__support/FPUtil/BasicOperations.h
+++ b/libc/src/__support/FPUtil/BasicOperations.h
@@ -176,10 +176,8 @@ LIBC_INLINE T fdim(T x, T y) {
 
   return (x > y ? x - y : 0);
 }
-  
-LIBC_INLINE float fmul(double x, double y) {
-  return static_cast<float>(x * y);
-}
+
+LIBC_INLINE float fmul(double x, double y) { return static_cast<float>(x * y); }
 
 LIBC_INLINE float fmull(long double x, long double y) {
   return static_cast<float>(x * y);
diff --git a/libc/test/src/math/smoke/FMulTest.h b/libc/test/src/math/smoke/FMulTest.h
index d217dc369ef81..412d9eedf1f7c 100644
--- a/libc/test/src/math/smoke/FMulTest.h
+++ b/libc/test/src/math/smoke/FMulTest.h
@@ -19,14 +19,11 @@ class FmulTest : public LIBC_NAMESPACE::testing::FEnvSafeTest {
 public:
   typedef T (*FMulFunc)(R, R);
 
-  void testMul(FMulFunc func) {
-    EXPECT_FP_EQ(T(1.0), func(1.0, 1.0));
-  }
-  
+  void testMul(FMulFunc func) { EXPECT_FP_EQ(T(1.0), func(1.0, 1.0)); }
 };
 
-#define LIST_FMUL_TESTS(T, R, func)				       \
-  using LlvmLibcFmulTest = FmulTest<T,R>;                                \
-  TEST_F(LlvmLibcFmulTest, Mul) { testMul(&func); }                
+#define LIST_FMUL_TESTS(T, R, func)                                            \
+  using LlvmLibcFmulTest = FmulTest<T, R>;                                     \
+  TEST_F(LlvmLibcFmulTest, Mul) { testMul(&func); }
 
 #endif // LLVM_LIBC_TEST_SRC_MATH_SMOKE_FMULTEST_H

>From f00264caf86d1ba451f669d2c51c81a681df8f3e Mon Sep 17 00:00:00 2001
From: Job Hernandez <hj93 at protonmail.com>
Date: Sat, 1 Jun 2024 09:19:25 -0700
Subject: [PATCH 3/8] make actual llvm impl.

---
 libc/src/__support/FPUtil/BasicOperations.h | 154 +++++++++++++++++++-
 libc/test/src/math/smoke/FMulTest.h         |   6 +-
 2 files changed, 158 insertions(+), 2 deletions(-)

diff --git a/libc/src/__support/FPUtil/BasicOperations.h b/libc/src/__support/FPUtil/BasicOperations.h
index e3581910e9799..58dc5c7276edc 100644
--- a/libc/src/__support/FPUtil/BasicOperations.h
+++ b/libc/src/__support/FPUtil/BasicOperations.h
@@ -13,6 +13,7 @@
 #include "FPBits.h"
 
 #include "FEnvImpl.h"
+#include "src/__support/CPP/bit.h"
 #include "src/__support/CPP/type_traits.h"
 #include "src/__support/common.h"
 #include "src/__support/macros/optimization.h" // LIBC_UNLIKELY
@@ -177,8 +178,159 @@ LIBC_INLINE T fdim(T x, T y) {
   return (x > y ? x - y : 0);
 }
 
-LIBC_INLINE float fmul(double x, double y) { return static_cast<float>(x * y); }
+// helpers for fmul
 
+uint64_t maxu(uint64_t A, uint64_t B) {
+  return A > B ? A : B;
+}
+
+uint64_t mul(uint64_t a, uint64_t b) {
+    __uint128_t product = static_cast<__uint128_t>(a) * static_cast<__uint128_t>(b);
+    return static_cast<uint64_t>(product >> 64);
+}
+
+
+uint64_t mullow(uint64_t a, uint64_t b) {
+    __uint128_t product = static_cast<__uint128_t>(a) * static_cast<__uint128_t>(b);
+    return static_cast<uint64_t>(product);
+}
+
+
+uint64_t nlz(uint64_t x) {
+  uint64_t z = 0;
+
+  if (x == 0) return 64;
+  if (x <= 0x00000000FFFFFFFF) {
+    z = z + 32;
+    x = x << 32;
+  }
+  if (x <= 0x0000FFFFFFFFFFFF) {
+    z = z + 16;
+    x = x << 16;
+  }
+  if (x <= 0x00FFFFFFFFFFFFFF) {
+    z = z + 8;
+    x = x << 8;
+  }
+  if (x <= 0x0FFFFFFFFFFFFFFF) {
+    z = z + 4;
+    x = x << 4;
+  }
+  if (x <= 0x3FFFFFFFFFFFFFFF) {
+    z = z + 2;
+    x = x << 2;
+  }
+  if (x <= 0x7FFFFFFFFFFFFFFF) {
+    z = z + 1;
+  }
+  return z;
+}
+  
+LIBC_INLINE float fmul(double x, double y) {
+  
+  
+  auto x_bits = FPBits<double>(x);
+  uint64_t x_u = x_bits.uintval();
+
+  auto y_bits = FPBits<double>(y);
+  uint64_t y_u = y_bits.uintval();
+
+  uint64_t absx = x_u & 0x7FFFFFFFFFFFFFFF;
+  uint64_t absy = y_u & 0x7FFFFFFFFFFFFFFF;
+
+  uint64_t exponent_x = absx >> 52;
+  uint64_t exponent_y = absy >> 52;
+
+  //uint64_t x_normal_bit = absx >= 0x10000000000000;
+  //uint64_t y_normal_bit = absy >= 0x10000000000000;
+
+  uint64_t mx, my;
+  
+  mx = maxu(nlz(absx), 11); 
+
+  //lambdax = mx - 11;  
+
+  my = maxu(nlz(absy), 11);
+
+  //lambday = my - 11;
+
+  int32_t dm1;
+  uint64_t mpx, mpy, highs,lows, b;
+  uint64_t  g, hight, lowt, morlowt, c, m;
+  mpx = (x_u << mx) | 0x8000000000000000;
+  mpy = (y_u << my) |  0x8000000000000000;
+  highs = mul(mpx, mpy); 
+  c = highs >= 0x8000000000000000; 
+  lows = mullow(mpx, mpy);
+
+  lowt = (lows != 0);
+
+  g = (highs >> (38 + c)) & 1;
+  hight = (highs << (55 - c)) != 0;
+
+  int32_t exint = static_cast<uint32_t>(exponent_x);
+  int32_t eyint = static_cast<uint32_t>(exponent_y);
+  int32_t cint = static_cast<uint32_t>(c);
+  dm1 = ((exint + eyint) - 1919) + cint;
+
+  if (dm1 >= 255) {
+    dm1 = 255;
+    
+    m = 0;
+  }
+  else if (dm1<=0) {
+    m = static_cast<uint32_t>((highs >> (39 + c)) >> (1 - dm1));
+    dm1 = 0;
+    morlowt = m | lowt;
+    b = g & (morlowt | hight);
+  }
+  else  {
+    m = static_cast<uint32_t>(highs >> (39 + c));
+    morlowt = m | lowt;
+    b = g & (morlowt | hight);
+  }
+
+  uint32_t sr = static_cast<uint32_t>((x_u ^ y_u) & 0x8000000000000000);
+  
+
+  uint32_t exp16 = sr | (static_cast<uint32_t>(dm1) << 23);
+  if (dm1 == 0) {
+    uint32_t m2 = static_cast<uint32_t>(m);
+    uint32_t result = (static_cast<uint32_t>(exp16) + m2) + static_cast<uint32_t>(b);
+   
+    
+    float result32 = cpp::bit_cast<float>(result);
+    
+    return result32;
+      
+ } else {
+   constexpr uint32_t FLOAT32_MANTISSA_MASK = 0b00000000011111111111111111111111;
+   uint32_t m2 = static_cast<uint32_t>(m) & FLOAT32_MANTISSA_MASK;
+   
+   uint32_t result = (static_cast<uint32_t>(exp16) + m2) + static_cast<uint32_t>(b);
+  
+  float result16 = cpp::bit_cast<float>(result);
+  // std::memcpy(&result16, &result, sizeof(result16));  
+  //result = *reinterpret_cast<_Float16*>(&result);
+  
+  return result16;
+ }
+}
+  /*
+LIBC_INLINE float fmul(double x, double y){
+  FPBits<double> bitx(x), bity(y);
+
+  long long int product;
+  long long int p;
+  
+  if (bitx.is_normal() && bity.is_normal()) {
+    product = bitx.get_mantissa() * bity.get_mantissa();
+    p = product & ((1ULL << 24) -1);
+    sr = bitx.sign() ^ bity.sign();
+    return pow(-1, sr) * pow(p, p.biased_exponent());
+  }
+}
+  */
 LIBC_INLINE float fmull(long double x, long double y) {
   return static_cast<float>(x * y);
 }
diff --git a/libc/test/src/math/smoke/FMulTest.h b/libc/test/src/math/smoke/FMulTest.h
index 412d9eedf1f7c..afaca64d866c1 100644
--- a/libc/test/src/math/smoke/FMulTest.h
+++ b/libc/test/src/math/smoke/FMulTest.h
@@ -19,7 +19,11 @@ class FmulTest : public LIBC_NAMESPACE::testing::FEnvSafeTest {
 public:
   typedef T (*FMulFunc)(R, R);
 
-  void testMul(FMulFunc func) { EXPECT_FP_EQ(T(1.0), func(1.0, 1.0)); }
+  void testMul(FMulFunc func) {
+    EXPECT_FP_EQ(T(15.0), func(3.0, 5.0));
+    EXPECT_FP_EQ(T(7.34684e-40), func(0x1.0p1,0x1.0p-131));
+    EXPECT_FP_EQ(T(5.87747e-39), func(0x1.0p2, 0x1.0p-129));
+}
 };
 
 #define LIST_FMUL_TESTS(T, R, func)                                            \

>From beaea4eca07ab2ffb4ea19fbdbdb39a35072dec4 Mon Sep 17 00:00:00 2001
From: Job Hernandez <hj93 at protonmail.com>
Date: Sat, 1 Jun 2024 09:20:02 -0700
Subject: [PATCH 4/8] format code

---
 libc/src/__support/FPUtil/BasicOperations.h | 117 ++++++++++----------
 libc/test/src/math/smoke/FMulTest.h         |   4 +-
 2 files changed, 59 insertions(+), 62 deletions(-)

diff --git a/libc/src/__support/FPUtil/BasicOperations.h b/libc/src/__support/FPUtil/BasicOperations.h
index 58dc5c7276edc..1926659e5f78f 100644
--- a/libc/src/__support/FPUtil/BasicOperations.h
+++ b/libc/src/__support/FPUtil/BasicOperations.h
@@ -180,26 +180,25 @@ LIBC_INLINE T fdim(T x, T y) {
 
 // helpers for fmul
 
-uint64_t maxu(uint64_t A, uint64_t B) {
-  return A > B ? A : B;
-}
+uint64_t maxu(uint64_t A, uint64_t B) { return A > B ? A : B; }
 
 uint64_t mul(uint64_t a, uint64_t b) {
-    __uint128_t product = static_cast<__uint128_t>(a) * static_cast<__uint128_t>(b);
-    return static_cast<uint64_t>(product >> 64);
+  __uint128_t product =
+      static_cast<__uint128_t>(a) * static_cast<__uint128_t>(b);
+  return static_cast<uint64_t>(product >> 64);
 }
 
-
 uint64_t mullow(uint64_t a, uint64_t b) {
-    __uint128_t product = static_cast<__uint128_t>(a) * static_cast<__uint128_t>(b);
-    return static_cast<uint64_t>(product);
+  __uint128_t product =
+      static_cast<__uint128_t>(a) * static_cast<__uint128_t>(b);
+  return static_cast<uint64_t>(product);
 }
 
-
 uint64_t nlz(uint64_t x) {
   uint64_t z = 0;
 
-  if (x == 0) return 64;
+  if (x == 0)
+    return 64;
   if (x <= 0x00000000FFFFFFFF) {
     z = z + 32;
     x = x << 32;
@@ -225,10 +224,9 @@ uint64_t nlz(uint64_t x) {
   }
   return z;
 }
-  
+
 LIBC_INLINE float fmul(double x, double y) {
-  
-  
+
   auto x_bits = FPBits<double>(x);
   uint64_t x_u = x_bits.uintval();
 
@@ -241,26 +239,26 @@ LIBC_INLINE float fmul(double x, double y) {
   uint64_t exponent_x = absx >> 52;
   uint64_t exponent_y = absy >> 52;
 
-  //uint64_t x_normal_bit = absx >= 0x10000000000000;
-  //uint64_t y_normal_bit = absy >= 0x10000000000000;
+  // uint64_t x_normal_bit = absx >= 0x10000000000000;
+  // uint64_t y_normal_bit = absy >= 0x10000000000000;
 
   uint64_t mx, my;
-  
-  mx = maxu(nlz(absx), 11); 
 
-  //lambdax = mx - 11;  
+  mx = maxu(nlz(absx), 11);
+
+  // lambdax = mx - 11;
 
   my = maxu(nlz(absy), 11);
 
-  //lambday = my - 11;
+  // lambday = my - 11;
 
   int32_t dm1;
-  uint64_t mpx, mpy, highs,lows, b;
-  uint64_t  g, hight, lowt, morlowt, c, m;
+  uint64_t mpx, mpy, highs, lows, b;
+  uint64_t g, hight, lowt, morlowt, c, m;
   mpx = (x_u << mx) | 0x8000000000000000;
-  mpy = (y_u << my) |  0x8000000000000000;
-  highs = mul(mpx, mpy); 
-  c = highs >= 0x8000000000000000; 
+  mpy = (y_u << my) | 0x8000000000000000;
+  highs = mul(mpx, mpy);
+  c = highs >= 0x8000000000000000;
   lows = mullow(mpx, mpy);
 
   lowt = (lows != 0);
@@ -275,62 +273,61 @@ LIBC_INLINE float fmul(double x, double y) {
 
   if (dm1 >= 255) {
     dm1 = 255;
-    
+
     m = 0;
-  }
-  else if (dm1<=0) {
+  } else if (dm1 <= 0) {
     m = static_cast<uint32_t>((highs >> (39 + c)) >> (1 - dm1));
     dm1 = 0;
     morlowt = m | lowt;
     b = g & (morlowt | hight);
-  }
-  else  {
+  } else {
     m = static_cast<uint32_t>(highs >> (39 + c));
     morlowt = m | lowt;
     b = g & (morlowt | hight);
   }
 
   uint32_t sr = static_cast<uint32_t>((x_u ^ y_u) & 0x8000000000000000);
-  
 
   uint32_t exp16 = sr | (static_cast<uint32_t>(dm1) << 23);
   if (dm1 == 0) {
     uint32_t m2 = static_cast<uint32_t>(m);
-    uint32_t result = (static_cast<uint32_t>(exp16) + m2) + static_cast<uint32_t>(b);
-   
-    
+    uint32_t result =
+        (static_cast<uint32_t>(exp16) + m2) + static_cast<uint32_t>(b);
+
     float result32 = cpp::bit_cast<float>(result);
-    
+
     return result32;
-      
- } else {
-   constexpr uint32_t FLOAT32_MANTISSA_MASK = 0b00000000011111111111111111111111;
-   uint32_t m2 = static_cast<uint32_t>(m) & FLOAT32_MANTISSA_MASK;
-   
-   uint32_t result = (static_cast<uint32_t>(exp16) + m2) + static_cast<uint32_t>(b);
-  
-  float result16 = cpp::bit_cast<float>(result);
-  // std::memcpy(&result16, &result, sizeof(result16));  
-  //result = *reinterpret_cast<_Float16*>(&result);
-  
-  return result16;
- }
+
+  } else {
+    constexpr uint32_t FLOAT32_MANTISSA_MASK =
+        0b00000000011111111111111111111111;
+    uint32_t m2 = static_cast<uint32_t>(m) & FLOAT32_MANTISSA_MASK;
+
+    uint32_t result =
+        (static_cast<uint32_t>(exp16) + m2) + static_cast<uint32_t>(b);
+
+    float result16 = cpp::bit_cast<float>(result);
+    // std::memcpy(&result16, &result, sizeof(result16));
+    // result = *reinterpret_cast<_Float16*>(&result);
+
+    return result16;
+  }
 }
-  /*
+/*
 LIBC_INLINE float fmul(double x, double y){
-  FPBits<double> bitx(x), bity(y);
-
-  long long int product;
-  long long int p;
-  
-  if (bitx.is_normal() && bity.is_normal()) {
-    product = bitx.get_mantissa() * bity.get_mantissa();
-    p = product & ((1ULL << 24) -1);
-    sr = bitx.sign() ^ bity.sign();
-    return pow(-1, sr) * pow(p, p.biased_exponent());
-  }
+FPBits<double> bitx(x), bity(y);
+
+long long int product;
+long long int p;
+
+if (bitx.is_normal() && bity.is_normal()) {
+  product = bitx.get_mantissa() * bity.get_mantissa();
+  p = product & ((1ULL << 24) -1);
+  sr = bitx.sign() ^ bity.sign();
+  return pow(-1, sr) * pow(p, p.biased_exponent());
+}
 }
-  */
+*/
 LIBC_INLINE float fmull(long double x, long double y) {
   return static_cast<float>(x * y);
 }
diff --git a/libc/test/src/math/smoke/FMulTest.h b/libc/test/src/math/smoke/FMulTest.h
index afaca64d866c1..44cab1e2560a7 100644
--- a/libc/test/src/math/smoke/FMulTest.h
+++ b/libc/test/src/math/smoke/FMulTest.h
@@ -21,9 +21,9 @@ class FmulTest : public LIBC_NAMESPACE::testing::FEnvSafeTest {
 
   void testMul(FMulFunc func) {
     EXPECT_FP_EQ(T(15.0), func(3.0, 5.0));
-    EXPECT_FP_EQ(T(7.34684e-40), func(0x1.0p1,0x1.0p-131));
+    EXPECT_FP_EQ(T(7.34684e-40), func(0x1.0p1, 0x1.0p-131));
     EXPECT_FP_EQ(T(5.87747e-39), func(0x1.0p2, 0x1.0p-129));
-}
+  }
 };
 
 #define LIST_FMUL_TESTS(T, R, func)                                            \

>From 54dee21e4a904ae331b45dd41e66ae24ed458c79 Mon Sep 17 00:00:00 2001
From: Job Hernandez <hj93 at protonmail.com>
Date: Sun, 2 Jun 2024 12:28:39 -0700
Subject: [PATCH 5/8] move fmul function to generic/fmul.cpp

---
 libc/src/__support/FPUtil/BasicOperations.h | 150 --------------------
 libc/src/math/generic/fmul.cpp              | 130 ++++++++++++++++-
 2 files changed, 129 insertions(+), 151 deletions(-)

diff --git a/libc/src/__support/FPUtil/BasicOperations.h b/libc/src/__support/FPUtil/BasicOperations.h
index 1926659e5f78f..c63f4d35eaa61 100644
--- a/libc/src/__support/FPUtil/BasicOperations.h
+++ b/libc/src/__support/FPUtil/BasicOperations.h
@@ -178,156 +178,6 @@ LIBC_INLINE T fdim(T x, T y) {
   return (x > y ? x - y : 0);
 }
 
-// helpers for fmul
-
-uint64_t maxu(uint64_t A, uint64_t B) { return A > B ? A : B; }
-
-uint64_t mul(uint64_t a, uint64_t b) {
-  __uint128_t product =
-      static_cast<__uint128_t>(a) * static_cast<__uint128_t>(b);
-  return static_cast<uint64_t>(product >> 64);
-}
-
-uint64_t mullow(uint64_t a, uint64_t b) {
-  __uint128_t product =
-      static_cast<__uint128_t>(a) * static_cast<__uint128_t>(b);
-  return static_cast<uint64_t>(product);
-}
-
-uint64_t nlz(uint64_t x) {
-  uint64_t z = 0;
-
-  if (x == 0)
-    return 64;
-  if (x <= 0x00000000FFFFFFFF) {
-    z = z + 32;
-    x = x << 32;
-  }
-  if (x <= 0x0000FFFFFFFFFFFF) {
-    z = z + 16;
-    x = x << 16;
-  }
-  if (x <= 0x00FFFFFFFFFFFFFF) {
-    z = z + 8;
-    x = x << 8;
-  }
-  if (x <= 0x0FFFFFFFFFFFFFFF) {
-    z = z + 4;
-    x = x << 4;
-  }
-  if (x <= 0x3FFFFFFFFFFFFFFF) {
-    z = z + 2;
-    x = x << 2;
-  }
-  if (x <= 0x7FFFFFFFFFFFFFFF) {
-    z = z + 1;
-  }
-  return z;
-}
-
-LIBC_INLINE float fmul(double x, double y) {
-
-  auto x_bits = FPBits<double>(x);
-  uint64_t x_u = x_bits.uintval();
-
-  auto y_bits = FPBits<double>(y);
-  uint64_t y_u = y_bits.uintval();
-
-  uint64_t absx = x_u & 0x7FFFFFFFFFFFFFFF;
-  uint64_t absy = y_u & 0x7FFFFFFFFFFFFFFF;
-
-  uint64_t exponent_x = absx >> 52;
-  uint64_t exponent_y = absy >> 52;
-
-  // uint64_t x_normal_bit = absx >= 0x10000000000000;
-  // uint64_t y_normal_bit = absy >= 0x10000000000000;
-
-  uint64_t mx, my;
-
-  mx = maxu(nlz(absx), 11);
-
-  // lambdax = mx - 11;
-
-  my = maxu(nlz(absy), 11);
-
-  // lambday = my - 11;
-
-  int32_t dm1;
-  uint64_t mpx, mpy, highs, lows, b;
-  uint64_t g, hight, lowt, morlowt, c, m;
-  mpx = (x_u << mx) | 0x8000000000000000;
-  mpy = (y_u << my) | 0x8000000000000000;
-  highs = mul(mpx, mpy);
-  c = highs >= 0x8000000000000000;
-  lows = mullow(mpx, mpy);
-
-  lowt = (lows != 0);
-
-  g = (highs >> (38 + c)) & 1;
-  hight = (highs << (55 - c)) != 0;
-
-  int32_t exint = static_cast<uint32_t>(exponent_x);
-  int32_t eyint = static_cast<uint32_t>(exponent_y);
-  int32_t cint = static_cast<uint32_t>(c);
-  dm1 = ((exint + eyint) - 1919) + cint;
-
-  if (dm1 >= 255) {
-    dm1 = 255;
-
-    m = 0;
-  } else if (dm1 <= 0) {
-    m = static_cast<uint32_t>((highs >> (39 + c)) >> (1 - dm1));
-    dm1 = 0;
-    morlowt = m | lowt;
-    b = g & (morlowt | hight);
-  } else {
-    m = static_cast<uint32_t>(highs >> (39 + c));
-    morlowt = m | lowt;
-    b = g & (morlowt | hight);
-  }
-
-  uint32_t sr = static_cast<uint32_t>((x_u ^ y_u) & 0x8000000000000000);
-
-  uint32_t exp16 = sr | (static_cast<uint32_t>(dm1) << 23);
-  if (dm1 == 0) {
-    uint32_t m2 = static_cast<uint32_t>(m);
-    uint32_t result =
-        (static_cast<uint32_t>(exp16) + m2) + static_cast<uint32_t>(b);
-
-    float result32 = cpp::bit_cast<float>(result);
-
-    return result32;
-
-  } else {
-    constexpr uint32_t FLOAT32_MANTISSA_MASK =
-        0b00000000011111111111111111111111;
-    uint32_t m2 = static_cast<uint32_t>(m) & FLOAT32_MANTISSA_MASK;
-
-    uint32_t result =
-        (static_cast<uint32_t>(exp16) + m2) + static_cast<uint32_t>(b);
-
-    float result16 = cpp::bit_cast<float>(result);
-    // std::memcpy(&result16, &result, sizeof(result16));
-    // result = *reinterpret_cast<_Float16*>(&result);
-
-    return result16;
-  }
-}
-/*
-LIBC_INLINE float fmul(double x, double y){
-FPBits<double> bitx(x), bity(y);
-
-long long int product;
-long long int p;
-
-if (bitx.is_normal() && bity.is_normal()) {
-  product = bitx.get_mantissa() * bity.get_mantissa();
-  p = product & ((1ULL << 24) -1);
-  sr = bitx.sign() ^ bity.sign();
-  return pow(-1, sr) * pow(p, p.biased_exponent());
-}
-}
-*/
 LIBC_INLINE float fmull(long double x, long double y) {
   return static_cast<float>(x * y);
 }
diff --git a/libc/src/math/generic/fmul.cpp b/libc/src/math/generic/fmul.cpp
index 77452013951d7..f618a168bb42b 100644
--- a/libc/src/math/generic/fmul.cpp
+++ b/libc/src/math/generic/fmul.cpp
@@ -8,12 +8,140 @@
 
 #include "src/math/fmul.h"
 #include "src/__support/FPUtil/BasicOperations.h"
+#include "src/__support/FPUtil/FPBits.h"
 #include "src/__support/common.h"
 
 namespace LIBC_NAMESPACE {
 
+namespace Fmul  {
+uint64_t maxu(uint64_t A, uint64_t B) { return A > B ? A : B; }
+
+uint64_t mul(uint64_t a, uint64_t b) {
+  __uint128_t product =
+      static_cast<__uint128_t>(a) * static_cast<__uint128_t>(b);
+  return static_cast<uint64_t>(product >> 64);
+}
+
+uint64_t mullow(uint64_t a, uint64_t b) {
+  __uint128_t product =
+      static_cast<__uint128_t>(a) * static_cast<__uint128_t>(b);
+  return static_cast<uint64_t>(product);
+}
+
+uint64_t nlz(uint64_t x) {
+  uint64_t z = 0;
+
+  if (x == 0)
+    return 64;
+  if (x <= 0x00000000FFFFFFFF) {
+    z = z + 32;
+    x = x << 32;
+  }
+  if (x <= 0x0000FFFFFFFFFFFF) {
+    z = z + 16;
+    x = x << 16;
+  }
+  if (x <= 0x00FFFFFFFFFFFFFF) {
+    z = z + 8;
+    x = x << 8;
+  }
+  if (x <= 0x0FFFFFFFFFFFFFFF) {
+    z = z + 4;
+    x = x << 4;
+  }
+  if (x <= 0x3FFFFFFFFFFFFFFF) {
+    z = z + 2;
+    x = x << 2;
+  }
+  if (x <= 0x7FFFFFFFFFFFFFFF) {
+    z = z + 1;
+  }
+  return z;
+}
+
+float fmul(double x, double y) {
+
+  auto x_bits = fputil::FPBits<double>(x);
+  uint64_t x_u = x_bits.uintval();
+
+  auto y_bits = fputil::FPBits<double>(y);
+  uint64_t y_u = y_bits.uintval();
+
+  uint64_t absx = x_u & 0x7FFFFFFFFFFFFFFF;
+  uint64_t absy = y_u & 0x7FFFFFFFFFFFFFFF;
+
+  uint64_t exponent_x = absx >> 52;
+  uint64_t exponent_y = absy >> 52;
+
+  uint64_t mx, my;
+
+  mx = maxu(nlz(absx), 11);
+
+  my = maxu(nlz(absy), 11);
+
+  int32_t dm1;
+  uint64_t mpx, mpy, highs, lows, b;
+  uint64_t g, hight, lowt, morlowt, c, m;
+  mpx = (x_u << mx) | 0x8000000000000000;
+  mpy = (y_u << my) | 0x8000000000000000;
+  highs = mul(mpx, mpy);
+  c = highs >= 0x8000000000000000;
+  lows = mullow(mpx, mpy);
+
+  lowt = (lows != 0);
+
+  g = (highs >> (38 + c)) & 1;
+  hight = (highs << (55 - c)) != 0;
+
+  int32_t exint = static_cast<int32_t>(exponent_x);
+  int32_t eyint = static_cast<int32_t>(exponent_y);
+  int32_t cint = static_cast<int32_t>(c);
+  dm1 = ((exint + eyint) - 1919) + cint;
+
+  if (dm1 >= 255) {
+    dm1 = 255;
+
+    m = 0;
+  } else if (dm1 <= 0) {
+    m = static_cast<uint32_t>((highs >> (39 + c)) >> (1 - dm1));
+    dm1 = 0;
+    morlowt = m | lowt;
+    b = g & (morlowt | hight);
+  } else {
+    m = static_cast<uint32_t>(highs >> (39 + c));
+    morlowt = m | lowt;
+    b = g & (morlowt | hight);
+  }
+
+  uint32_t sr = static_cast<uint32_t>((x_u ^ y_u) & 0x8000000000000000);
+
+  uint32_t exp16 = sr | (dm1 << 23);
+  if (dm1 == 0) {
+    uint32_t m2 = static_cast<uint32_t>(m);
+    uint32_t result =
+        (static_cast<uint32_t>(exp16) + m2) + static_cast<uint32_t>(b);
+
+    float result32 = cpp::bit_cast<float>(result);
+
+    return result32;
+
+  } else {
+    constexpr uint32_t FLOAT32_MANTISSA_MASK =
+        0b00000000011111111111111111111111;
+    uint32_t m2 = static_cast<uint32_t>(m) & FLOAT32_MANTISSA_MASK;
+
+    uint32_t result =
+        (static_cast<uint32_t>(exp16) + m2) + static_cast<uint32_t>(b);
+
+    float result16 = cpp::bit_cast<float>(result);
+
+    return result16;
+  }
+}
+}
+
 LLVM_LIBC_FUNCTION(float, fmul, (double x, double y)) {
-  return fputil::fmul(x, y);
+  return Fmul::fmul(x, y);
 }
 
 } // namespace LIBC_NAMESPACE

>From 3bf640d06c04698ce5b7f98d6e1add86d76d50a8 Mon Sep 17 00:00:00 2001
From: Job Hernandez <hj93 at protonmail.com>
Date: Sun, 2 Jun 2024 12:28:54 -0700
Subject: [PATCH 6/8] fix tests

---
 libc/test/src/math/smoke/FMulTest.h | 5 +++--
 1 file changed, 3 insertions(+), 2 deletions(-)

diff --git a/libc/test/src/math/smoke/FMulTest.h b/libc/test/src/math/smoke/FMulTest.h
index 44cab1e2560a7..881ed00c8bd67 100644
--- a/libc/test/src/math/smoke/FMulTest.h
+++ b/libc/test/src/math/smoke/FMulTest.h
@@ -21,8 +21,9 @@ class FmulTest : public LIBC_NAMESPACE::testing::FEnvSafeTest {
 
   void testMul(FMulFunc func) {
     EXPECT_FP_EQ(T(15.0), func(3.0, 5.0));
-    EXPECT_FP_EQ(T(7.34684e-40), func(0x1.0p1, 0x1.0p-131));
-    EXPECT_FP_EQ(T(5.87747e-39), func(0x1.0p2, 0x1.0p-129));
+    EXPECT_FP_EQ(T(0x1.0p-130), func(0x1.0p1, 0x1.0p-131));
+    EXPECT_FP_EQ(T(0x1.0p-127), func(0x1.0p2, 0x1.0p-129));
+    
   }
 };
 

>From 0c290cef9a0024d27885c81582da76e254c2ebb0 Mon Sep 17 00:00:00 2001
From: Job Hernandez <hj93 at protonmail.com>
Date: Sun, 2 Jun 2024 12:29:39 -0700
Subject: [PATCH 7/8] clang format

---
 libc/src/math/generic/fmul.cpp      | 4 ++--
 libc/test/src/math/smoke/FMulTest.h | 1 -
 2 files changed, 2 insertions(+), 3 deletions(-)

diff --git a/libc/src/math/generic/fmul.cpp b/libc/src/math/generic/fmul.cpp
index f618a168bb42b..e543888b0dc93 100644
--- a/libc/src/math/generic/fmul.cpp
+++ b/libc/src/math/generic/fmul.cpp
@@ -13,7 +13,7 @@
 
 namespace LIBC_NAMESPACE {
 
-namespace Fmul  {
+namespace Fmul {
 uint64_t maxu(uint64_t A, uint64_t B) { return A > B ? A : B; }
 
 uint64_t mul(uint64_t a, uint64_t b) {
@@ -138,7 +138,7 @@ float fmul(double x, double y) {
     return result16;
   }
 }
-}
+} // namespace Fmul
 
 LLVM_LIBC_FUNCTION(float, fmul, (double x, double y)) {
   return Fmul::fmul(x, y);
diff --git a/libc/test/src/math/smoke/FMulTest.h b/libc/test/src/math/smoke/FMulTest.h
index 881ed00c8bd67..6bfb66eaaf869 100644
--- a/libc/test/src/math/smoke/FMulTest.h
+++ b/libc/test/src/math/smoke/FMulTest.h
@@ -23,7 +23,6 @@ class FmulTest : public LIBC_NAMESPACE::testing::FEnvSafeTest {
     EXPECT_FP_EQ(T(15.0), func(3.0, 5.0));
     EXPECT_FP_EQ(T(0x1.0p-130), func(0x1.0p1, 0x1.0p-131));
     EXPECT_FP_EQ(T(0x1.0p-127), func(0x1.0p2, 0x1.0p-129));
-    
   }
 };
 

>From c33af75640da42b54d588e2b83af8ce16d5b23a9 Mon Sep 17 00:00:00 2001
From: Job Hernandez <hj93 at protonmail.com>
Date: Sun, 2 Jun 2024 17:16:53 -0700
Subject: [PATCH 8/8] add more test cases

---
 libc/src/math/generic/fmul.cpp      | 23 +++++++++++++++++++++++
 libc/test/src/math/smoke/FMulTest.h |  8 ++++++++
 2 files changed, 31 insertions(+)

diff --git a/libc/src/math/generic/fmul.cpp b/libc/src/math/generic/fmul.cpp
index e543888b0dc93..6aaeb1b32141c 100644
--- a/libc/src/math/generic/fmul.cpp
+++ b/libc/src/math/generic/fmul.cpp
@@ -61,12 +61,35 @@ uint64_t nlz(uint64_t x) {
 
 float fmul(double x, double y) {
 
+  
   auto x_bits = fputil::FPBits<double>(x);
   uint64_t x_u = x_bits.uintval();
 
   auto y_bits = fputil::FPBits<double>(y);
   uint64_t y_u = y_bits.uintval();
 
+  if ((x_bits.is_subnormal() || x_bits.is_normal()) && y_bits.is_inf())
+    return y_bits.inf().get_val();
+  if (x_bits.is_inf() && (y_bits.is_subnormal() || y_bits.is_normal()))
+    return x_bits.inf().get_val();
+  if ((x_bits.is_subnormal() || x_bits.is_normal()) && y_bits.is_zero())
+    return y_bits.zero().get_val();
+  if (x_bits.is_zero() && (y_bits.is_subnormal() || y_bits.is_normal()))
+    return x_bits.zero().get_val();
+  if (x_bits.is_zero() && y_bits.is_zero())
+    return x_bits.zero().get_val();
+  if (x_bits.is_inf() && y_bits.is_zero())
+    return fputil::FPBits<float>::quiet_nan().get_val();
+  if (y_bits.is_inf() && x_bits.is_zero())
+    return fputil::FPBits<float>::quiet_nan().get_val();
+  if ((x_bits.is_zero() || x_bits.is_normal() || x_bits.is_subnormal() || x_bits.is_inf()) && y_bits.is_nan()) {
+    fputil::raise_except_if_required(FE_INVALID);
+    return fputil::FPBits<float>::quiet_nan().get_val();
+  }
+  if ((y_bits.is_zero() || x_bits.is_normal() || x_bits.is_subnormal() || x_bits.is_inf()) &&  x_bits.is_nan()) {
+    fputil::raise_except_if_required(FE_INVALID);
+    return fputil::FPBits<float>::quiet_nan().get_val();
+  }
   uint64_t absx = x_u & 0x7FFFFFFFFFFFFFFF;
   uint64_t absy = y_u & 0x7FFFFFFFFFFFFFFF;
 
diff --git a/libc/test/src/math/smoke/FMulTest.h b/libc/test/src/math/smoke/FMulTest.h
index 6bfb66eaaf869..c0e5650f53dc4 100644
--- a/libc/test/src/math/smoke/FMulTest.h
+++ b/libc/test/src/math/smoke/FMulTest.h
@@ -16,6 +16,8 @@
 template <typename T, typename R>
 class FmulTest : public LIBC_NAMESPACE::testing::FEnvSafeTest {
 
+  DECLARE_SPECIAL_CONSTANTS(T)
+
 public:
   typedef T (*FMulFunc)(R, R);
 
@@ -23,6 +25,12 @@ class FmulTest : public LIBC_NAMESPACE::testing::FEnvSafeTest {
     EXPECT_FP_EQ(T(15.0), func(3.0, 5.0));
     EXPECT_FP_EQ(T(0x1.0p-130), func(0x1.0p1, 0x1.0p-131));
     EXPECT_FP_EQ(T(0x1.0p-127), func(0x1.0p2, 0x1.0p-129));
+    EXPECT_FP_EQ(inf, func(inf, 0x1.0p-129));
+    EXPECT_FP_EQ(inf, func(0x1.0p-129, inf));
+    EXPECT_FP_EQ(inf, func(inf, 2.0));
+    EXPECT_FP_EQ(inf, func(3.0, inf));
+    EXPECT_FP_EQ(0.0, func(0.0, 0.0));
+    EXPECT_FP_EQ(0.0, func(0.0, 0x1.0p-129));
   }
 };
 



More information about the libc-commits mailing list