[libc] Move implementations of expf and exp2f from the AOR to src/math.

[sivachandra]

Summary

Depends on D79150, D79278.

[sivachandra]

This change is not quite ready but feel free to share any opinions/comments you might have on the changes I am making. I started with moving exp and all its friends but soon ended with a large change which even I was getting confused about. Also, I made a mistake with how LLVM libc results are compared with MPFR results. I have to fix that first before I can add more tests to this patch.

[abrachet]

libc/src/math/exp2f.cpp

17–19

These are only used once or twice and we just expand these manually?

26

We probably want these to be //

libc/src/math/exp2f.h

1–2

Classic clang-format

10–11

MATH_EXP2F

libc/src/math/exp_utils.cpp

1

Missing comment

libc/src/math/expf.h

9–11

MATH_EXPF

[sivachandra]

• Rebase and add tests
• Address comments from @abrachet

[sivachandra]

I am mostly done with everything I wanted to do in this patch but I am still leaving it without reviewers as I want to do one more pass tomorrow before I add reviewers. Like with the other math change, I want to keep the general look of the code similar to how it is in the AOR directory. So, I have not made widespread stylistic changes which would change the look and feel of the actual code. That said, as I start moving the double precision functions, it might be more appropriate to consolidate the now dispersed floating point utils into a single, more C++ like util library. I will start doing it in the following patches as I learn more about the double precision implementations.

[sivachandra]

Added reviewers now. Just pointing out again: I have not done any stylistic changes to the main code. The idea is to keep it as familiar as possible to the AOR developers when they start working with LLVM-libc.

[phosek]

LGTM

libc/utils/MPFRWrapper/MPFRUtils.h

38–39

This could be enum class and we could avoid the OP_ prefix but this is fine as well.

libc/config/linux/api.td

let Macros = [
     IsNanMacro,
   ];
   let TypeDeclarations = [
     DoubleT,
     FloatT,
   ];
   let Functions = [
    "cosf",
+   "expf",
+   "exp2f",
    "round",
    "sincosf",
    "sinf",
   ];
 }
 
 def StringAPI : PublicAPI<"string.h"> {
   let Functions = [

libc/spec/stdc.td

HeaderSpec Math = HeaderSpec<
           NamedType<"float_t">,
           NamedType<"double_t">,
       ],
       [], // Enumerations
       [
           FunctionSpec<"cosf", RetValSpec<FloatType>, [ArgSpec<FloatType>]>,
           FunctionSpec<"sinf", RetValSpec<FloatType>, [ArgSpec<FloatType>]>,
 
+          FunctionSpec<"expf", RetValSpec<FloatType>, [ArgSpec<FloatType>]>,
+          FunctionSpec<"exp2f", RetValSpec<FloatType>, [ArgSpec<FloatType>]>,
+
           FunctionSpec<"round", RetValSpec<DoubleType>, [ArgSpec<DoubleType>]>,
       ]
   >;
 
   HeaderSpec StdIO = HeaderSpec<
       "stdio.h",
       [], // Macros
       [ // Types

libc/src/math/CMakeLists.txt

SRCS
     sincosf.cpp
   HDRS
     sincosf.h
   DEPENDS
     .sincosf_utils
     libc.include.math
     libc.src.errno.__errno_location
 )
+
+add_object_library(
+  exp_utils
+  HDRS
+    exp_utils.h
+  SRCS
+    exp_utils.cpp
+)
+
+add_entrypoint_object(
+  expf
+  SRCS
+    expf.cpp
+  HDRS
+    expf.h
+  DEPENDS
+    .exp_utils
+    .math_utils
+    libc.include.math
+)
+
+add_entrypoint_object(
+  exp2f
+  SRCS
+    exp2f.cpp
+  HDRS
+    exp2f.h
+  DEPENDS
+    .exp_utils
+    .math_utils
+    libc.include.math
+)

libc/src/math/cosf.cpp

if (abstop12(y) < abstop12(pio4)) {
     s = p->sign[(n + sign) & 3];
 
     if ((n + sign) & 2)
       p = &__sincosf_table[1];
 
     return sinf_poly(x * s, x * x, p, n ^ 1);
   }
 
-  return invalidf(y);
+  return invalid(y);
 }
 
 } // namespace __llvm_libc

libc/src/math/exp2f.h

+//===-- Implementation header for exp2f --------------------------*- 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_SINF_H
+#define LLVM_LIBC_SRC_MATH_SINF_H
+
+namespace __llvm_libc {
+
+float exp2f(float x);
+
+} // namespace __llvm_libc
+
+#endif // LLVM_LIBC_SRC_MATH_SINF_H

libc/src/math/exp2f.cpp

+//===-- Single-precision 2^x 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 "exp_utils.h"
+#include "math_utils.h"
+
+#include "include/math.h"
+#include "src/__support/common.h"
+
+#include <stdint.h>
+
+#define T exp2f_data.tab
+#define C exp2f_data.poly
+#define SHIFT exp2f_data.shift_scaled
+
+namespace __llvm_libc {
+
+float LLVM_LIBC_ENTRYPOINT(exp2f)(float x) {
+  uint32_t abstop;
+  uint64_t ki, t;
+  /* double_t for better performance on targets with FLT_EVAL_METHOD==2.  */
+  double_t kd, xd, z, r, r2, y, s;
+
+  xd = static_cast<double_t>(x);
+  abstop = top12_bits(x) & 0x7ff;
+  if (unlikely(abstop >= top12_bits(128.0f))) {
+    // |x| >= 128 or x is nan.
+    if (as_uint32_bits(x) == as_uint32_bits(-INFINITY))
+      return 0.0f;
+    if (abstop >= top12_bits(INFINITY))
+      return x + x;
+    if (x > 0.0f)
+      return overflow<float>(0);
+    if (x <= -150.0f)
+      return underflow<float>(0);
+    if (x < -149.0f)
+      return may_underflow<float>(0);
+  }
+
+  /* x = k/N + r with r in [-1/(2N), 1/(2N)] and int k.  */
+  kd = static_cast<double>(xd + SHIFT);
+  ki = as_uint64_bits(kd);
+  kd -= SHIFT; /* k/N for int k.  */
+  r = xd - kd;
+
+  /* exp2(x) = 2^(k/N) * 2^r ~= s * (C0*r^3 + C1*r^2 + C2*r + 1) */
+  t = T[ki % N];
+  t += ki << (52 - EXP2F_TABLE_BITS);
+  s = as_double(t);
+  z = C[0] * r + C[1];
+  r2 = r * r;
+  y = C[2] * r + 1;
+  y = z * r2 + y;
+  y = y * s;
+  return static_cast<float>(y);
+}
+
+} // namespace __llvm_libc

libc/src/math/exp_utils.h

+//===-- Collection of utils for exp and friends -----------------*- 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_EXP_UTILS_H
+#define LLVM_LIBC_SRC_MATH_EXP_UTILS_H
+
+#include "math_utils.h"
+#include <stdint.h>
+
+#define EXP2F_TABLE_BITS 5
+#define EXP2F_POLY_ORDER 3
+#define N (1 << EXP2F_TABLE_BITS)
+
+namespace __llvm_libc {
+
+struct Exp2fDataTable {
+  uint64_t tab[1 << EXP2F_TABLE_BITS];
+  double shift_scaled;
+  double poly[EXP2F_POLY_ORDER];
+  double shift;
+  double invln2_scaled;
+  double poly_scaled[EXP2F_POLY_ORDER];
+};
+
+extern const Exp2fDataTable exp2f_data;
+
+} // namespace __llvm_libc
+
+#endif // LLVM_LIBC_SRC_MATH_EXP_UTILS_H

libc/src/math/exp_utils.cpp

+//===----------------------------------------------------------------------===//
+//
+// 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 "exp_utils.h"
+
+namespace __llvm_libc {
+
+constexpr Exp2fDataTable exp2f_data = {
+    /* :tab[i] = uint(2^(i/N)) - (i << 52-BITS)
+        used for computing 2^(k/N) for an int |k| < 150 N as
+        double(tab[k%N] + (k << 52-BITS)) */
+    .tab =
+        {
+#if N == 8
+            0x3ff0000000000000,
+            0x3fef72b83c7d517b,
+            0x3fef06fe0a31b715,
+            0x3feebfdad5362a27,
+            0x3feea09e667f3bcd,
+            0x3feeace5422aa0db,
+            0x3feee89f995ad3ad,
+            0x3fef5818dcfba487,
+#elif N == 16
+            0x3ff0000000000000,
+            0x3fefb5586cf9890f,
+            0x3fef72b83c7d517b,
+            0x3fef387a6e756238,
+            0x3fef06fe0a31b715,
+            0x3feedea64c123422,
+            0x3feebfdad5362a27,
+            0x3feeab07dd485429,
+            0x3feea09e667f3bcd,
+            0x3feea11473eb0187,
+            0x3feeace5422aa0db,
+            0x3feec49182a3f090,
+            0x3feee89f995ad3ad,
+            0x3fef199bdd85529c,
+            0x3fef5818dcfba487,
+            0x3fefa4afa2a490da,
+#elif N == 32
+            0x3ff0000000000000, 0x3fefd9b0d3158574, 0x3fefb5586cf9890f,
+            0x3fef9301d0125b51, 0x3fef72b83c7d517b, 0x3fef54873168b9aa,
+            0x3fef387a6e756238, 0x3fef1e9df51fdee1, 0x3fef06fe0a31b715,
+            0x3feef1a7373aa9cb, 0x3feedea64c123422, 0x3feece086061892d,
+            0x3feebfdad5362a27, 0x3feeb42b569d4f82, 0x3feeab07dd485429,
+            0x3feea47eb03a5585, 0x3feea09e667f3bcd, 0x3fee9f75e8ec5f74,
+            0x3feea11473eb0187, 0x3feea589994cce13, 0x3feeace5422aa0db,
+            0x3feeb737b0cdc5e5, 0x3feec49182a3f090, 0x3feed503b23e255d,
+            0x3feee89f995ad3ad, 0x3feeff76f2fb5e47, 0x3fef199bdd85529c,
+            0x3fef3720dcef9069, 0x3fef5818dcfba487, 0x3fef7c97337b9b5f,
+            0x3fefa4afa2a490da, 0x3fefd0765b6e4540,
+#elif N == 64
+            0x3ff0000000000000, 0x3fefec9a3e778061, 0x3fefd9b0d3158574,
+            0x3fefc74518759bc8, 0x3fefb5586cf9890f, 0x3fefa3ec32d3d1a2,
+            0x3fef9301d0125b51, 0x3fef829aaea92de0, 0x3fef72b83c7d517b,
+            0x3fef635beb6fcb75, 0x3fef54873168b9aa, 0x3fef463b88628cd6,
+            0x3fef387a6e756238, 0x3fef2b4565e27cdd, 0x3fef1e9df51fdee1,
+            0x3fef1285a6e4030b, 0x3fef06fe0a31b715, 0x3feefc08b26416ff,
+            0x3feef1a7373aa9cb, 0x3feee7db34e59ff7, 0x3feedea64c123422,
+            0x3feed60a21f72e2a, 0x3feece086061892d, 0x3feec6a2b5c13cd0,
+            0x3feebfdad5362a27, 0x3feeb9b2769d2ca7, 0x3feeb42b569d4f82,
+            0x3feeaf4736b527da, 0x3feeab07dd485429, 0x3feea76f15ad2148,
+            0x3feea47eb03a5585, 0x3feea23882552225, 0x3feea09e667f3bcd,
+            0x3fee9fb23c651a2f, 0x3fee9f75e8ec5f74, 0x3fee9feb564267c9,
+            0x3feea11473eb0187, 0x3feea2f336cf4e62, 0x3feea589994cce13,
+            0x3feea8d99b4492ed, 0x3feeace5422aa0db, 0x3feeb1ae99157736,
+            0x3feeb737b0cdc5e5, 0x3feebd829fde4e50, 0x3feec49182a3f090,
+            0x3feecc667b5de565, 0x3feed503b23e255d, 0x3feede6b5579fdbf,
+            0x3feee89f995ad3ad, 0x3feef3a2b84f15fb, 0x3feeff76f2fb5e47,
+            0x3fef0c1e904bc1d2, 0x3fef199bdd85529c, 0x3fef27f12e57d14b,
+            0x3fef3720dcef9069, 0x3fef472d4a07897c, 0x3fef5818dcfba487,
+            0x3fef69e603db3285, 0x3fef7c97337b9b5f, 0x3fef902ee78b3ff6,
+            0x3fefa4afa2a490da, 0x3fefba1bee615a27, 0x3fefd0765b6e4540,
+            0x3fefe7c1819e90d8,
+#endif
+        },
+    .shift_scaled = 0x1.8p+52 / N,
+    .poly =
+        {
+#if N == 8
+            0x1.c6a00335106e2p-5,
+            0x1.ec0c313449f55p-3,
+            0x1.62e431111f69fp-1,
+#elif N == 16
+            0x1.c6ac6aa313963p-5,
+            0x1.ebfff4532d9bap-3,
+            0x1.62e43001bc49fp-1,
+#elif N == 32
+            0x1.c6af84b912394p-5,
+            0x1.ebfce50fac4f3p-3,
+            0x1.62e42ff0c52d6p-1,
+#elif N == 64
+            0x1.c6b04b4221b2ap-5,
+            0x1.ebfc213e184d7p-3,
+            0x1.62e42fefb5b7fp-1,
+#endif
+        },
+    .shift = 0x1.8p+52,
+    .invln2_scaled = 0x1.71547652b82fep+0 * N,
+    .poly_scaled =
+        {
+#if N == 8
+            0x1.c6a00335106e2p-5 / N / N / N,
+            0x1.ec0c313449f55p-3 / N / N,
+            0x1.62e431111f69fp-1 / N,
+#elif N == 16
+            0x1.c6ac6aa313963p-5 / N / N / N,
+            0x1.ebfff4532d9bap-3 / N / N,
+            0x1.62e43001bc49fp-1 / N,
+#elif N == 32
+            0x1.c6af84b912394p-5 / N / N / N,
+            0x1.ebfce50fac4f3p-3 / N / N,
+            0x1.62e42ff0c52d6p-1 / N,
+#elif N == 64
+            0x1.c6b04b4221b2ap-5 / N / N / N,
+            0x1.ebfc213e184d7p-3 / N / N,
+            0x1.62e42fefb5b7fp-1 / N,
+#endif
+        },
+};
+
+} // namespace __llvm_libc

libc/src/math/expf.h

+//===-- Implementation header for expf --------------------------*- 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_SINF_H
+#define LLVM_LIBC_SRC_MATH_SINF_H
+
+namespace __llvm_libc {
+
+float expf(float x);
+
+} // namespace __llvm_libc
+
+#endif // LLVM_LIBC_SRC_MATH_SINF_H

libc/src/math/expf.cpp

+//===-- Single-precision e^x 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 "exp_utils.h"
+#include "math_utils.h"
+
+#include "include/math.h"
+#include "src/__support/common.h"
+
+#include <stdint.h>
+
+#define InvLn2N exp2f_data.invln2_scaled
+#define T exp2f_data.tab
+#define C exp2f_data.poly_scaled
+#define SHIFT exp2f_data.shift
+
+namespace __llvm_libc {
+
+float LLVM_LIBC_ENTRYPOINT(expf)(float x) {
+  uint32_t abstop;
+  uint64_t ki, t;
+  // double_t for better performance on targets with FLT_EVAL_METHOD == 2.
+  double_t kd, xd, z, r, r2, y, s;
+
+  xd = static_cast<double_t>(x);
+  abstop = top12_bits(x) & 0x7ff;
+  if (unlikely(abstop >= top12_bits(88.0f))) {
+    // |x| >= 88 or x is nan.
+    if (as_uint32_bits(x) == as_uint32_bits(-INFINITY))
+      return 0.0f;
+    if (abstop >= top12_bits(INFINITY))
+      return x + x;
+    if (x > 0x1.62e42ep6f) // x > log(0x1p128) ~= 88.72
+      return overflow<float>(0);
+    if (x < -0x1.9fe368p6f) // x < log(0x1p-150) ~= -103.97
+      return underflow<float>(0);
+    if (x < -0x1.9d1d9ep6f) // x < log(0x1p-149) ~= -103.28
+      return may_underflow<float>(0);
+  }
+
+  // x*N/Ln2 = k + r with r in [-1/2, 1/2] and int k.
+  z = InvLn2N * xd;
+
+  // Round and convert z to int, the result is in [-150*N, 128*N] and
+  // ideally nearest int is used, otherwise the magnitude of r can be
+  // bigger which gives larger approximation error.
+  kd = static_cast<double>(z + SHIFT);
+  ki = as_uint64_bits(kd);
+  kd -= SHIFT;
+  r = z - kd;
+
+  // exp(x) = 2^(k/N) * 2^(r/N) ~= s *(C0*r^3 + C1*r^2 + C2*r + 1)
+  t = T[ki % N];
+  t += ki << (52 - EXP2F_TABLE_BITS);
+  s = as_double(t);
+  z = C[0] * r + C[1];
+  r2 = r * r;
+  y = C[2] * r + 1;
+  y = z * r2 + y;
+  y = y * s;
+  return static_cast<float>(y);
+}
+
+} // namespace __llvm_libc

libc/src/math/math_utils.h

 //===-- Collection of utils for implementing math functions -----*- 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_MATH_UTILS_H
 #define LLVM_LIBC_SRC_MATH_MATH_UTILS_H
 
 #include "include/errno.h"
 #include "include/math.h"
-
 #include "src/__support/common.h"
 #include "src/errno/llvmlibc_errno.h"
+#include "utils/CPP/TypeTraits.h"
 
 #include <stdint.h>
 
 namespace __llvm_libc {
 
-static inline float with_errnof(float x, int err) {
+template <typename T> struct IsFloatOrDouble {
+  static constexpr bool Value =
+      cpp::IsSame<T, float>::Value || cpp::IsSame<T, double>::Value;
+};
+
+template <typename T>
+using EnableIfFloatOrDouble = cpp::EnableIfType<IsFloatOrDouble<T>::Value, int>;
+
+template <typename T> struct XFlowValues;
+
+// Single-precission floating point values to trigger underflow or overflow.
+template <> struct XFlowValues<float> {
+  static constexpr float overflow_value = 0x1p97f;
+  static constexpr float underflow_value = 0x1p-95f;
+  static constexpr float may_underflow_value = 0x1.4p-75f;
+};
+
+// Double-precission floating point values to trigger underflow or overflow.
+template <> struct XFlowValues<double> {
+  static constexpr double overflow_value = 0x1p769;
+  static constexpr double underflow_value = 0x1p-767;
+  static constexpr double may_underflow_value = 0x1.8p-538;
+};
+
+template <typename T> static inline T with_errno(T x, int err) {
   if (math_errhandling & MATH_ERRNO)
     llvmlibc_errno = err;
   return x;
 }
 
 static inline uint32_t as_uint32_bits(float x) {
   return *reinterpret_cast<uint32_t *>(&x);
 }
 
+static inline uint64_t as_uint64_bits(double x) {
+  return *reinterpret_cast<uint64_t *>(&x);
+}
+
 static inline float as_float(uint32_t x) {
   return *reinterpret_cast<float *>(&x);
 }
 
 static inline double as_double(uint64_t x) {
   return *reinterpret_cast<double *>(&x);
 }
 
-static inline constexpr float invalidf(float x) {
-  float y = (x - x) / (x - x);
-  return isnan(x) ? y : with_errnof(y, EDOM);
+template <typename T, EnableIfFloatOrDouble<T> = 0>
+static inline constexpr float invalid(T x) {
+  T y = (x - x) / (x - x);
+  return isnan(x) ? y : with_errno(y, EDOM);
 }
 
-static inline void force_eval_float(float x) { volatile float y UNUSED = x; }
+template <typename T> static inline void force_eval(T x) {
+  volatile T y UNUSED = x;
+}
+
+template <typename T> static inline T opt_barrier(T x) {
+  volatile T y = x;
+  return y;
+}
+
+template <typename T, EnableIfFloatOrDouble<T> = 0>
+T xflow(uint32_t sign, T y) {
+  // Underflow happens when two extremely small values are multiplied.
+  // Likewise, overflow happens when two large values are mulitplied.
+  y = opt_barrier(sign ? -y : y) * y;
+  return with_errno(y, ERANGE);
+}
+
+template <typename T, EnableIfFloatOrDouble<T> = 0> T overflow(uint32_t sign) {
+  return xflow(sign, XFlowValues<T>::overflow_value);
+}
+
+template <typename T, EnableIfFloatOrDouble<T> = 0> T underflow(uint32_t sign) {
+  return xflow(sign, XFlowValues<T>::underflow_value);
+}
+
+template <typename T, EnableIfFloatOrDouble<T> = 0>
+T may_underflow(uint32_t sign) {
+  return xflow(sign, XFlowValues<T>::may_underflow_value);
+}
+
+static inline uint32_t top12_bits(float x) { return as_uint32_bits(x) >> 20; }
+
+static inline uint32_t top12_bits(double x) { return as_uint64_bits(x) >> 52; }
 
 } // namespace __llvm_libc
 
 #endif // LLVM_LIBC_SRC_MATH_MATH_UTILS_H

libc/src/math/sincosf.cpp

void LLVM_LIBC_ENTRYPOINT(sincosf)(float y, float *sinp, float *cosp) {
   const sincos_t *p = &__sincosf_table[0];
 
   if (abstop12(y) < abstop12(pio4)) {
     double x2 = x * x;
 
     if (unlikely(abstop12(y) < abstop12(as_float(0x39800000)))) {
       if (unlikely(abstop12(y) < abstop12(as_float(0x800000))))
         // Force underflow for tiny y.
-        force_eval_float(x2);
+        force_eval<float>(x2);
       *sinp = y;
       *cosp = 1.0f;
       return;
     }
 
     sincosf_poly(x, x2, p, 0, sinp, cosp);
   } else if (abstop12(y) < abstop12(120.0f)) {
     x = reduce_fast(x, p, &n);
if (abstop12(y) < abstop12(pio4)) {
     sincosf_poly(x * s, x * x, p, n, sinp, cosp);
   } else {
     // Return NaN if Inf or NaN for both sin and cos.
     *sinp = *cosp = y - y;
 
     // Needed to set errno for +-Inf, the add is a hack to work
     // around a gcc register allocation issue: just passing y
     // affects code generation in the fast path.
-    invalidf(y + y);
+    invalid(y + y);
   }
 }
 
 } // namespace __llvm_libc

libc/src/math/sinf.cpp

float LLVM_LIBC_ENTRYPOINT(sinf)(float y) {
   const sincos_t *p = &__sincosf_table[0];
 
   if (abstop12(y) < abstop12(pio4)) {
     s = x * x;
 
     if (unlikely(abstop12(y) < abstop12(as_float(0x39800000)))) {
       if (unlikely(abstop12(y) < abstop12(as_float(0x800000))))
         // Force underflow for tiny y.
-        force_eval_float(s);
+        force_eval<float>(s);
       return y;
     }
 
     return sinf_poly(x, s, p, 0);
   } else if (likely(abstop12(y) < abstop12(120.0f))) {
     x = reduce_fast(x, p, &n);
 
     // Setup the signs for sin and cos.
if (abstop12(y) < abstop12(pio4)) {
     s = p->sign[(n + sign) & 3];
 
     if ((n + sign) & 2)
       p = &__sincosf_table[1];
 
     return sinf_poly(x * s, x * x, p, n);
   }
 
-  return invalidf(y);
+  return invalid(y);
 }
 
 } // namespace __llvm_libc

libc/test/src/math/CMakeLists.txt

SRCS
     sincosf_test.cpp
   HDRS
     sdcomp26094.h
   DEPENDS
     .float_utils
     libc.src.math.sincosf
     libc.utils.CPP.standalone_cpp
 )
+
+add_math_unittest(
+  expf_test
+  NEED_MPFR
+  SUITE
+    libc_math_unittests
+  SRCS
+    expf_test.cpp
+  DEPENDS
+    .float_utils
+    libc.include.math
+    libc.src.math.expf
+)

libc/test/src/math/expf_test.cpp

+//===-- Unittests for expf ------------------------------------------------===//
+//
+// 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 "include/math.h"
+#include "src/math/expf.h"
+#include "test/src/math/float.h"
+#include "utils/MPFRWrapper/MPFRUtils.h"
+#include "utils/UnitTest/Test.h"
+
+#include <stdint.h>
+
+using __llvm_libc::as_float;
+namespace mpfr = __llvm_libc::testing::mpfr;
+
+// 12 additional bits of precision over the base precision of a |float|
+// value.
+static constexpr mpfr::Tolerance tolerance{mpfr::Tolerance::floatPrecision, 12,
+                                           3 * 0x1000 / 4};
+
+TEST(ExpfTest, InFloatRange) {
+  constexpr uint32_t count = 1000000;
+  constexpr uint32_t step = UINT32_MAX / count;
+  for (uint32_t i = 0, v = 0; i <= count; ++i, v += step) {
+    float x = as_float(v);
+    if (isnan(x) || isinf(x))
+      continue;
+    ASSERT_TRUE(mpfr::equalsExp(x, __llvm_libc::expf(x), tolerance));
+  }
+}

libc/utils/MPFRWrapper/MPFRUtils.h

struct Tolerance {
 
   unsigned int width; // Number of valid LSB bits in |value|.
 
   // The bits in the tolerance value. The tolerance value will be
   // sum(bits[width - i] * 2 ^ (- basePrecision - i)) for |i| in
   // range [1, width].
   uint32_t bits;
 };
 
 // Return true if |libcOutput| is within the tolerance |t| of the cos(x)
 // value as evaluated by MPFR.
 bool equalsCos(float x, float libcOutput, const Tolerance &t);
 
 // Return true if |libcOutput| is within the tolerance |t| of the sin(x)
 // value as evaluated by MPFR.
 bool equalsSin(float x, float libcOutput, const Tolerance &t);
 
+// Return true if |libcOutput| is within the tolerance |t| of the sin(x)
+// value as evaluated by MPFR.
+bool equalsExp(float x, float libcOutput, const Tolerance &t);
+
 } // namespace mpfr
 } // namespace testing
 } // namespace __llvm_libc
 
 #endif // LLVM_LIBC_UTILS_TESTUTILS_MPFRUTILS_H

libc/utils/MPFRWrapper/MPFRUtils.cpp

 //===-- Utils which wrap MPFR ---------------------------------------------===//
 //
 // 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 "MPFRUtils.h"
 
 #include <iostream>
+
+#include <math.h>
 #include <mpfr.h>
 
 namespace __llvm_libc {
 namespace testing {
 namespace mpfr {
 
 class MPFRNumber {
   // A precision value which allows sufficiently large additional
public:
   }
 
   static MPFRNumber sin(float x) {
     MPFRNumber result;
     MPFRNumber mpfrX(x);
     mpfr_sin(result.value, mpfrX.value, MPFR_RNDN);
     return result;
   }
+
+  static MPFRNumber exp(float x) {
+    MPFRNumber result;
+    MPFRNumber mpfrX(x);
+    mpfr_exp(result.value, mpfrX.value, MPFR_RNDN);
+    return result;
+  }
 };
 
 bool equalsCos(float input, float libcOutput, const Tolerance &t) {
   MPFRNumber mpfrResult = MPFRNumber::cos(input);
   MPFRNumber libcResult(libcOutput);
   return mpfrResult.isEqual(libcResult, t);
 }
 
 bool equalsSin(float input, float libcOutput, const Tolerance &t) {
   MPFRNumber mpfrResult = MPFRNumber::sin(input);
   MPFRNumber libcResult(libcOutput);
   return mpfrResult.isEqual(libcResult, t);
 }
 
+bool equalsExp(float input, float libcOutput, const Tolerance &t) {
+  MPFRNumber mpfrResult = MPFRNumber::exp(input);
+  MPFRNumber libcResult(libcOutput);
+  return mpfrResult.isEqual(libcResult, t);
+}
+
 } // namespace mpfr
 } // namespace testing
 } // namespace __llvm_libc