diff --git a/libc/src/__support/FPUtil/CMakeLists.txt b/libc/src/__support/FPUtil/CMakeLists.txt
--- a/libc/src/__support/FPUtil/CMakeLists.txt
+++ b/libc/src/__support/FPUtil/CMakeLists.txt
@@ -13,6 +13,7 @@
     NormalFloat.h
     PlatformDefs.h
     UInt.h
+    builtin_wrappers.h
   DEPENDS
     libc.include.math
     libc.include.errno
diff --git a/libc/src/__support/FPUtil/FPBits.h b/libc/src/__support/FPUtil/FPBits.h
--- a/libc/src/__support/FPUtil/FPBits.h
+++ b/libc/src/__support/FPUtil/FPBits.h
@@ -13,6 +13,8 @@
 
 #include "src/__support/CPP/Bit.h"
 #include "src/__support/CPP/TypeTraits.h"
+#include "src/__support/FPUtil/builtin_wrappers.h"
+#include "src/__support/common.h"
 
 #include "FloatProperties.h"
 #include <stdint.h>
@@ -28,32 +30,6 @@
   static constexpr unsigned VALUE = FloatProperties<T>::EXPONENT_WIDTH;
 };
 
-// The following overloads are matched based on what is accepted by
-// __builtin_clz/ctz* rather than using the exactly-sized aliases from stdint.h.
-// This way, we can avoid making any assumptions about integer sizes and let the
-// compiler match for us.
-template <typename T> static inline int clz(T val);
-template <> inline int clz<unsigned int>(unsigned int val) {
-  return __builtin_clz(val);
-}
-template <> inline int clz<unsigned long int>(unsigned long int val) {
-  return __builtin_clzl(val);
-}
-template <> inline int clz<unsigned long long int>(unsigned long long int val) {
-  return __builtin_clzll(val);
-}
-
-template <typename T> static inline int ctz(T val);
-template <> inline int ctz<unsigned int>(unsigned int val) {
-  return __builtin_ctz(val);
-}
-template <> inline int ctz<unsigned long int>(unsigned long int val) {
-  return __builtin_ctzl(val);
-}
-template <> inline int ctz<unsigned long long int>(unsigned long long int val) {
-  return __builtin_ctzll(val);
-}
-
 // A generic class to represent single precision, double precision, and quad
 // precision IEEE 754 floating point formats.
 // On most platforms, the 'float' type corresponds to single precision floating
@@ -83,6 +59,11 @@
 
   UIntType get_mantissa() const { return bits & FloatProp::MANTISSA_MASK; }
 
+  // The function return mantissa with implicit bit set for normal values.
+  constexpr UIntType get_explicit_mantissa() {
+    return (FloatProp::MANTISSA_MASK + 1) | (FloatProp::MANTISSA_MASK & bits);
+  }
+
   void set_unbiased_exponent(UIntType expVal) {
     expVal = (expVal << (FloatProp::MANTISSA_WIDTH)) & FloatProp::EXPONENT_MASK;
     bits &= ~(FloatProp::EXPONENT_MASK);
@@ -182,6 +163,27 @@
     bits.set_mantissa(v);
     return T(bits);
   }
+
+  // The function convert integer number and unbiased exponent to proper Float.
+  // Result = number * 2^(ep+1)
+  // Be careful! The function add to +1 to ep for seamless normalized to
+  // denormalized transition.
+  inline static constexpr FPBits<T> make_value(UIntType number, int ep) {
+    FPBits<T> result;
+    // offset: +1 for sign, but -1 for implicit first bit
+    int lz = fputil::clz(number) - FloatProp::EXPONENT_WIDTH;
+    number <<= lz;
+    ep -= lz;
+
+    if (likely(ep >= 0)) {
+      // Implicit number bit will be removed by mask
+      result.set_mantissa(number);
+      result.set_unbiased_exponent(ep + 1);
+    } else {
+      result.set_mantissa(number >> -ep);
+    }
+    return result;
+  }
 };
 
 } // namespace fputil
diff --git a/libc/src/__support/FPUtil/builtin_wrappers.h b/libc/src/__support/FPUtil/builtin_wrappers.h
new file mode 100644
--- /dev/null
+++ b/libc/src/__support/FPUtil/builtin_wrappers.h
@@ -0,0 +1,61 @@
+//===--Convenient template for builtins -------------------------*- C++ -*-===//
+//             (Count Lead Zeroes) and (Count Trailing Zeros)
+//
+// 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_SUPPORT_FPUTIL_CLZ_H
+#define LLVM_LIBC_SRC_SUPPORT_FPUTIL_CLZ_H
+
+namespace __llvm_libc {
+namespace fputil {
+
+// The following overloads are matched based on what is accepted by
+// __builtin_clz/ctz* rather than using the exactly-sized aliases from stdint.h.
+// This way, we can avoid making any assumptions about integer sizes and let the
+// compiler match for us.
+template <typename T> static inline int clz(T val);
+template <> inline int clz<unsigned int>(unsigned int val) {
+  return __builtin_clz(val);
+}
+template <> inline int clz<unsigned long int>(unsigned long int val) {
+  return __builtin_clzl(val);
+}
+template <> inline int clz<unsigned long long int>(unsigned long long int val) {
+  return __builtin_clzll(val);
+}
+
+template <typename T> static inline int ctz(T val);
+template <> inline int ctz<unsigned int>(unsigned int val) {
+  return __builtin_ctz(val);
+}
+template <> inline int ctz<unsigned long int>(unsigned long int val) {
+  return __builtin_ctzl(val);
+}
+template <> inline int ctz<unsigned long long int>(unsigned long long int val) {
+  return __builtin_ctzll(val);
+}
+
+template <typename T> static inline bool isnan(T val) {
+  return __builtin_isnan(val);
+}
+
+template <typename T> static inline bool isinf(T val) {
+  return __builtin_isinf(val);
+}
+
+template <typename T> static inline bool isfinite(T val) {
+  return __builtin_isfinite(val);
+}
+
+inline float quiet_NaN(float) { return __builtin_nanf(""); }
+inline double quiet_NaN(double) { return __builtin_nan(""); }
+inline long double quiet_NaN(long double) { return __builtin_nanl(""); }
+
+} // namespace fputil
+} // namespace __llvm_libc
+
+#endif // LLVM_LIBC_SRC_SUPPORT_FPUTIL_CLZ_H
diff --git a/libc/src/__support/FPUtil/generic/CMakeLists.txt b/libc/src/__support/FPUtil/generic/CMakeLists.txt
--- a/libc/src/__support/FPUtil/generic/CMakeLists.txt
+++ b/libc/src/__support/FPUtil/generic/CMakeLists.txt
@@ -15,4 +15,4 @@
   fmod
   HDRS
     FMod.h
-)
\ No newline at end of file
+)
diff --git a/libc/src/__support/FPUtil/generic/FMod.h b/libc/src/__support/FPUtil/generic/FMod.h
--- a/libc/src/__support/FPUtil/generic/FMod.h
+++ b/libc/src/__support/FPUtil/generic/FMod.h
@@ -9,17 +9,13 @@
 #ifndef LLVM_LIBC_SRC_SUPPORT_FPUTIL_GENERIC_FMOD_H
 #define LLVM_LIBC_SRC_SUPPORT_FPUTIL_GENERIC_FMOD_H
 
+#include "src/__support/CPP/TypeTraits.h"
 #include "src/__support/FPUtil/FPBits.h"
+#include "src/__support/FPUtil/builtin_wrappers.h"
 #include "src/__support/common.h"
 #include "src/math/generic/math_utils.h"
 
-#include <algorithm>
 #include <cerrno>
-#include <cmath>
-#include <limits>
-#include <optional>
-#include <stdint.h>
-#include <type_traits>
 
 namespace __llvm_libc {
 namespace fputil {
@@ -75,44 +71,42 @@
 template <typename T> class FModGCCInspiredWrapper {
 
   static_assert(cpp::IsFloatingPointType<T>::Value,
-                "FPBits instantiated with invalid type.");
+                "FModGCCInspiredWrapper instantiated with invalid type.");
 
 public:
-  static std::optional<T> PreCheck(T x, T y) {
-    if (unlikely(y == 0 || std::isnan(y) || !std::isfinite(x))) {
+  static bool PreCheck(T x, T y, T &out) {
+    if (unlikely(y == 0 || fputil::isnan(y) || !fputil::isfinite(x))) {
       force_eval<T>((x * y) / (y * x));
       // y=0,or x not finite or y is NaN
-      if (std::isnan(x) || std::isnan(y))
-        return std::numeric_limits<T>::quiet_NaN();
-      else
-        return with_errno(std::numeric_limits<T>::quiet_NaN(), EDOM);
+      if (fputil::isnan(x) || fputil::isnan(y)) {
+        out = fputil::quiet_NaN(T(0));
+      } else {
+        out = with_errno(quiet_NaN(T(0)), EDOM);
+      }
+      return true;
     }
-    return std::nullopt;
+    return false;
   }
 };
 
 template <typename T> class FModFastMathWrapper {
 
   static_assert(cpp::IsFloatingPointType<T>::Value,
-                "FPBits instantiated with invalid type.");
+                "FModFastMathWrapper instantiated with invalid type.");
 
 public:
-  static std::optional<T> PreCheck(T, T) { return std::nullopt; }
+  static bool PreCheck(T, T, T &) { return false; }
 };
 
 template <typename T, class Wrapper = FModGCCInspiredWrapper<T>,
           bool InverseMultiplication = false>
 class FMod {
   static_assert(cpp::IsFloatingPointType<T>::Value,
-                "FPBits instantiated with invalid type.");
+                "FMod instantiated with invalid type.");
 
 private:
   using FPB = FPBits<T>;
   using intU_t = typename FPB::UIntType;
-  static constexpr intU_t set_normal_value(intU_t a) {
-    return (FPB::FloatProp::MANTISSA_MASK + 1) |
-           (FPB::FloatProp::MANTISSA_MASK & a);
-  }
 
   template <bool B = InverseMultiplication>
   typename std::enable_if<
@@ -158,23 +152,6 @@
     return hx;
   }
 
-  inline static constexpr FPB make_value(intU_t number, int ep) {
-    FPB result;
-    // offset: +1 for sign, but -1 for implicit first bit
-    int lz = fputil::clz(number) - FPB::FloatProp::EXPONENT_WIDTH;
-    number <<= lz;
-    ep -= lz;
-
-    if (likely(ep >= 0)) {
-      // Implicit number bit will be removed by mask
-      result.set_mantissa(number);
-      result.set_unbiased_exponent(ep + 1);
-    } else {
-      result.set_mantissa(number >> -ep);
-    }
-    return result;
-  }
-
   inline static constexpr FPB make_internal(FPB sx, FPB sy) {
 
     if (likely(sx.uintval() <= sy.uintval())) {
@@ -183,38 +160,37 @@
       return FPB::zero(); // |x|=|y| return 0.0
     }
 
-    intU_t hx = sx.uintval();
     int ix = sx.get_unbiased_exponent();
-
-    intU_t hy = sy.uintval();
     int iy = sy.get_unbiased_exponent();
 
     // Most common case where |y| is "very normal" and |x/y| < 2^EXPONENT_WIDTH
     if (likely(iy > int(FPB::FloatProp::MANTISSA_WIDTH) &&
                ix - iy <= int(FPB::FloatProp::EXPONENT_WIDTH))) {
-      hx = set_normal_value(hx);
-      hy = set_normal_value(hy);
+      intU_t hx = sx.get_explicit_mantissa();
+      intU_t hy = sy.get_explicit_mantissa();
       intU_t d = (ix == iy) ? (hx - hy) : (hx << (ix - iy)) % hy;
       if (d == 0)
         return FPB::zero();
       // iy - 1 because of "zero power" for number with power 1
-      return make_value(d, iy - 1);
+      return FPB::make_value(d, iy - 1);
     }
     /* Both subnormal special case. */
     if (unlikely(ix == 0 && iy == 0)) {
       FPB d;
-      d.set_mantissa(hx % hy);
+      d.set_mantissa(sx.uintval() % sy.uintval());
       return d;
     }
 
     // Note that hx is not subnormal by conditions above.
-    hx = set_normal_value(hx);
+    intU_t hx = sx.get_explicit_mantissa();
     ix--;
+
+    intU_t hy = sy.get_explicit_mantissa();
     int lzhy = FPB::FloatProp::EXPONENT_WIDTH;
     if (likely(iy > 0)) {
-      hy = set_normal_value(hy);
       iy--;
     } else {
+      hy = sy.get_mantissa();
       lzhy = clz(hy);
     }
 
@@ -226,7 +202,7 @@
 
     {
       // Shift hy right until the end or n = 0
-      int right_shift = std::min(n, tzhy);
+      int right_shift = n < tzhy ? n : tzhy;
       hy >>= right_shift;
       n -= right_shift;
       iy += right_shift;
@@ -234,7 +210,9 @@
 
     {
       // Shift hx left until the end or n = 0
-      int left_shift = std::min<int>(n, FPB::FloatProp::EXPONENT_WIDTH);
+      int left_shift = n < int(FPB::FloatProp::EXPONENT_WIDTH)
+                           ? n
+                           : FPB::FloatProp::EXPONENT_WIDTH;
       hx <<= left_shift;
       n -= left_shift;
     }
@@ -244,18 +222,17 @@
       return FPB::zero();
 
     if (n == 0)
-      return make_value(hx, iy);
+      return FPB::make_value(hx, iy);
 
     /* hx next can't be 0, because hx < hy, hy % 2 == 1 hx * 2^i % hy != 0 */
     hx = division_loop(n, hyltzeroes, hx, hy);
-    return make_value(hx, iy);
+    return FPB::make_value(hx, iy);
   }
 
 public:
   static inline T make(T x, T y) {
-    std::optional<T> out = Wrapper::PreCheck(x, y);
-    if (out.has_value())
-      return out.value();
+    if (T out; Wrapper::PreCheck(x, y, out))
+      return out;
     FPB sx(x), sy(y);
     bool sign = sx.get_sign();
     sx.set_sign(false);
diff --git a/libc/src/__support/FPUtil/generic/sqrt.h b/libc/src/__support/FPUtil/generic/sqrt.h
--- a/libc/src/__support/FPUtil/generic/sqrt.h
+++ b/libc/src/__support/FPUtil/generic/sqrt.h
@@ -15,6 +15,7 @@
 #include "src/__support/FPUtil/FEnvImpl.h"
 #include "src/__support/FPUtil/FPBits.h"
 #include "src/__support/FPUtil/PlatformDefs.h"
+#include "src/__support/FPUtil/builtin_wrappers.h"
 
 namespace __llvm_libc {
 namespace fputil {
diff --git a/libc/test/src/math/differential_testing/CMakeLists.txt b/libc/test/src/math/differential_testing/CMakeLists.txt
--- a/libc/test/src/math/differential_testing/CMakeLists.txt
+++ b/libc/test/src/math/differential_testing/CMakeLists.txt
@@ -487,4 +487,4 @@
     libc.src.math.fmod
   COMPILE_OPTIONS
     -fno-builtin
-)
\ No newline at end of file
+)
diff --git a/libc/test/src/math/exhaustive/FMod_test.cpp b/libc/test/src/math/exhaustive/FMod_test.cpp
--- a/libc/test/src/math/exhaustive/FMod_test.cpp
+++ b/libc/test/src/math/exhaustive/FMod_test.cpp
@@ -11,17 +11,17 @@
 #include "utils/UnitTest/Test.h"
 
 #include <array>
-// #include <cmath>
 #include <limits>
 
 namespace mpfr = __llvm_libc::testing::mpfr;
 
 template <typename T, bool InverseMultiplication>
 class LlvmLibcFModTest : public __llvm_libc::testing::Test {
-  static constexpr std::array<T, 10> test_bases = {
+  static constexpr std::array<T, 11> test_bases = {
       T(0.0),
       T(1.0),
       T(3.0),
+      T(27.0),
       T(11.0 / 8.0),
       T(2.764443),
       T(1.0) - std::numeric_limits<T>::epsilon(),
@@ -40,12 +40,12 @@
     int max2 = nl::max_exponent + 3;
     for (T by : test_bases) {
       for (int iy = min2; iy < max2; iy++) {
-        T y = by * std::pow<T>(2, iy);
+        T y = by * std::ldexp(2, iy);
         if (y == 0 || !std::isfinite(y))
           continue;
         for (T bx : test_bases) {
           for (int ix = min2; ix < max2; ix++) {
-            T x = bx * std::pow<T>(2, ix);
+            T x = bx * std::ldexp(2, ix);
             if (!std::isfinite(x))
               continue;
             T result = FMod::make(x, y);