diff --git a/library/compiler-builtins/libm/crates/libm-macros/src/shared.rs b/library/compiler-builtins/libm/crates/libm-macros/src/shared.rs
index d17bc6ffccb7..e7d3d18d9a1b 100644
--- a/library/compiler-builtins/libm/crates/libm-macros/src/shared.rs
+++ b/library/compiler-builtins/libm/crates/libm-macros/src/shared.rs
@@ -9,7 +9,7 @@ const ALL_OPERATIONS_NESTED: &[(FloatTy, Signature, Option<Signature>, &[&str])]
         FloatTy::F16,
         Signature { args: &[Ty::F16], returns: &[Ty::F16] },
         None,
-        &["fabsf16", "sqrtf16", "truncf16"],
+        &["ceilf16", "fabsf16", "sqrtf16", "truncf16"],
     ),
     (
         // `fn(f32) -> f32`
@@ -40,7 +40,7 @@ const ALL_OPERATIONS_NESTED: &[(FloatTy, Signature, Option<Signature>, &[&str])]
         FloatTy::F128,
         Signature { args: &[Ty::F128], returns: &[Ty::F128] },
         None,
-        &["fabsf128", "sqrtf128", "truncf128"],
+        &["ceilf128", "fabsf128", "sqrtf128", "truncf128"],
     ),
     (
         // `(f16, f16) -> f16`
diff --git a/library/compiler-builtins/libm/crates/libm-test/benches/icount.rs b/library/compiler-builtins/libm/crates/libm-test/benches/icount.rs
index c8451f88cd38..84be3d5245b6 100644
--- a/library/compiler-builtins/libm/crates/libm-test/benches/icount.rs
+++ b/library/compiler-builtins/libm/crates/libm-test/benches/icount.rs
@@ -69,6 +69,8 @@ main!(
     icount_bench_cbrt_group,
     icount_bench_cbrtf_group,
     icount_bench_ceil_group,
+    icount_bench_ceilf128_group,
+    icount_bench_ceilf16_group,
     icount_bench_ceilf_group,
     icount_bench_copysign_group,
     icount_bench_copysignf128_group,
diff --git a/library/compiler-builtins/libm/crates/libm-test/benches/random.rs b/library/compiler-builtins/libm/crates/libm-test/benches/random.rs
index 026841202532..511e26d91267 100644
--- a/library/compiler-builtins/libm/crates/libm-test/benches/random.rs
+++ b/library/compiler-builtins/libm/crates/libm-test/benches/random.rs
@@ -117,7 +117,9 @@ libm_macros::for_each_function! {
         exp10 | exp10f | exp2 | exp2f => (true, Some(musl_math_sys::MACRO_FN_NAME)),
 
         // Musl does not provide `f16` and `f128` functions
-        copysignf128
+        ceilf128
+        | ceilf16
+        | copysignf128
         | copysignf16
         | fabsf128
         | fabsf16
diff --git a/library/compiler-builtins/libm/crates/libm-test/src/f8_impl.rs b/library/compiler-builtins/libm/crates/libm-test/src/f8_impl.rs
index 299553d202ad..96b78392453f 100644
--- a/library/compiler-builtins/libm/crates/libm-test/src/f8_impl.rs
+++ b/library/compiler-builtins/libm/crates/libm-test/src/f8_impl.rs
@@ -30,6 +30,8 @@ impl Float for f8 {
     const INFINITY: Self = Self(0b0_1111_000);
     const NEG_INFINITY: Self = Self(0b1_1111_000);
     const NAN: Self = Self(0b0_1111_100);
+    // FIXME: incorrect values
+    const EPSILON: Self = Self::ZERO;
     const PI: Self = Self::ZERO;
     const NEG_PI: Self = Self::ZERO;
     const FRAC_PI_2: Self = Self::ZERO;
diff --git a/library/compiler-builtins/libm/crates/libm-test/src/mpfloat.rs b/library/compiler-builtins/libm/crates/libm-test/src/mpfloat.rs
index 9d95356d39b3..bbd19dbb07ba 100644
--- a/library/compiler-builtins/libm/crates/libm-test/src/mpfloat.rs
+++ b/library/compiler-builtins/libm/crates/libm-test/src/mpfloat.rs
@@ -137,6 +137,8 @@ libm_macros::for_each_function! {
         // Most of these need a manual implementation
         ceil,
         ceilf,
+        ceilf128,
+        ceilf16,
         copysign,
         copysignf,
         copysignf128,
@@ -237,12 +239,14 @@ impl_no_round! {
 #[cfg(f16_enabled)]
 impl_no_round! {
     fabsf16 => abs_mut;
+    ceilf16 => ceil_mut;
     truncf16 => trunc_mut;
 }
 
 #[cfg(f128_enabled)]
 impl_no_round! {
     fabsf128 => abs_mut;
+    ceilf128 => ceil_mut;
     truncf128 => trunc_mut;
 }
 
diff --git a/library/compiler-builtins/libm/crates/libm-test/tests/compare_built_musl.rs b/library/compiler-builtins/libm/crates/libm-test/tests/compare_built_musl.rs
index 24703f273168..e13acf3de216 100644
--- a/library/compiler-builtins/libm/crates/libm-test/tests/compare_built_musl.rs
+++ b/library/compiler-builtins/libm/crates/libm-test/tests/compare_built_musl.rs
@@ -79,6 +79,8 @@ libm_macros::for_each_function! {
         ynf,
 
         // Not provided by musl
+        ceilf128,
+        ceilf16,
         copysignf128,
         copysignf16,
         fabsf128,
diff --git a/library/compiler-builtins/libm/crates/util/src/main.rs b/library/compiler-builtins/libm/crates/util/src/main.rs
index cd68d9afd1a2..810919339ec6 100644
--- a/library/compiler-builtins/libm/crates/util/src/main.rs
+++ b/library/compiler-builtins/libm/crates/util/src/main.rs
@@ -84,7 +84,9 @@ fn do_eval(basis: &str, op: &str, inputs: &[&str]) {
         emit_types: [CFn, RustFn, RustArgs],
         extra: (basis, op, inputs),
         fn_extra: match MACRO_FN_NAME {
-            copysignf128
+            ceilf128
+            | ceilf16
+            | copysignf128
             | copysignf16
             | fabsf128
             | fabsf16
diff --git a/library/compiler-builtins/libm/etc/function-definitions.json b/library/compiler-builtins/libm/etc/function-definitions.json
index 2d0af3bcff02..c75152f63019 100644
--- a/library/compiler-builtins/libm/etc/function-definitions.json
+++ b/library/compiler-builtins/libm/etc/function-definitions.json
@@ -109,17 +109,33 @@
             "src/libm_helper.rs",
             "src/math/arch/i586.rs",
             "src/math/arch/wasm32.rs",
-            "src/math/ceil.rs"
+            "src/math/ceil.rs",
+            "src/math/generic/ceil.rs"
         ],
         "type": "f64"
     },
     "ceilf": {
         "sources": [
             "src/math/arch/wasm32.rs",
-            "src/math/ceilf.rs"
+            "src/math/ceilf.rs",
+            "src/math/generic/ceil.rs"
         ],
         "type": "f32"
     },
+    "ceilf128": {
+        "sources": [
+            "src/math/ceilf128.rs",
+            "src/math/generic/ceil.rs"
+        ],
+        "type": "f128"
+    },
+    "ceilf16": {
+        "sources": [
+            "src/math/ceilf16.rs",
+            "src/math/generic/ceil.rs"
+        ],
+        "type": "f16"
+    },
     "copysign": {
         "sources": [
             "src/libm_helper.rs",
diff --git a/library/compiler-builtins/libm/etc/function-list.txt b/library/compiler-builtins/libm/etc/function-list.txt
index 47c34ab90521..337e7e434e5f 100644
--- a/library/compiler-builtins/libm/etc/function-list.txt
+++ b/library/compiler-builtins/libm/etc/function-list.txt
@@ -17,6 +17,8 @@ cbrt
 cbrtf
 ceil
 ceilf
+ceilf128
+ceilf16
 copysign
 copysignf
 copysignf128
diff --git a/library/compiler-builtins/libm/src/math/ceil.rs b/library/compiler-builtins/libm/src/math/ceil.rs
index 398bfee47a21..535f434ac243 100644
--- a/library/compiler-builtins/libm/src/math/ceil.rs
+++ b/library/compiler-builtins/libm/src/math/ceil.rs
@@ -1,8 +1,3 @@
-#![allow(unreachable_code)]
-use core::f64;
-
-const TOINT: f64 = 1. / f64::EPSILON;
-
 /// Ceil (f64)
 ///
 /// Finds the nearest integer greater than or equal to `x`.
@@ -15,40 +10,5 @@ pub fn ceil(x: f64) -> f64 {
         args: x,
     }
 
-    let u: u64 = x.to_bits();
-    let e: i64 = ((u >> 52) & 0x7ff) as i64;
-    let y: f64;
-
-    if e >= 0x3ff + 52 || x == 0. {
-        return x;
-    }
-    // y = int(x) - x, where int(x) is an integer neighbor of x
-    y = if (u >> 63) != 0 { x - TOINT + TOINT - x } else { x + TOINT - TOINT - x };
-    // special case because of non-nearest rounding modes
-    if e < 0x3ff {
-        force_eval!(y);
-        return if (u >> 63) != 0 { -0. } else { 1. };
-    }
-    if y < 0. { x + y + 1. } else { x + y }
-}
-
-#[cfg(test)]
-mod tests {
-    use super::*;
-
-    #[test]
-    fn sanity_check() {
-        assert_eq!(ceil(1.1), 2.0);
-        assert_eq!(ceil(2.9), 3.0);
-    }
-
-    /// The spec: https://en.cppreference.com/w/cpp/numeric/math/ceil
-    #[test]
-    fn spec_tests() {
-        // Not Asserted: that the current rounding mode has no effect.
-        assert!(ceil(f64::NAN).is_nan());
-        for f in [0.0, -0.0, f64::INFINITY, f64::NEG_INFINITY].iter().copied() {
-            assert_eq!(ceil(f), f);
-        }
-    }
+    super::generic::ceil(x)
 }
diff --git a/library/compiler-builtins/libm/src/math/ceilf.rs b/library/compiler-builtins/libm/src/math/ceilf.rs
index 9e8e78e3e532..66d44189c921 100644
--- a/library/compiler-builtins/libm/src/math/ceilf.rs
+++ b/library/compiler-builtins/libm/src/math/ceilf.rs
@@ -1,5 +1,3 @@
-use core::f32;
-
 /// Ceil (f32)
 ///
 /// Finds the nearest integer greater than or equal to `x`.
@@ -11,52 +9,5 @@ pub fn ceilf(x: f32) -> f32 {
         args: x,
     }
 
-    let mut ui = x.to_bits();
-    let e = (((ui >> 23) & 0xff).wrapping_sub(0x7f)) as i32;
-
-    if e >= 23 {
-        return x;
-    }
-    if e >= 0 {
-        let m = 0x007fffff >> e;
-        if (ui & m) == 0 {
-            return x;
-        }
-        force_eval!(x + f32::from_bits(0x7b800000));
-        if ui >> 31 == 0 {
-            ui += m;
-        }
-        ui &= !m;
-    } else {
-        force_eval!(x + f32::from_bits(0x7b800000));
-        if ui >> 31 != 0 {
-            return -0.0;
-        } else if ui << 1 != 0 {
-            return 1.0;
-        }
-    }
-    f32::from_bits(ui)
-}
-
-// PowerPC tests are failing on LLVM 13: https://github.com/rust-lang/rust/issues/88520
-#[cfg(not(target_arch = "powerpc64"))]
-#[cfg(test)]
-mod tests {
-    use super::*;
-
-    #[test]
-    fn sanity_check() {
-        assert_eq!(ceilf(1.1), 2.0);
-        assert_eq!(ceilf(2.9), 3.0);
-    }
-
-    /// The spec: https://en.cppreference.com/w/cpp/numeric/math/ceil
-    #[test]
-    fn spec_tests() {
-        // Not Asserted: that the current rounding mode has no effect.
-        assert!(ceilf(f32::NAN).is_nan());
-        for f in [0.0, -0.0, f32::INFINITY, f32::NEG_INFINITY].iter().copied() {
-            assert_eq!(ceilf(f), f);
-        }
-    }
+    super::generic::ceil(x)
 }
diff --git a/library/compiler-builtins/libm/src/math/ceilf128.rs b/library/compiler-builtins/libm/src/math/ceilf128.rs
new file mode 100644
index 000000000000..89980858eb5b
--- /dev/null
+++ b/library/compiler-builtins/libm/src/math/ceilf128.rs
@@ -0,0 +1,7 @@
+/// Ceil (f128)
+///
+/// Finds the nearest integer greater than or equal to `x`.
+#[cfg_attr(all(test, assert_no_panic), no_panic::no_panic)]
+pub fn ceilf128(x: f128) -> f128 {
+    super::generic::ceil(x)
+}
diff --git a/library/compiler-builtins/libm/src/math/ceilf16.rs b/library/compiler-builtins/libm/src/math/ceilf16.rs
new file mode 100644
index 000000000000..2af67eff0d2a
--- /dev/null
+++ b/library/compiler-builtins/libm/src/math/ceilf16.rs
@@ -0,0 +1,7 @@
+/// Ceil (f16)
+///
+/// Finds the nearest integer greater than or equal to `x`.
+#[cfg_attr(all(test, assert_no_panic), no_panic::no_panic)]
+pub fn ceilf16(x: f16) -> f16 {
+    super::generic::ceil(x)
+}
diff --git a/library/compiler-builtins/libm/src/math/generic/ceil.rs b/library/compiler-builtins/libm/src/math/generic/ceil.rs
new file mode 100644
index 000000000000..34261faf782f
--- /dev/null
+++ b/library/compiler-builtins/libm/src/math/generic/ceil.rs
@@ -0,0 +1,101 @@
+/* SPDX-License-Identifier: MIT */
+/* origin: musl src/math/ceilf.c */
+
+//! Generic `ceil` algorithm.
+//!
+//! Note that this uses the algorithm from musl's `ceilf` rather than `ceil` or `ceill` because
+//! performance seems to be better (based on icount) and it does not seem to experience rounding
+//! errors on i386.
+
+use super::super::{Float, Int, IntTy, MinInt};
+
+pub fn ceil<F: Float>(x: F) -> F {
+    let zero = IntTy::<F>::ZERO;
+
+    let mut ix = x.to_bits();
+    let e = x.exp_unbiased();
+
+    // If the represented value has no fractional part, no truncation is needed.
+    if e >= F::SIG_BITS as i32 {
+        return x;
+    }
+
+    if e >= 0 {
+        // |x| >= 1.0
+
+        let m = F::SIG_MASK >> e.unsigned();
+        if (ix & m) == zero {
+            // Portion to be masked is already zero; no adjustment needed.
+            return x;
+        }
+
+        // Otherwise, raise an inexact exception.
+        force_eval!(x + F::MAX);
+        if x.is_sign_positive() {
+            ix += m;
+        }
+        ix &= !m;
+    } else {
+        // |x| < 1.0, raise an inexact exception since truncation will happen (unless x == 0).
+        force_eval!(x + F::MAX);
+
+        if x.is_sign_negative() {
+            // -1.0 < x <= -0.0; rounding up goes toward -0.0.
+            return F::NEG_ZERO;
+        } else if ix << 1 != zero {
+            // 0.0 < x < 1.0; rounding up goes toward +1.0.
+            return F::ONE;
+        }
+    }
+
+    F::from_bits(ix)
+}
+
+#[cfg(test)]
+mod tests {
+    use super::*;
+
+    /// Test against https://en.cppreference.com/w/cpp/numeric/math/ceil
+    fn spec_test<F: Float>() {
+        // Not Asserted: that the current rounding mode has no effect.
+        for f in [F::ZERO, F::NEG_ZERO, F::INFINITY, F::NEG_INFINITY].iter().copied() {
+            assert_biteq!(ceil(f), f);
+        }
+    }
+
+    /* Skipping f16 / f128 "sanity_check"s due to rejected literal lexing at MSRV */
+
+    #[test]
+    #[cfg(f16_enabled)]
+    fn spec_tests_f16() {
+        spec_test::<f16>();
+    }
+
+    #[test]
+    fn sanity_check_f32() {
+        assert_eq!(ceil(1.1f32), 2.0);
+        assert_eq!(ceil(2.9f32), 3.0);
+    }
+
+    #[test]
+    fn spec_tests_f32() {
+        spec_test::<f32>();
+    }
+
+    #[test]
+    fn sanity_check_f64() {
+        assert_eq!(ceil(1.1f64), 2.0);
+        assert_eq!(ceil(2.9f64), 3.0);
+    }
+
+    #[test]
+    fn spec_tests_f64() {
+        spec_test::<f64>();
+    }
+
+    #[test]
+    #[cfg(f128_enabled)]
+    fn spec_tests_f128() {
+        spec_test::<f128>();
+    }
+}
diff --git a/library/compiler-builtins/libm/src/math/generic/mod.rs b/library/compiler-builtins/libm/src/math/generic/mod.rs
index 3b5a2c3effd7..f8bb9fa6adb1 100644
--- a/library/compiler-builtins/libm/src/math/generic/mod.rs
+++ b/library/compiler-builtins/libm/src/math/generic/mod.rs
@@ -1,9 +1,11 @@
+mod ceil;
 mod copysign;
 mod fabs;
 mod fdim;
 mod sqrt;
 mod trunc;
 
+pub use ceil::ceil;
 pub use copysign::copysign;
 pub use fabs::fabs;
 pub use fdim::fdim;
diff --git a/library/compiler-builtins/libm/src/math/generic/sqrt.rs b/library/compiler-builtins/libm/src/math/generic/sqrt.rs
index c892f9997241..22ee93f57a3c 100644
--- a/library/compiler-builtins/libm/src/math/generic/sqrt.rs
+++ b/library/compiler-builtins/libm/src/math/generic/sqrt.rs
@@ -96,7 +96,7 @@ where
         ix = scaled.to_bits();
         match top {
             Exp::Shifted(ref mut v) => {
-                *v = scaled.exp().unsigned();
+                *v = scaled.exp();
                 *v = (*v).wrapping_sub(F::SIG_BITS);
             }
             Exp::NoShift(()) => {
diff --git a/library/compiler-builtins/libm/src/math/mod.rs b/library/compiler-builtins/libm/src/math/mod.rs
index 3684025a6c10..5228e78b7f33 100644
--- a/library/compiler-builtins/libm/src/math/mod.rs
+++ b/library/compiler-builtins/libm/src/math/mod.rs
@@ -341,12 +341,14 @@ pub use self::truncf::truncf;
 
 cfg_if! {
     if #[cfg(f16_enabled)] {
+        mod ceilf16;
         mod copysignf16;
         mod fabsf16;
         mod fdimf16;
         mod sqrtf16;
         mod truncf16;
 
+        pub use self::ceilf16::ceilf16;
         pub use self::copysignf16::copysignf16;
         pub use self::fabsf16::fabsf16;
         pub use self::fdimf16::fdimf16;
@@ -357,12 +359,14 @@ cfg_if! {
 
 cfg_if! {
     if #[cfg(f128_enabled)] {
+        mod ceilf128;
         mod copysignf128;
         mod fabsf128;
         mod fdimf128;
         mod sqrtf128;
         mod truncf128;
 
+        pub use self::ceilf128::ceilf128;
         pub use self::copysignf128::copysignf128;
         pub use self::fabsf128::fabsf128;
         pub use self::fdimf128::fdimf128;
diff --git a/library/compiler-builtins/libm/src/math/support/float_traits.rs b/library/compiler-builtins/libm/src/math/support/float_traits.rs
index 647f4f5e251b..1abb7c4ded84 100644
--- a/library/compiler-builtins/libm/src/math/support/float_traits.rs
+++ b/library/compiler-builtins/libm/src/math/support/float_traits.rs
@@ -1,6 +1,6 @@
 use core::{fmt, mem, ops};
 
-use super::int_traits::{CastFrom, CastInto, Int, MinInt};
+use super::int_traits::{CastFrom, Int, MinInt};
 
 /// Trait for some basic operations on floats
 #[allow(dead_code)]
@@ -34,6 +34,7 @@ pub trait Float:
     const NAN: Self;
     const MAX: Self;
     const MIN: Self;
+    const EPSILON: Self;
     const PI: Self;
     const NEG_PI: Self;
     const FRAC_PI_2: Self;
@@ -107,13 +108,13 @@ pub trait Float:
     }
 
     /// Returns the exponent, not adjusting for bias, not accounting for subnormals or zero.
-    fn exp(self) -> i32 {
-        ((self.to_bits() & Self::EXP_MASK) >> Self::SIG_BITS).cast()
+    fn exp(self) -> u32 {
+        u32::cast_from(self.to_bits() >> Self::SIG_BITS) & Self::EXP_MAX
     }
 
     /// Extract the exponent and adjust it for bias, not accounting for subnormals or zero.
     fn exp_unbiased(self) -> i32 {
-        self.exp() - (Self::EXP_BIAS as i32)
+        self.exp().signed() - (Self::EXP_BIAS as i32)
     }
 
     /// Returns the significand with no implicit bit (or the "fractional" part)
@@ -180,6 +181,7 @@ macro_rules! float_impl {
             const MAX: Self = -Self::MIN;
             // Sign bit set, saturated mantissa, saturated exponent with last bit zeroed
             const MIN: Self = $from_bits(Self::Int::MAX & !(1 << Self::SIG_BITS));
+            const EPSILON: Self = <$ty>::EPSILON;
 
             const PI: Self = core::$ty::consts::PI;
             const NEG_PI: Self = -Self::PI;