Add math functions for `f16` and `f128`

This adds missing functions for math operations on the new float types.

Platform support is pretty spotty at this point, since even platforms
with generally good support can be missing math functions.
`std/build.rs` is updated to reflect this.

6 files changed, 3415 insertions, 92 deletions

diff --git a/library/std/src/f128.rs b/library/std/src/f128.rs
index a5b00d57cef..f6df6259137 100644
--- a/library/std/src/f128.rs
+++ b/library/std/src/f128.rs
@@ -12,25 +12,180 @@ pub use core::f128::consts;
 
 #[cfg(not(test))]
 use crate::intrinsics;
+#[cfg(not(test))]
+use crate::sys::cmath;
 
 #[cfg(not(test))]
 impl f128 {
-    /// Raises a number to an integer power.
+    /// Returns the largest integer less than or equal to `self`.
     ///
-    /// Using this function is generally faster than using `powf`.
-    /// It might have a different sequence of rounding operations than `powf`,
-    /// so the results are not guaranteed to agree.
+    /// This function always returns the precise result.
     ///
-    /// # Unspecified precision
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let f = 3.7_f128;
+    /// let g = 3.0_f128;
+    /// let h = -3.7_f128;
     ///
-    /// The precision of this function is non-deterministic. This means it varies by platform, Rust version, and
-    /// can even differ within the same execution from one invocation to the next.
+    /// assert_eq!(f.floor(), 3.0);
+    /// assert_eq!(g.floor(), 3.0);
+    /// assert_eq!(h.floor(), -4.0);
+    /// # }
+    /// ```
     #[inline]
     #[rustc_allow_incoherent_impl]
     #[unstable(feature = "f128", issue = "116909")]
     #[must_use = "method returns a new number and does not mutate the original value"]
-    pub fn powi(self, n: i32) -> f128 {
-        unsafe { intrinsics::powif128(self, n) }
+    pub fn floor(self) -> f128 {
+        unsafe { intrinsics::floorf128(self) }
+    }
+
+    /// Returns the smallest integer greater than or equal to `self`.
+    ///
+    /// This function always returns the precise result.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let f = 3.01_f128;
+    /// let g = 4.0_f128;
+    ///
+    /// assert_eq!(f.ceil(), 4.0);
+    /// assert_eq!(g.ceil(), 4.0);
+    /// # }
+    /// ```
+    #[inline]
+    #[doc(alias = "ceiling")]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn ceil(self) -> f128 {
+        unsafe { intrinsics::ceilf128(self) }
+    }
+
+    /// Returns the nearest integer to `self`. If a value is half-way between two
+    /// integers, round away from `0.0`.
+    ///
+    /// This function always returns the precise result.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let f = 3.3_f128;
+    /// let g = -3.3_f128;
+    /// let h = -3.7_f128;
+    /// let i = 3.5_f128;
+    /// let j = 4.5_f128;
+    ///
+    /// assert_eq!(f.round(), 3.0);
+    /// assert_eq!(g.round(), -3.0);
+    /// assert_eq!(h.round(), -4.0);
+    /// assert_eq!(i.round(), 4.0);
+    /// assert_eq!(j.round(), 5.0);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn round(self) -> f128 {
+        unsafe { intrinsics::roundf128(self) }
+    }
+
+    /// Returns the nearest integer to a number. Rounds half-way cases to the number
+    /// with an even least significant digit.
+    ///
+    /// This function always returns the precise result.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let f = 3.3_f128;
+    /// let g = -3.3_f128;
+    /// let h = 3.5_f128;
+    /// let i = 4.5_f128;
+    ///
+    /// assert_eq!(f.round_ties_even(), 3.0);
+    /// assert_eq!(g.round_ties_even(), -3.0);
+    /// assert_eq!(h.round_ties_even(), 4.0);
+    /// assert_eq!(i.round_ties_even(), 4.0);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn round_ties_even(self) -> f128 {
+        unsafe { intrinsics::rintf128(self) }
+    }
+
+    /// Returns the integer part of `self`.
+    /// This means that non-integer numbers are always truncated towards zero.
+    ///
+    /// This function always returns the precise result.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let f = 3.7_f128;
+    /// let g = 3.0_f128;
+    /// let h = -3.7_f128;
+    ///
+    /// assert_eq!(f.trunc(), 3.0);
+    /// assert_eq!(g.trunc(), 3.0);
+    /// assert_eq!(h.trunc(), -3.0);
+    /// # }
+    /// ```
+    #[inline]
+    #[doc(alias = "truncate")]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn trunc(self) -> f128 {
+        unsafe { intrinsics::truncf128(self) }
+    }
+
+    /// Returns the fractional part of `self`.
+    ///
+    /// This function always returns the precise result.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let x = 3.6_f128;
+    /// let y = -3.6_f128;
+    /// let abs_difference_x = (x.fract() - 0.6).abs();
+    /// let abs_difference_y = (y.fract() - (-0.6)).abs();
+    ///
+    /// assert!(abs_difference_x <= f128::EPSILON);
+    /// assert!(abs_difference_y <= f128::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn fract(self) -> f128 {
+        self - self.trunc()
     }
 
     /// Computes the absolute value of `self`.
@@ -41,7 +196,7 @@ impl f128 {
     ///
     /// ```
     /// #![feature(f128)]
-    /// # #[cfg(reliable_f128)] { // FIXME(f16_f128): reliable_f128
+    /// # #[cfg(reliable_f128)] {
     ///
     /// let x = 3.5_f128;
     /// let y = -3.5_f128;
@@ -61,4 +216,1129 @@ impl f128 {
         // We don't do this now because LLVM has lowering bugs for f128 math.
         Self::from_bits(self.to_bits() & !(1 << 127))
     }
+
+    /// Returns a number that represents the sign of `self`.
+    ///
+    /// - `1.0` if the number is positive, `+0.0` or `INFINITY`
+    /// - `-1.0` if the number is negative, `-0.0` or `NEG_INFINITY`
+    /// - NaN if the number is NaN
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let f = 3.5_f128;
+    ///
+    /// assert_eq!(f.signum(), 1.0);
+    /// assert_eq!(f128::NEG_INFINITY.signum(), -1.0);
+    ///
+    /// assert!(f128::NAN.signum().is_nan());
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn signum(self) -> f128 {
+        if self.is_nan() { Self::NAN } else { 1.0_f128.copysign(self) }
+    }
+
+    /// Returns a number composed of the magnitude of `self` and the sign of
+    /// `sign`.
+    ///
+    /// Equal to `self` if the sign of `self` and `sign` are the same, otherwise
+    /// equal to `-self`. If `self` is a NaN, then a NaN with the sign bit of
+    /// `sign` is returned. Note, however, that conserving the sign bit on NaN
+    /// across arithmetical operations is not generally guaranteed.
+    /// See [explanation of NaN as a special value](primitive@f128) for more info.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let f = 3.5_f128;
+    ///
+    /// assert_eq!(f.copysign(0.42), 3.5_f128);
+    /// assert_eq!(f.copysign(-0.42), -3.5_f128);
+    /// assert_eq!((-f).copysign(0.42), 3.5_f128);
+    /// assert_eq!((-f).copysign(-0.42), -3.5_f128);
+    ///
+    /// assert!(f128::NAN.copysign(1.0).is_nan());
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn copysign(self, sign: f128) -> f128 {
+        unsafe { intrinsics::copysignf128(self, sign) }
+    }
+
+    /// Fused multiply-add. Computes `(self * a) + b` with only one rounding
+    /// error, yielding a more accurate result than an unfused multiply-add.
+    ///
+    /// Using `mul_add` *may* be more performant than an unfused multiply-add if
+    /// the target architecture has a dedicated `fma` CPU instruction. However,
+    /// this is not always true, and will be heavily dependant on designing
+    /// algorithms with specific target hardware in mind.
+    ///
+    /// # Precision
+    ///
+    /// The result of this operation is guaranteed to be the rounded
+    /// infinite-precision result. It is specified by IEEE 754 as
+    /// `fusedMultiplyAdd` and guaranteed not to change.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let m = 10.0_f128;
+    /// let x = 4.0_f128;
+    /// let b = 60.0_f128;
+    ///
+    /// assert_eq!(m.mul_add(x, b), 100.0);
+    /// assert_eq!(m * x + b, 100.0);
+    ///
+    /// let one_plus_eps = 1.0_f128 + f128::EPSILON;
+    /// let one_minus_eps = 1.0_f128 - f128::EPSILON;
+    /// let minus_one = -1.0_f128;
+    ///
+    /// // The exact result (1 + eps) * (1 - eps) = 1 - eps * eps.
+    /// assert_eq!(one_plus_eps.mul_add(one_minus_eps, minus_one), -f128::EPSILON * f128::EPSILON);
+    /// // Different rounding with the non-fused multiply and add.
+    /// assert_eq!(one_plus_eps * one_minus_eps + minus_one, 0.0);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn mul_add(self, a: f128, b: f128) -> f128 {
+        unsafe { intrinsics::fmaf128(self, a, b) }
+    }
+
+    /// Calculates Euclidean division, the matching method for `rem_euclid`.
+    ///
+    /// This computes the integer `n` such that
+    /// `self = n * rhs + self.rem_euclid(rhs)`.
+    /// In other words, the result is `self / rhs` rounded to the integer `n`
+    /// such that `self >= n * rhs`.
+    ///
+    /// # Precision
+    ///
+    /// The result of this operation is guaranteed to be the rounded
+    /// infinite-precision result.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let a: f128 = 7.0;
+    /// let b = 4.0;
+    /// assert_eq!(a.div_euclid(b), 1.0); // 7.0 > 4.0 * 1.0
+    /// assert_eq!((-a).div_euclid(b), -2.0); // -7.0 >= 4.0 * -2.0
+    /// assert_eq!(a.div_euclid(-b), -1.0); // 7.0 >= -4.0 * -1.0
+    /// assert_eq!((-a).div_euclid(-b), 2.0); // -7.0 >= -4.0 * 2.0
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn div_euclid(self, rhs: f128) -> f128 {
+        let q = (self / rhs).trunc();
+        if self % rhs < 0.0 {
+            return if rhs > 0.0 { q - 1.0 } else { q + 1.0 };
+        }
+        q
+    }
+
+    /// Calculates the least nonnegative remainder of `self (mod rhs)`.
+    ///
+    /// In particular, the return value `r` satisfies `0.0 <= r < rhs.abs()` in
+    /// most cases. However, due to a floating point round-off error it can
+    /// result in `r == rhs.abs()`, violating the mathematical definition, if
+    /// `self` is much smaller than `rhs.abs()` in magnitude and `self < 0.0`.
+    /// This result is not an element of the function's codomain, but it is the
+    /// closest floating point number in the real numbers and thus fulfills the
+    /// property `self == self.div_euclid(rhs) * rhs + self.rem_euclid(rhs)`
+    /// approximately.
+    ///
+    /// # Precision
+    ///
+    /// The result of this operation is guaranteed to be the rounded
+    /// infinite-precision result.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let a: f128 = 7.0;
+    /// let b = 4.0;
+    /// assert_eq!(a.rem_euclid(b), 3.0);
+    /// assert_eq!((-a).rem_euclid(b), 1.0);
+    /// assert_eq!(a.rem_euclid(-b), 3.0);
+    /// assert_eq!((-a).rem_euclid(-b), 1.0);
+    /// // limitation due to round-off error
+    /// assert!((-f128::EPSILON).rem_euclid(3.0) != 0.0);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[doc(alias = "modulo", alias = "mod")]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn rem_euclid(self, rhs: f128) -> f128 {
+        let r = self % rhs;
+        if r < 0.0 { r + rhs.abs() } else { r }
+    }
+
+    /// Raises a number to an integer power.
+    ///
+    /// Using this function is generally faster than using `powf`.
+    /// It might have a different sequence of rounding operations than `powf`,
+    /// so the results are not guaranteed to agree.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn powi(self, n: i32) -> f128 {
+        unsafe { intrinsics::powif128(self, n) }
+    }
+
+    /// Raises a number to a floating point power.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let x = 2.0_f128;
+    /// let abs_difference = (x.powf(2.0) - (x * x)).abs();
+    ///
+    /// assert!(abs_difference <= f128::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn powf(self, n: f128) -> f128 {
+        unsafe { intrinsics::powf128(self, n) }
+    }
+
+    /// Returns the square root of a number.
+    ///
+    /// Returns NaN if `self` is a negative number other than `-0.0`.
+    ///
+    /// # Precision
+    ///
+    /// The result of this operation is guaranteed to be the rounded
+    /// infinite-precision result. It is specified by IEEE 754 as `squareRoot`
+    /// and guaranteed not to change.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let positive = 4.0_f128;
+    /// let negative = -4.0_f128;
+    /// let negative_zero = -0.0_f128;
+    ///
+    /// assert_eq!(positive.sqrt(), 2.0);
+    /// assert!(negative.sqrt().is_nan());
+    /// assert!(negative_zero.sqrt() == negative_zero);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn sqrt(self) -> f128 {
+        unsafe { intrinsics::sqrtf128(self) }
+    }
+
+    /// Returns `e^(self)`, (the exponential function).
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let one = 1.0f128;
+    /// // e^1
+    /// let e = one.exp();
+    ///
+    /// // ln(e) - 1 == 0
+    /// let abs_difference = (e.ln() - 1.0).abs();
+    ///
+    /// assert!(abs_difference <= f128::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn exp(self) -> f128 {
+        unsafe { intrinsics::expf128(self) }
+    }
+
+    /// Returns `2^(self)`.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let f = 2.0f128;
+    ///
+    /// // 2^2 - 4 == 0
+    /// let abs_difference = (f.exp2() - 4.0).abs();
+    ///
+    /// assert!(abs_difference <= f128::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn exp2(self) -> f128 {
+        unsafe { intrinsics::exp2f128(self) }
+    }
+
+    /// Returns the natural logarithm of the number.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let one = 1.0f128;
+    /// // e^1
+    /// let e = one.exp();
+    ///
+    /// // ln(e) - 1 == 0
+    /// let abs_difference = (e.ln() - 1.0).abs();
+    ///
+    /// assert!(abs_difference <= f128::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn ln(self) -> f128 {
+        unsafe { intrinsics::logf128(self) }
+    }
+
+    /// Returns the logarithm of the number with respect to an arbitrary base.
+    ///
+    /// The result might not be correctly rounded owing to implementation details;
+    /// `self.log2()` can produce more accurate results for base 2, and
+    /// `self.log10()` can produce more accurate results for base 10.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let five = 5.0f128;
+    ///
+    /// // log5(5) - 1 == 0
+    /// let abs_difference = (five.log(5.0) - 1.0).abs();
+    ///
+    /// assert!(abs_difference <= f128::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn log(self, base: f128) -> f128 {
+        self.ln() / base.ln()
+    }
+
+    /// Returns the base 2 logarithm of the number.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let two = 2.0f128;
+    ///
+    /// // log2(2) - 1 == 0
+    /// let abs_difference = (two.log2() - 1.0).abs();
+    ///
+    /// assert!(abs_difference <= f128::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn log2(self) -> f128 {
+        unsafe { intrinsics::log2f128(self) }
+    }
+
+    /// Returns the base 10 logarithm of the number.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let ten = 10.0f128;
+    ///
+    /// // log10(10) - 1 == 0
+    /// let abs_difference = (ten.log10() - 1.0).abs();
+    ///
+    /// assert!(abs_difference <= f128::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn log10(self) -> f128 {
+        unsafe { intrinsics::log10f128(self) }
+    }
+
+    /// Returns the cube root of a number.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    ///
+    /// This function currently corresponds to the `cbrtf128` from libc on Unix
+    /// and Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let x = 8.0f128;
+    ///
+    /// // x^(1/3) - 2 == 0
+    /// let abs_difference = (x.cbrt() - 2.0).abs();
+    ///
+    /// assert!(abs_difference <= f128::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn cbrt(self) -> f128 {
+        unsafe { cmath::cbrtf128(self) }
+    }
+
+    /// Compute the distance between the origin and a point (`x`, `y`) on the
+    /// Euclidean plane. Equivalently, compute the length of the hypotenuse of a
+    /// right-angle triangle with other sides having length `x.abs()` and
+    /// `y.abs()`.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    ///
+    /// This function currently corresponds to the `hypotf128` from libc on Unix
+    /// and Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let x = 2.0f128;
+    /// let y = 3.0f128;
+    ///
+    /// // sqrt(x^2 + y^2)
+    /// let abs_difference = (x.hypot(y) - (x.powi(2) + y.powi(2)).sqrt()).abs();
+    ///
+    /// assert!(abs_difference <= f128::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn hypot(self, other: f128) -> f128 {
+        unsafe { cmath::hypotf128(self, other) }
+    }
+
+    /// Computes the sine of a number (in radians).
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let x = std::f128::consts::FRAC_PI_2;
+    ///
+    /// let abs_difference = (x.sin() - 1.0).abs();
+    ///
+    /// assert!(abs_difference <= f128::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn sin(self) -> f128 {
+        unsafe { intrinsics::sinf128(self) }
+    }
+
+    /// Computes the cosine of a number (in radians).
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let x = 2.0 * std::f128::consts::PI;
+    ///
+    /// let abs_difference = (x.cos() - 1.0).abs();
+    ///
+    /// assert!(abs_difference <= f128::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn cos(self) -> f128 {
+        unsafe { intrinsics::cosf128(self) }
+    }
+
+    /// Computes the tangent of a number (in radians).
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `tanf128` from libc on Unix and
+    /// Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let x = std::f128::consts::FRAC_PI_4;
+    /// let abs_difference = (x.tan() - 1.0).abs();
+    ///
+    /// assert!(abs_difference <= f128::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn tan(self) -> f128 {
+        unsafe { cmath::tanf128(self) }
+    }
+
+    /// Computes the arcsine of a number. Return value is in radians in
+    /// the range [-pi/2, pi/2] or NaN if the number is outside the range
+    /// [-1, 1].
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `asinf128` from libc on Unix
+    /// and Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let f = std::f128::consts::FRAC_PI_2;
+    ///
+    /// // asin(sin(pi/2))
+    /// let abs_difference = (f.sin().asin() - std::f128::consts::FRAC_PI_2).abs();
+    ///
+    /// assert!(abs_difference <= f128::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[doc(alias = "arcsin")]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn asin(self) -> f128 {
+        unsafe { cmath::asinf128(self) }
+    }
+
+    /// Computes the arccosine of a number. Return value is in radians in
+    /// the range [0, pi] or NaN if the number is outside the range
+    /// [-1, 1].
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `acosf128` from libc on Unix
+    /// and Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let f = std::f128::consts::FRAC_PI_4;
+    ///
+    /// // acos(cos(pi/4))
+    /// let abs_difference = (f.cos().acos() - std::f128::consts::FRAC_PI_4).abs();
+    ///
+    /// assert!(abs_difference <= f128::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[doc(alias = "arccos")]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn acos(self) -> f128 {
+        unsafe { cmath::acosf128(self) }
+    }
+
+    /// Computes the arctangent of a number. Return value is in radians in the
+    /// range [-pi/2, pi/2];
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `atanf128` from libc on Unix
+    /// and Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let f = 1.0f128;
+    ///
+    /// // atan(tan(1))
+    /// let abs_difference = (f.tan().atan() - 1.0).abs();
+    ///
+    /// assert!(abs_difference <= f128::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[doc(alias = "arctan")]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn atan(self) -> f128 {
+        unsafe { cmath::atanf128(self) }
+    }
+
+    /// Computes the four quadrant arctangent of `self` (`y`) and `other` (`x`) in radians.
+    ///
+    /// * `x = 0`, `y = 0`: `0`
+    /// * `x >= 0`: `arctan(y/x)` -> `[-pi/2, pi/2]`
+    /// * `y >= 0`: `arctan(y/x) + pi` -> `(pi/2, pi]`
+    /// * `y < 0`: `arctan(y/x) - pi` -> `(-pi, -pi/2)`
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `atan2f128` from libc on Unix
+    /// and Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// // Positive angles measured counter-clockwise
+    /// // from positive x axis
+    /// // -pi/4 radians (45 deg clockwise)
+    /// let x1 = 3.0f128;
+    /// let y1 = -3.0f128;
+    ///
+    /// // 3pi/4 radians (135 deg counter-clockwise)
+    /// let x2 = -3.0f128;
+    /// let y2 = 3.0f128;
+    ///
+    /// let abs_difference_1 = (y1.atan2(x1) - (-std::f128::consts::FRAC_PI_4)).abs();
+    /// let abs_difference_2 = (y2.atan2(x2) - (3.0 * std::f128::consts::FRAC_PI_4)).abs();
+    ///
+    /// assert!(abs_difference_1 <= f128::EPSILON);
+    /// assert!(abs_difference_2 <= f128::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn atan2(self, other: f128) -> f128 {
+        unsafe { cmath::atan2f128(self, other) }
+    }
+
+    /// Simultaneously computes the sine and cosine of the number, `x`. Returns
+    /// `(sin(x), cos(x))`.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `(f128::sin(x),
+    /// f128::cos(x))`. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let x = std::f128::consts::FRAC_PI_4;
+    /// let f = x.sin_cos();
+    ///
+    /// let abs_difference_0 = (f.0 - x.sin()).abs();
+    /// let abs_difference_1 = (f.1 - x.cos()).abs();
+    ///
+    /// assert!(abs_difference_0 <= f128::EPSILON);
+    /// assert!(abs_difference_1 <= f128::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[doc(alias = "sincos")]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    pub fn sin_cos(self) -> (f128, f128) {
+        (self.sin(), self.cos())
+    }
+
+    /// Returns `e^(self) - 1` in a way that is accurate even if the
+    /// number is close to zero.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `expm1f128` from libc on Unix
+    /// and Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let x = 1e-8_f128;
+    ///
+    /// // for very small x, e^x is approximately 1 + x + x^2 / 2
+    /// let approx = x + x * x / 2.0;
+    /// let abs_difference = (x.exp_m1() - approx).abs();
+    ///
+    /// assert!(abs_difference < 1e-10);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn exp_m1(self) -> f128 {
+        unsafe { cmath::expm1f128(self) }
+    }
+
+    /// Returns `ln(1+n)` (natural logarithm) more accurately than if
+    /// the operations were performed separately.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `log1pf128` from libc on Unix
+    /// and Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let x = 1e-8_f128;
+    ///
+    /// // for very small x, ln(1 + x) is approximately x - x^2 / 2
+    /// let approx = x - x * x / 2.0;
+    /// let abs_difference = (x.ln_1p() - approx).abs();
+    ///
+    /// assert!(abs_difference < 1e-10);
+    /// # }
+    /// ```
+    #[inline]
+    #[doc(alias = "log1p")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    pub fn ln_1p(self) -> f128 {
+        unsafe { cmath::log1pf128(self) }
+    }
+
+    /// Hyperbolic sine function.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `sinhf128` from libc on Unix
+    /// and Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let e = std::f128::consts::E;
+    /// let x = 1.0f128;
+    ///
+    /// let f = x.sinh();
+    /// // Solving sinh() at 1 gives `(e^2-1)/(2e)`
+    /// let g = ((e * e) - 1.0) / (2.0 * e);
+    /// let abs_difference = (f - g).abs();
+    ///
+    /// assert!(abs_difference <= f128::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn sinh(self) -> f128 {
+        unsafe { cmath::sinhf128(self) }
+    }
+
+    /// Hyperbolic cosine function.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `coshf128` from libc on Unix
+    /// and Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let e = std::f128::consts::E;
+    /// let x = 1.0f128;
+    /// let f = x.cosh();
+    /// // Solving cosh() at 1 gives this result
+    /// let g = ((e * e) + 1.0) / (2.0 * e);
+    /// let abs_difference = (f - g).abs();
+    ///
+    /// // Same result
+    /// assert!(abs_difference <= f128::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn cosh(self) -> f128 {
+        unsafe { cmath::coshf128(self) }
+    }
+
+    /// Hyperbolic tangent function.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `tanhf128` from libc on Unix
+    /// and Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let e = std::f128::consts::E;
+    /// let x = 1.0f128;
+    ///
+    /// let f = x.tanh();
+    /// // Solving tanh() at 1 gives `(1 - e^(-2))/(1 + e^(-2))`
+    /// let g = (1.0 - e.powi(-2)) / (1.0 + e.powi(-2));
+    /// let abs_difference = (f - g).abs();
+    ///
+    /// assert!(abs_difference <= f128::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn tanh(self) -> f128 {
+        unsafe { cmath::tanhf128(self) }
+    }
+
+    /// Inverse hyperbolic sine function.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let x = 1.0f128;
+    /// let f = x.sinh().asinh();
+    ///
+    /// let abs_difference = (f - x).abs();
+    ///
+    /// assert!(abs_difference <= f128::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[doc(alias = "arcsinh")]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn asinh(self) -> f128 {
+        let ax = self.abs();
+        let ix = 1.0 / ax;
+        (ax + (ax / (Self::hypot(1.0, ix) + ix))).ln_1p().copysign(self)
+    }
+
+    /// Inverse hyperbolic cosine function.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let x = 1.0f128;
+    /// let f = x.cosh().acosh();
+    ///
+    /// let abs_difference = (f - x).abs();
+    ///
+    /// assert!(abs_difference <= f128::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[doc(alias = "arccosh")]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn acosh(self) -> f128 {
+        if self < 1.0 {
+            Self::NAN
+        } else {
+            (self + ((self - 1.0).sqrt() * (self + 1.0).sqrt())).ln()
+        }
+    }
+
+    /// Inverse hyperbolic tangent function.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let e = std::f128::consts::E;
+    /// let f = e.tanh().atanh();
+    ///
+    /// let abs_difference = (f - e).abs();
+    ///
+    /// assert!(abs_difference <= 1e-5);
+    /// # }
+    /// ```
+    #[inline]
+    #[doc(alias = "arctanh")]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn atanh(self) -> f128 {
+        0.5 * ((2.0 * self) / (1.0 - self)).ln_1p()
+    }
+
+    /// Gamma function.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `tgammaf128` from libc on Unix
+    /// and Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// #![feature(float_gamma)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let x = 5.0f128;
+    ///
+    /// let abs_difference = (x.gamma() - 24.0).abs();
+    ///
+    /// assert!(abs_difference <= f128::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn gamma(self) -> f128 {
+        unsafe { cmath::tgammaf128(self) }
+    }
+
+    /// Natural logarithm of the absolute value of the gamma function
+    ///
+    /// The integer part of the tuple indicates the sign of the gamma function.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `lgammaf128_r` from libc on Unix
+    /// and Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f128)]
+    /// #![feature(float_gamma)]
+    /// # #[cfg(reliable_f128_math)] {
+    ///
+    /// let x = 2.0f128;
+    ///
+    /// let abs_difference = (x.ln_gamma().0 - 0.0).abs();
+    ///
+    /// assert!(abs_difference <= f128::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f128", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn ln_gamma(self) -> (f128, i32) {
+        let mut signgamp: i32 = 0;
+        let x = unsafe { cmath::lgammaf128_r(self, &mut signgamp) };
+        (x, signgamp)
+    }
 }
diff --git a/library/std/src/f128/tests.rs b/library/std/src/f128/tests.rs
index 162c8dbad81..df806a639f6 100644
--- a/library/std/src/f128/tests.rs
+++ b/library/std/src/f128/tests.rs
@@ -4,6 +4,21 @@
 use crate::f128::consts;
 use crate::num::{FpCategory as Fp, *};
 
+// Note these tolerances make sense around zero, but not for more extreme exponents.
+
+/// For operations that are near exact, usually not involving math of different
+/// signs.
+const TOL_PRECISE: f128 = 1e-28;
+
+/// Default tolerances. Works for values that should be near precise but not exact. Roughly
+/// the precision carried by `100 * 100`.
+const TOL: f128 = 1e-12;
+
+/// Tolerances for math that is allowed to be imprecise, usually due to multiple chained
+/// operations.
+#[cfg(reliable_f128_math)]
+const TOL_IMPR: f128 = 1e-10;
+
 /// Smallest number
 const TINY_BITS: u128 = 0x1;
 
@@ -191,9 +206,100 @@ fn test_classify() {
     assert_eq!(1e-4932f128.classify(), Fp::Subnormal);
 }
 
-// FIXME(f16_f128): add missing math functions when available
+#[test]
+#[cfg(reliable_f128_math)]
+fn test_floor() {
+    assert_approx_eq!(1.0f128.floor(), 1.0f128, TOL_PRECISE);
+    assert_approx_eq!(1.3f128.floor(), 1.0f128, TOL_PRECISE);
+    assert_approx_eq!(1.5f128.floor(), 1.0f128, TOL_PRECISE);
+    assert_approx_eq!(1.7f128.floor(), 1.0f128, TOL_PRECISE);
+    assert_approx_eq!(0.0f128.floor(), 0.0f128, TOL_PRECISE);
+    assert_approx_eq!((-0.0f128).floor(), -0.0f128, TOL_PRECISE);
+    assert_approx_eq!((-1.0f128).floor(), -1.0f128, TOL_PRECISE);
+    assert_approx_eq!((-1.3f128).floor(), -2.0f128, TOL_PRECISE);
+    assert_approx_eq!((-1.5f128).floor(), -2.0f128, TOL_PRECISE);
+    assert_approx_eq!((-1.7f128).floor(), -2.0f128, TOL_PRECISE);
+}
+
+#[test]
+#[cfg(reliable_f128_math)]
+fn test_ceil() {
+    assert_approx_eq!(1.0f128.ceil(), 1.0f128, TOL_PRECISE);
+    assert_approx_eq!(1.3f128.ceil(), 2.0f128, TOL_PRECISE);
+    assert_approx_eq!(1.5f128.ceil(), 2.0f128, TOL_PRECISE);
+    assert_approx_eq!(1.7f128.ceil(), 2.0f128, TOL_PRECISE);
+    assert_approx_eq!(0.0f128.ceil(), 0.0f128, TOL_PRECISE);
+    assert_approx_eq!((-0.0f128).ceil(), -0.0f128, TOL_PRECISE);
+    assert_approx_eq!((-1.0f128).ceil(), -1.0f128, TOL_PRECISE);
+    assert_approx_eq!((-1.3f128).ceil(), -1.0f128, TOL_PRECISE);
+    assert_approx_eq!((-1.5f128).ceil(), -1.0f128, TOL_PRECISE);
+    assert_approx_eq!((-1.7f128).ceil(), -1.0f128, TOL_PRECISE);
+}
+
+#[test]
+#[cfg(reliable_f128_math)]
+fn test_round() {
+    assert_approx_eq!(2.5f128.round(), 3.0f128, TOL_PRECISE);
+    assert_approx_eq!(1.0f128.round(), 1.0f128, TOL_PRECISE);
+    assert_approx_eq!(1.3f128.round(), 1.0f128, TOL_PRECISE);
+    assert_approx_eq!(1.5f128.round(), 2.0f128, TOL_PRECISE);
+    assert_approx_eq!(1.7f128.round(), 2.0f128, TOL_PRECISE);
+    assert_approx_eq!(0.0f128.round(), 0.0f128, TOL_PRECISE);
+    assert_approx_eq!((-0.0f128).round(), -0.0f128, TOL_PRECISE);
+    assert_approx_eq!((-1.0f128).round(), -1.0f128, TOL_PRECISE);
+    assert_approx_eq!((-1.3f128).round(), -1.0f128, TOL_PRECISE);
+    assert_approx_eq!((-1.5f128).round(), -2.0f128, TOL_PRECISE);
+    assert_approx_eq!((-1.7f128).round(), -2.0f128, TOL_PRECISE);
+}
+
+#[test]
+#[cfg(reliable_f128_math)]
+fn test_round_ties_even() {
+    assert_approx_eq!(2.5f128.round_ties_even(), 2.0f128, TOL_PRECISE);
+    assert_approx_eq!(1.0f128.round_ties_even(), 1.0f128, TOL_PRECISE);
+    assert_approx_eq!(1.3f128.round_ties_even(), 1.0f128, TOL_PRECISE);
+    assert_approx_eq!(1.5f128.round_ties_even(), 2.0f128, TOL_PRECISE);
+    assert_approx_eq!(1.7f128.round_ties_even(), 2.0f128, TOL_PRECISE);
+    assert_approx_eq!(0.0f128.round_ties_even(), 0.0f128, TOL_PRECISE);
+    assert_approx_eq!((-0.0f128).round_ties_even(), -0.0f128, TOL_PRECISE);
+    assert_approx_eq!((-1.0f128).round_ties_even(), -1.0f128, TOL_PRECISE);
+    assert_approx_eq!((-1.3f128).round_ties_even(), -1.0f128, TOL_PRECISE);
+    assert_approx_eq!((-1.5f128).round_ties_even(), -2.0f128, TOL_PRECISE);
+    assert_approx_eq!((-1.7f128).round_ties_even(), -2.0f128, TOL_PRECISE);
+}
+
+#[test]
+#[cfg(reliable_f128_math)]
+fn test_trunc() {
+    assert_approx_eq!(1.0f128.trunc(), 1.0f128, TOL_PRECISE);
+    assert_approx_eq!(1.3f128.trunc(), 1.0f128, TOL_PRECISE);
+    assert_approx_eq!(1.5f128.trunc(), 1.0f128, TOL_PRECISE);
+    assert_approx_eq!(1.7f128.trunc(), 1.0f128, TOL_PRECISE);
+    assert_approx_eq!(0.0f128.trunc(), 0.0f128, TOL_PRECISE);
+    assert_approx_eq!((-0.0f128).trunc(), -0.0f128, TOL_PRECISE);
+    assert_approx_eq!((-1.0f128).trunc(), -1.0f128, TOL_PRECISE);
+    assert_approx_eq!((-1.3f128).trunc(), -1.0f128, TOL_PRECISE);
+    assert_approx_eq!((-1.5f128).trunc(), -1.0f128, TOL_PRECISE);
+    assert_approx_eq!((-1.7f128).trunc(), -1.0f128, TOL_PRECISE);
+}
+
+#[test]
+#[cfg(reliable_f128_math)]
+fn test_fract() {
+    assert_approx_eq!(1.0f128.fract(), 0.0f128, TOL_PRECISE);
+    assert_approx_eq!(1.3f128.fract(), 0.3f128, TOL_PRECISE);
+    assert_approx_eq!(1.5f128.fract(), 0.5f128, TOL_PRECISE);
+    assert_approx_eq!(1.7f128.fract(), 0.7f128, TOL_PRECISE);
+    assert_approx_eq!(0.0f128.fract(), 0.0f128, TOL_PRECISE);
+    assert_approx_eq!((-0.0f128).fract(), -0.0f128, TOL_PRECISE);
+    assert_approx_eq!((-1.0f128).fract(), -0.0f128, TOL_PRECISE);
+    assert_approx_eq!((-1.3f128).fract(), -0.3f128, TOL_PRECISE);
+    assert_approx_eq!((-1.5f128).fract(), -0.5f128, TOL_PRECISE);
+    assert_approx_eq!((-1.7f128).fract(), -0.7f128, TOL_PRECISE);
+}
 
 #[test]
+#[cfg(reliable_f128_math)]
 fn test_abs() {
     assert_eq!(f128::INFINITY.abs(), f128::INFINITY);
     assert_eq!(1f128.abs(), 1f128);
@@ -293,6 +399,24 @@ fn test_next_down() {
 }
 
 #[test]
+#[cfg(reliable_f128_math)]
+fn test_mul_add() {
+    let nan: f128 = f128::NAN;
+    let inf: f128 = f128::INFINITY;
+    let neg_inf: f128 = f128::NEG_INFINITY;
+    assert_approx_eq!(12.3f128.mul_add(4.5, 6.7), 62.05, TOL_PRECISE);
+    assert_approx_eq!((-12.3f128).mul_add(-4.5, -6.7), 48.65, TOL_PRECISE);
+    assert_approx_eq!(0.0f128.mul_add(8.9, 1.2), 1.2, TOL_PRECISE);
+    assert_approx_eq!(3.4f128.mul_add(-0.0, 5.6), 5.6, TOL_PRECISE);
+    assert!(nan.mul_add(7.8, 9.0).is_nan());
+    assert_eq!(inf.mul_add(7.8, 9.0), inf);
+    assert_eq!(neg_inf.mul_add(7.8, 9.0), neg_inf);
+    assert_eq!(8.9f128.mul_add(inf, 3.2), inf);
+    assert_eq!((-3.2f128).mul_add(2.4, neg_inf), neg_inf);
+}
+
+#[test]
+#[cfg(reliable_f16_math)]
 fn test_recip() {
     let nan: f128 = f128::NAN;
     let inf: f128 = f128::INFINITY;
@@ -301,11 +425,161 @@ fn test_recip() {
     assert_eq!(2.0f128.recip(), 0.5);
     assert_eq!((-0.4f128).recip(), -2.5);
     assert_eq!(0.0f128.recip(), inf);
+    assert_approx_eq!(
+        f128::MAX.recip(),
+        8.40525785778023376565669454330438228902076605e-4933,
+        1e-4900
+    );
     assert!(nan.recip().is_nan());
     assert_eq!(inf.recip(), 0.0);
     assert_eq!(neg_inf.recip(), 0.0);
 }
 
+// Many math functions allow for less accurate results, so the next tolerance up is used
+
+#[test]
+#[cfg(reliable_f128_math)]
+fn test_powi() {
+    let nan: f128 = f128::NAN;
+    let inf: f128 = f128::INFINITY;
+    let neg_inf: f128 = f128::NEG_INFINITY;
+    assert_eq!(1.0f128.powi(1), 1.0);
+    assert_approx_eq!((-3.1f128).powi(2), 9.6100000000000005506706202140776519387, TOL);
+    assert_approx_eq!(5.9f128.powi(-2), 0.028727377190462507313100483690639638451, TOL);
+    assert_eq!(8.3f128.powi(0), 1.0);
+    assert!(nan.powi(2).is_nan());
+    assert_eq!(inf.powi(3), inf);
+    assert_eq!(neg_inf.powi(2), inf);
+}
+
+#[test]
+#[cfg(reliable_f128_math)]
+fn test_powf() {
+    let nan: f128 = f128::NAN;
+    let inf: f128 = f128::INFINITY;
+    let neg_inf: f128 = f128::NEG_INFINITY;
+    assert_eq!(1.0f128.powf(1.0), 1.0);
+    assert_approx_eq!(3.4f128.powf(4.5), 246.40818323761892815995637964326426756, TOL_IMPR);
+    assert_approx_eq!(2.7f128.powf(-3.2), 0.041652009108526178281070304373500889273, TOL_IMPR);
+    assert_approx_eq!((-3.1f128).powf(2.0), 9.6100000000000005506706202140776519387, TOL_IMPR);
+    assert_approx_eq!(5.9f128.powf(-2.0), 0.028727377190462507313100483690639638451, TOL_IMPR);
+    assert_eq!(8.3f128.powf(0.0), 1.0);
+    assert!(nan.powf(2.0).is_nan());
+    assert_eq!(inf.powf(2.0), inf);
+    assert_eq!(neg_inf.powf(3.0), neg_inf);
+}
+
+#[test]
+#[cfg(reliable_f128_math)]
+fn test_sqrt_domain() {
+    assert!(f128::NAN.sqrt().is_nan());
+    assert!(f128::NEG_INFINITY.sqrt().is_nan());
+    assert!((-1.0f128).sqrt().is_nan());
+    assert_eq!((-0.0f128).sqrt(), -0.0);
+    assert_eq!(0.0f128.sqrt(), 0.0);
+    assert_eq!(1.0f128.sqrt(), 1.0);
+    assert_eq!(f128::INFINITY.sqrt(), f128::INFINITY);
+}
+
+#[test]
+#[cfg(reliable_f128_math)]
+fn test_exp() {
+    assert_eq!(1.0, 0.0f128.exp());
+    assert_approx_eq!(consts::E, 1.0f128.exp(), TOL);
+    assert_approx_eq!(148.41315910257660342111558004055227962348775, 5.0f128.exp(), TOL);
+
+    let inf: f128 = f128::INFINITY;
+    let neg_inf: f128 = f128::NEG_INFINITY;
+    let nan: f128 = f128::NAN;
+    assert_eq!(inf, inf.exp());
+    assert_eq!(0.0, neg_inf.exp());
+    assert!(nan.exp().is_nan());
+}
+
+#[test]
+#[cfg(reliable_f128_math)]
+fn test_exp2() {
+    assert_eq!(32.0, 5.0f128.exp2());
+    assert_eq!(1.0, 0.0f128.exp2());
+
+    let inf: f128 = f128::INFINITY;
+    let neg_inf: f128 = f128::NEG_INFINITY;
+    let nan: f128 = f128::NAN;
+    assert_eq!(inf, inf.exp2());
+    assert_eq!(0.0, neg_inf.exp2());
+    assert!(nan.exp2().is_nan());
+}
+
+#[test]
+#[cfg(reliable_f128_math)]
+fn test_ln() {
+    let nan: f128 = f128::NAN;
+    let inf: f128 = f128::INFINITY;
+    let neg_inf: f128 = f128::NEG_INFINITY;
+    assert_approx_eq!(1.0f128.exp().ln(), 1.0, TOL);
+    assert!(nan.ln().is_nan());
+    assert_eq!(inf.ln(), inf);
+    assert!(neg_inf.ln().is_nan());
+    assert!((-2.3f128).ln().is_nan());
+    assert_eq!((-0.0f128).ln(), neg_inf);
+    assert_eq!(0.0f128.ln(), neg_inf);
+    assert_approx_eq!(4.0f128.ln(), 1.3862943611198906188344642429163531366, TOL);
+}
+
+#[test]
+#[cfg(reliable_f128_math)]
+fn test_log() {
+    let nan: f128 = f128::NAN;
+    let inf: f128 = f128::INFINITY;
+    let neg_inf: f128 = f128::NEG_INFINITY;
+    assert_eq!(10.0f128.log(10.0), 1.0);
+    assert_approx_eq!(2.3f128.log(3.5), 0.66485771361478710036766645911922010272, TOL);
+    assert_eq!(1.0f128.exp().log(1.0f128.exp()), 1.0);
+    assert!(1.0f128.log(1.0).is_nan());
+    assert!(1.0f128.log(-13.9).is_nan());
+    assert!(nan.log(2.3).is_nan());
+    assert_eq!(inf.log(10.0), inf);
+    assert!(neg_inf.log(8.8).is_nan());
+    assert!((-2.3f128).log(0.1).is_nan());
+    assert_eq!((-0.0f128).log(2.0), neg_inf);
+    assert_eq!(0.0f128.log(7.0), neg_inf);
+}
+
+#[test]
+#[cfg(reliable_f128_math)]
+fn test_log2() {
+    let nan: f128 = f128::NAN;
+    let inf: f128 = f128::INFINITY;
+    let neg_inf: f128 = f128::NEG_INFINITY;
+    assert_approx_eq!(10.0f128.log2(), 3.32192809488736234787031942948939017, TOL);
+    assert_approx_eq!(2.3f128.log2(), 1.2016338611696504130002982471978765921, TOL);
+    assert_approx_eq!(1.0f128.exp().log2(), 1.4426950408889634073599246810018921381, TOL);
+    assert!(nan.log2().is_nan());
+    assert_eq!(inf.log2(), inf);
+    assert!(neg_inf.log2().is_nan());
+    assert!((-2.3f128).log2().is_nan());
+    assert_eq!((-0.0f128).log2(), neg_inf);
+    assert_eq!(0.0f128.log2(), neg_inf);
+}
+
+#[test]
+#[cfg(reliable_f128_math)]
+fn test_log10() {
+    let nan: f128 = f128::NAN;
+    let inf: f128 = f128::INFINITY;
+    let neg_inf: f128 = f128::NEG_INFINITY;
+    assert_eq!(10.0f128.log10(), 1.0);
+    assert_approx_eq!(2.3f128.log10(), 0.36172783601759284532595218865859309898, TOL);
+    assert_approx_eq!(1.0f128.exp().log10(), 0.43429448190325182765112891891660508222, TOL);
+    assert_eq!(1.0f128.log10(), 0.0);
+    assert!(nan.log10().is_nan());
+    assert_eq!(inf.log10(), inf);
+    assert!(neg_inf.log10().is_nan());
+    assert!((-2.3f128).log10().is_nan());
+    assert_eq!((-0.0f128).log10(), neg_inf);
+    assert_eq!(0.0f128.log10(), neg_inf);
+}
+
 #[test]
 fn test_to_degrees() {
     let pi: f128 = consts::PI;
@@ -313,8 +587,8 @@ fn test_to_degrees() {
     let inf: f128 = f128::INFINITY;
     let neg_inf: f128 = f128::NEG_INFINITY;
     assert_eq!(0.0f128.to_degrees(), 0.0);
-    assert_approx_eq!((-5.8f128).to_degrees(), -332.315521);
-    assert_eq!(pi.to_degrees(), 180.0);
+    assert_approx_eq!((-5.8f128).to_degrees(), -332.31552117587745090765431723855668471, TOL);
+    assert_approx_eq!(pi.to_degrees(), 180.0, TOL);
     assert!(nan.to_degrees().is_nan());
     assert_eq!(inf.to_degrees(), inf);
     assert_eq!(neg_inf.to_degrees(), neg_inf);
@@ -328,19 +602,122 @@ fn test_to_radians() {
     let inf: f128 = f128::INFINITY;
     let neg_inf: f128 = f128::NEG_INFINITY;
     assert_eq!(0.0f128.to_radians(), 0.0);
-    assert_approx_eq!(154.6f128.to_radians(), 2.698279);
-    assert_approx_eq!((-332.31f128).to_radians(), -5.799903);
+    assert_approx_eq!(154.6f128.to_radians(), 2.6982790235832334267135442069489767804, TOL);
+    assert_approx_eq!((-332.31f128).to_radians(), -5.7999036373023566567593094812182763013, TOL);
     // check approx rather than exact because round trip for pi doesn't fall on an exactly
     // representable value (unlike `f32` and `f64`).
-    assert_approx_eq!(180.0f128.to_radians(), pi);
+    assert_approx_eq!(180.0f128.to_radians(), pi, TOL_PRECISE);
     assert!(nan.to_radians().is_nan());
     assert_eq!(inf.to_radians(), inf);
     assert_eq!(neg_inf.to_radians(), neg_inf);
 }
 
 #[test]
+#[cfg(reliable_f128_math)]
+fn test_asinh() {
+    // Lower accuracy results are allowed, use increased tolerances
+    assert_eq!(0.0f128.asinh(), 0.0f128);
+    assert_eq!((-0.0f128).asinh(), -0.0f128);
+
+    let inf: f128 = f128::INFINITY;
+    let neg_inf: f128 = f128::NEG_INFINITY;
+    let nan: f128 = f128::NAN;
+    assert_eq!(inf.asinh(), inf);
+    assert_eq!(neg_inf.asinh(), neg_inf);
+    assert!(nan.asinh().is_nan());
+    assert!((-0.0f128).asinh().is_sign_negative());
+
+    // issue 63271
+    assert_approx_eq!(2.0f128.asinh(), 1.443635475178810342493276740273105f128, TOL_IMPR);
+    assert_approx_eq!((-2.0f128).asinh(), -1.443635475178810342493276740273105f128, TOL_IMPR);
+    // regression test for the catastrophic cancellation fixed in 72486
+    assert_approx_eq!(
+        (-67452098.07139316f128).asinh(),
+        -18.720075426274544393985484294000831757220,
+        TOL_IMPR
+    );
+
+    // test for low accuracy from issue 104548
+    assert_approx_eq!(60.0f128, 60.0f128.sinh().asinh(), TOL_IMPR);
+    // mul needed for approximate comparison to be meaningful
+    assert_approx_eq!(1.0f128, 1e-15f128.sinh().asinh() * 1e15f128, TOL_IMPR);
+}
+
+#[test]
+#[cfg(reliable_f128_math)]
+fn test_acosh() {
+    assert_eq!(1.0f128.acosh(), 0.0f128);
+    assert!(0.999f128.acosh().is_nan());
+
+    let inf: f128 = f128::INFINITY;
+    let neg_inf: f128 = f128::NEG_INFINITY;
+    let nan: f128 = f128::NAN;
+    assert_eq!(inf.acosh(), inf);
+    assert!(neg_inf.acosh().is_nan());
+    assert!(nan.acosh().is_nan());
+    assert_approx_eq!(2.0f128.acosh(), 1.31695789692481670862504634730796844f128, TOL_IMPR);
+    assert_approx_eq!(3.0f128.acosh(), 1.76274717403908605046521864995958461f128, TOL_IMPR);
+
+    // test for low accuracy from issue 104548
+    assert_approx_eq!(60.0f128, 60.0f128.cosh().acosh(), TOL_IMPR);
+}
+
+#[test]
+#[cfg(reliable_f128_math)]
+fn test_atanh() {
+    assert_eq!(0.0f128.atanh(), 0.0f128);
+    assert_eq!((-0.0f128).atanh(), -0.0f128);
+
+    let inf: f128 = f128::INFINITY;
+    let neg_inf: f128 = f128::NEG_INFINITY;
+    let nan: f128 = f128::NAN;
+    assert_eq!(1.0f128.atanh(), inf);
+    assert_eq!((-1.0f128).atanh(), neg_inf);
+    assert!(2f128.atanh().atanh().is_nan());
+    assert!((-2f128).atanh().atanh().is_nan());
+    assert!(inf.atanh().is_nan());
+    assert!(neg_inf.atanh().is_nan());
+    assert!(nan.atanh().is_nan());
+    assert_approx_eq!(0.5f128.atanh(), 0.54930614433405484569762261846126285f128, TOL_IMPR);
+    assert_approx_eq!((-0.5f128).atanh(), -0.54930614433405484569762261846126285f128, TOL_IMPR);
+}
+
+#[test]
+#[cfg(reliable_f128_math)]
+fn test_gamma() {
+    // precision can differ among platforms
+    assert_approx_eq!(1.0f128.gamma(), 1.0f128, TOL_IMPR);
+    assert_approx_eq!(2.0f128.gamma(), 1.0f128, TOL_IMPR);
+    assert_approx_eq!(3.0f128.gamma(), 2.0f128, TOL_IMPR);
+    assert_approx_eq!(4.0f128.gamma(), 6.0f128, TOL_IMPR);
+    assert_approx_eq!(5.0f128.gamma(), 24.0f128, TOL_IMPR);
+    assert_approx_eq!(0.5f128.gamma(), consts::PI.sqrt(), TOL_IMPR);
+    assert_approx_eq!((-0.5f128).gamma(), -2.0 * consts::PI.sqrt(), TOL_IMPR);
+    assert_eq!(0.0f128.gamma(), f128::INFINITY);
+    assert_eq!((-0.0f128).gamma(), f128::NEG_INFINITY);
+    assert!((-1.0f128).gamma().is_nan());
+    assert!((-2.0f128).gamma().is_nan());
+    assert!(f128::NAN.gamma().is_nan());
+    assert!(f128::NEG_INFINITY.gamma().is_nan());
+    assert_eq!(f128::INFINITY.gamma(), f128::INFINITY);
+    assert_eq!(1760.9f128.gamma(), f128::INFINITY);
+}
+
+#[test]
+#[cfg(reliable_f128_math)]
+fn test_ln_gamma() {
+    assert_approx_eq!(1.0f128.ln_gamma().0, 0.0f128, TOL_IMPR);
+    assert_eq!(1.0f128.ln_gamma().1, 1);
+    assert_approx_eq!(2.0f128.ln_gamma().0, 0.0f128, TOL_IMPR);
+    assert_eq!(2.0f128.ln_gamma().1, 1);
+    assert_approx_eq!(3.0f128.ln_gamma().0, 2.0f128.ln(), TOL_IMPR);
+    assert_eq!(3.0f128.ln_gamma().1, 1);
+    assert_approx_eq!((-0.5f128).ln_gamma().0, (2.0 * consts::PI.sqrt()).ln(), TOL_IMPR);
+    assert_eq!((-0.5f128).ln_gamma().1, -1);
+}
+
+#[test]
 fn test_real_consts() {
-    // FIXME(f16_f128): add math tests when available
     use super::consts;
 
     let pi: f128 = consts::PI;
@@ -351,29 +728,34 @@ fn test_real_consts() {
     let frac_pi_8: f128 = consts::FRAC_PI_8;
     let frac_1_pi: f128 = consts::FRAC_1_PI;
     let frac_2_pi: f128 = consts::FRAC_2_PI;
-    // let frac_2_sqrtpi: f128 = consts::FRAC_2_SQRT_PI;
-    // let sqrt2: f128 = consts::SQRT_2;
-    // let frac_1_sqrt2: f128 = consts::FRAC_1_SQRT_2;
-    // let e: f128 = consts::E;
-    // let log2_e: f128 = consts::LOG2_E;
-    // let log10_e: f128 = consts::LOG10_E;
-    // let ln_2: f128 = consts::LN_2;
-    // let ln_10: f128 = consts::LN_10;
-
-    assert_approx_eq!(frac_pi_2, pi / 2f128);
-    assert_approx_eq!(frac_pi_3, pi / 3f128);
-    assert_approx_eq!(frac_pi_4, pi / 4f128);
-    assert_approx_eq!(frac_pi_6, pi / 6f128);
-    assert_approx_eq!(frac_pi_8, pi / 8f128);
-    assert_approx_eq!(frac_1_pi, 1f128 / pi);
-    assert_approx_eq!(frac_2_pi, 2f128 / pi);
-    // assert_approx_eq!(frac_2_sqrtpi, 2f128 / pi.sqrt());
-    // assert_approx_eq!(sqrt2, 2f128.sqrt());
-    // assert_approx_eq!(frac_1_sqrt2, 1f128 / 2f128.sqrt());
-    // assert_approx_eq!(log2_e, e.log2());
-    // assert_approx_eq!(log10_e, e.log10());
-    // assert_approx_eq!(ln_2, 2f128.ln());
-    // assert_approx_eq!(ln_10, 10f128.ln());
+
+    assert_approx_eq!(frac_pi_2, pi / 2f128, TOL_PRECISE);
+    assert_approx_eq!(frac_pi_3, pi / 3f128, TOL_PRECISE);
+    assert_approx_eq!(frac_pi_4, pi / 4f128, TOL_PRECISE);
+    assert_approx_eq!(frac_pi_6, pi / 6f128, TOL_PRECISE);
+    assert_approx_eq!(frac_pi_8, pi / 8f128, TOL_PRECISE);
+    assert_approx_eq!(frac_1_pi, 1f128 / pi, TOL_PRECISE);
+    assert_approx_eq!(frac_2_pi, 2f128 / pi, TOL_PRECISE);
+
+    #[cfg(reliable_f128_math)]
+    {
+        let frac_2_sqrtpi: f128 = consts::FRAC_2_SQRT_PI;
+        let sqrt2: f128 = consts::SQRT_2;
+        let frac_1_sqrt2: f128 = consts::FRAC_1_SQRT_2;
+        let e: f128 = consts::E;
+        let log2_e: f128 = consts::LOG2_E;
+        let log10_e: f128 = consts::LOG10_E;
+        let ln_2: f128 = consts::LN_2;
+        let ln_10: f128 = consts::LN_10;
+
+        assert_approx_eq!(frac_2_sqrtpi, 2f128 / pi.sqrt(), TOL_PRECISE);
+        assert_approx_eq!(sqrt2, 2f128.sqrt(), TOL_PRECISE);
+        assert_approx_eq!(frac_1_sqrt2, 1f128 / 2f128.sqrt(), TOL_PRECISE);
+        assert_approx_eq!(log2_e, e.log2(), TOL_PRECISE);
+        assert_approx_eq!(log10_e, e.log10(), TOL_PRECISE);
+        assert_approx_eq!(ln_2, 2f128.ln(), TOL_PRECISE);
+        assert_approx_eq!(ln_10, 10f128.ln(), TOL_PRECISE);
+    }
 }
 
 #[test]
@@ -382,10 +764,10 @@ fn test_float_bits_conv() {
     assert_eq!((12.5f128).to_bits(), 0x40029000000000000000000000000000);
     assert_eq!((1337f128).to_bits(), 0x40094e40000000000000000000000000);
     assert_eq!((-14.25f128).to_bits(), 0xc002c800000000000000000000000000);
-    assert_approx_eq!(f128::from_bits(0x3fff0000000000000000000000000000), 1.0);
-    assert_approx_eq!(f128::from_bits(0x40029000000000000000000000000000), 12.5);
-    assert_approx_eq!(f128::from_bits(0x40094e40000000000000000000000000), 1337.0);
-    assert_approx_eq!(f128::from_bits(0xc002c800000000000000000000000000), -14.25);
+    assert_approx_eq!(f128::from_bits(0x3fff0000000000000000000000000000), 1.0, TOL_PRECISE);
+    assert_approx_eq!(f128::from_bits(0x40029000000000000000000000000000), 12.5, TOL_PRECISE);
+    assert_approx_eq!(f128::from_bits(0x40094e40000000000000000000000000), 1337.0, TOL_PRECISE);
+    assert_approx_eq!(f128::from_bits(0xc002c800000000000000000000000000), -14.25, TOL_PRECISE);
 
     // Check that NaNs roundtrip their bits regardless of signaling-ness
     // 0xA is 0b1010; 0x5 is 0b0101 -- so these two together clobbers all the mantissa bits
diff --git a/library/std/src/f16.rs b/library/std/src/f16.rs
index e3024defed7..10908332762 100644
--- a/library/std/src/f16.rs
+++ b/library/std/src/f16.rs
@@ -12,25 +12,180 @@ pub use core::f16::consts;
 
 #[cfg(not(test))]
 use crate::intrinsics;
+#[cfg(not(test))]
+use crate::sys::cmath;
 
 #[cfg(not(test))]
 impl f16 {
-    /// Raises a number to an integer power.
+    /// Returns the largest integer less than or equal to `self`.
     ///
-    /// Using this function is generally faster than using `powf`.
-    /// It might have a different sequence of rounding operations than `powf`,
-    /// so the results are not guaranteed to agree.
+    /// This function always returns the precise result.
     ///
-    /// # Unspecified precision
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let f = 3.7_f16;
+    /// let g = 3.0_f16;
+    /// let h = -3.7_f16;
     ///
-    /// The precision of this function is non-deterministic. This means it varies by platform, Rust version, and
-    /// can even differ within the same execution from one invocation to the next.
+    /// assert_eq!(f.floor(), 3.0);
+    /// assert_eq!(g.floor(), 3.0);
+    /// assert_eq!(h.floor(), -4.0);
+    /// # }
+    /// ```
     #[inline]
     #[rustc_allow_incoherent_impl]
     #[unstable(feature = "f16", issue = "116909")]
     #[must_use = "method returns a new number and does not mutate the original value"]
-    pub fn powi(self, n: i32) -> f16 {
-        unsafe { intrinsics::powif16(self, n) }
+    pub fn floor(self) -> f16 {
+        unsafe { intrinsics::floorf16(self) }
+    }
+
+    /// Returns the smallest integer greater than or equal to `self`.
+    ///
+    /// This function always returns the precise result.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let f = 3.01_f16;
+    /// let g = 4.0_f16;
+    ///
+    /// assert_eq!(f.ceil(), 4.0);
+    /// assert_eq!(g.ceil(), 4.0);
+    /// # }
+    /// ```
+    #[inline]
+    #[doc(alias = "ceiling")]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn ceil(self) -> f16 {
+        unsafe { intrinsics::ceilf16(self) }
+    }
+
+    /// Returns the nearest integer to `self`. If a value is half-way between two
+    /// integers, round away from `0.0`.
+    ///
+    /// This function always returns the precise result.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let f = 3.3_f16;
+    /// let g = -3.3_f16;
+    /// let h = -3.7_f16;
+    /// let i = 3.5_f16;
+    /// let j = 4.5_f16;
+    ///
+    /// assert_eq!(f.round(), 3.0);
+    /// assert_eq!(g.round(), -3.0);
+    /// assert_eq!(h.round(), -4.0);
+    /// assert_eq!(i.round(), 4.0);
+    /// assert_eq!(j.round(), 5.0);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn round(self) -> f16 {
+        unsafe { intrinsics::roundf16(self) }
+    }
+
+    /// Returns the nearest integer to a number. Rounds half-way cases to the number
+    /// with an even least significant digit.
+    ///
+    /// This function always returns the precise result.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let f = 3.3_f16;
+    /// let g = -3.3_f16;
+    /// let h = 3.5_f16;
+    /// let i = 4.5_f16;
+    ///
+    /// assert_eq!(f.round_ties_even(), 3.0);
+    /// assert_eq!(g.round_ties_even(), -3.0);
+    /// assert_eq!(h.round_ties_even(), 4.0);
+    /// assert_eq!(i.round_ties_even(), 4.0);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn round_ties_even(self) -> f16 {
+        unsafe { intrinsics::rintf16(self) }
+    }
+
+    /// Returns the integer part of `self`.
+    /// This means that non-integer numbers are always truncated towards zero.
+    ///
+    /// This function always returns the precise result.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let f = 3.7_f16;
+    /// let g = 3.0_f16;
+    /// let h = -3.7_f16;
+    ///
+    /// assert_eq!(f.trunc(), 3.0);
+    /// assert_eq!(g.trunc(), 3.0);
+    /// assert_eq!(h.trunc(), -3.0);
+    /// # }
+    /// ```
+    #[inline]
+    #[doc(alias = "truncate")]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn trunc(self) -> f16 {
+        unsafe { intrinsics::truncf16(self) }
+    }
+
+    /// Returns the fractional part of `self`.
+    ///
+    /// This function always returns the precise result.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let x = 3.6_f16;
+    /// let y = -3.6_f16;
+    /// let abs_difference_x = (x.fract() - 0.6).abs();
+    /// let abs_difference_y = (y.fract() - (-0.6)).abs();
+    ///
+    /// assert!(abs_difference_x <= f16::EPSILON);
+    /// assert!(abs_difference_y <= f16::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn fract(self) -> f16 {
+        self - self.trunc()
     }
 
     /// Computes the absolute value of `self`.
@@ -60,4 +215,1127 @@ impl f16 {
         // FIXME(f16_f128): replace with `intrinsics::fabsf16` when available
         Self::from_bits(self.to_bits() & !(1 << 15))
     }
+
+    /// Returns a number that represents the sign of `self`.
+    ///
+    /// - `1.0` if the number is positive, `+0.0` or `INFINITY`
+    /// - `-1.0` if the number is negative, `-0.0` or `NEG_INFINITY`
+    /// - NaN if the number is NaN
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let f = 3.5_f16;
+    ///
+    /// assert_eq!(f.signum(), 1.0);
+    /// assert_eq!(f16::NEG_INFINITY.signum(), -1.0);
+    ///
+    /// assert!(f16::NAN.signum().is_nan());
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn signum(self) -> f16 {
+        if self.is_nan() { Self::NAN } else { 1.0_f16.copysign(self) }
+    }
+
+    /// Returns a number composed of the magnitude of `self` and the sign of
+    /// `sign`.
+    ///
+    /// Equal to `self` if the sign of `self` and `sign` are the same, otherwise
+    /// equal to `-self`. If `self` is a NaN, then a NaN with the sign bit of
+    /// `sign` is returned. Note, however, that conserving the sign bit on NaN
+    /// across arithmetical operations is not generally guaranteed.
+    /// See [explanation of NaN as a special value](primitive@f16) for more info.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let f = 3.5_f16;
+    ///
+    /// assert_eq!(f.copysign(0.42), 3.5_f16);
+    /// assert_eq!(f.copysign(-0.42), -3.5_f16);
+    /// assert_eq!((-f).copysign(0.42), 3.5_f16);
+    /// assert_eq!((-f).copysign(-0.42), -3.5_f16);
+    ///
+    /// assert!(f16::NAN.copysign(1.0).is_nan());
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn copysign(self, sign: f16) -> f16 {
+        unsafe { intrinsics::copysignf16(self, sign) }
+    }
+
+    /// Fused multiply-add. Computes `(self * a) + b` with only one rounding
+    /// error, yielding a more accurate result than an unfused multiply-add.
+    ///
+    /// Using `mul_add` *may* be more performant than an unfused multiply-add if
+    /// the target architecture has a dedicated `fma` CPU instruction. However,
+    /// this is not always true, and will be heavily dependant on designing
+    /// algorithms with specific target hardware in mind.
+    ///
+    /// # Precision
+    ///
+    /// The result of this operation is guaranteed to be the rounded
+    /// infinite-precision result. It is specified by IEEE 754 as
+    /// `fusedMultiplyAdd` and guaranteed not to change.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let m = 10.0_f16;
+    /// let x = 4.0_f16;
+    /// let b = 60.0_f16;
+    ///
+    /// assert_eq!(m.mul_add(x, b), 100.0);
+    /// assert_eq!(m * x + b, 100.0);
+    ///
+    /// let one_plus_eps = 1.0_f16 + f16::EPSILON;
+    /// let one_minus_eps = 1.0_f16 - f16::EPSILON;
+    /// let minus_one = -1.0_f16;
+    ///
+    /// // The exact result (1 + eps) * (1 - eps) = 1 - eps * eps.
+    /// assert_eq!(one_plus_eps.mul_add(one_minus_eps, minus_one), -f16::EPSILON * f16::EPSILON);
+    /// // Different rounding with the non-fused multiply and add.
+    /// assert_eq!(one_plus_eps * one_minus_eps + minus_one, 0.0);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn mul_add(self, a: f16, b: f16) -> f16 {
+        unsafe { intrinsics::fmaf16(self, a, b) }
+    }
+
+    /// Calculates Euclidean division, the matching method for `rem_euclid`.
+    ///
+    /// This computes the integer `n` such that
+    /// `self = n * rhs + self.rem_euclid(rhs)`.
+    /// In other words, the result is `self / rhs` rounded to the integer `n`
+    /// such that `self >= n * rhs`.
+    ///
+    /// # Precision
+    ///
+    /// The result of this operation is guaranteed to be the rounded
+    /// infinite-precision result.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let a: f16 = 7.0;
+    /// let b = 4.0;
+    /// assert_eq!(a.div_euclid(b), 1.0); // 7.0 > 4.0 * 1.0
+    /// assert_eq!((-a).div_euclid(b), -2.0); // -7.0 >= 4.0 * -2.0
+    /// assert_eq!(a.div_euclid(-b), -1.0); // 7.0 >= -4.0 * -1.0
+    /// assert_eq!((-a).div_euclid(-b), 2.0); // -7.0 >= -4.0 * 2.0
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn div_euclid(self, rhs: f16) -> f16 {
+        let q = (self / rhs).trunc();
+        if self % rhs < 0.0 {
+            return if rhs > 0.0 { q - 1.0 } else { q + 1.0 };
+        }
+        q
+    }
+
+    /// Calculates the least nonnegative remainder of `self (mod rhs)`.
+    ///
+    /// In particular, the return value `r` satisfies `0.0 <= r < rhs.abs()` in
+    /// most cases. However, due to a floating point round-off error it can
+    /// result in `r == rhs.abs()`, violating the mathematical definition, if
+    /// `self` is much smaller than `rhs.abs()` in magnitude and `self < 0.0`.
+    /// This result is not an element of the function's codomain, but it is the
+    /// closest floating point number in the real numbers and thus fulfills the
+    /// property `self == self.div_euclid(rhs) * rhs + self.rem_euclid(rhs)`
+    /// approximately.
+    ///
+    /// # Precision
+    ///
+    /// The result of this operation is guaranteed to be the rounded
+    /// infinite-precision result.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let a: f16 = 7.0;
+    /// let b = 4.0;
+    /// assert_eq!(a.rem_euclid(b), 3.0);
+    /// assert_eq!((-a).rem_euclid(b), 1.0);
+    /// assert_eq!(a.rem_euclid(-b), 3.0);
+    /// assert_eq!((-a).rem_euclid(-b), 1.0);
+    /// // limitation due to round-off error
+    /// assert!((-f16::EPSILON).rem_euclid(3.0) != 0.0);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[doc(alias = "modulo", alias = "mod")]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn rem_euclid(self, rhs: f16) -> f16 {
+        let r = self % rhs;
+        if r < 0.0 { r + rhs.abs() } else { r }
+    }
+
+    /// Raises a number to an integer power.
+    ///
+    /// Using this function is generally faster than using `powf`.
+    /// It might have a different sequence of rounding operations than `powf`,
+    /// so the results are not guaranteed to agree.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn powi(self, n: i32) -> f16 {
+        unsafe { intrinsics::powif16(self, n) }
+    }
+
+    /// Raises a number to a floating point power.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let x = 2.0_f16;
+    /// let abs_difference = (x.powf(2.0) - (x * x)).abs();
+    ///
+    /// assert!(abs_difference <= f16::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn powf(self, n: f16) -> f16 {
+        unsafe { intrinsics::powf16(self, n) }
+    }
+
+    /// Returns the square root of a number.
+    ///
+    /// Returns NaN if `self` is a negative number other than `-0.0`.
+    ///
+    /// # Precision
+    ///
+    /// The result of this operation is guaranteed to be the rounded
+    /// infinite-precision result. It is specified by IEEE 754 as `squareRoot`
+    /// and guaranteed not to change.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let positive = 4.0_f16;
+    /// let negative = -4.0_f16;
+    /// let negative_zero = -0.0_f16;
+    ///
+    /// assert_eq!(positive.sqrt(), 2.0);
+    /// assert!(negative.sqrt().is_nan());
+    /// assert!(negative_zero.sqrt() == negative_zero);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn sqrt(self) -> f16 {
+        unsafe { intrinsics::sqrtf16(self) }
+    }
+
+    /// Returns `e^(self)`, (the exponential function).
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let one = 1.0f16;
+    /// // e^1
+    /// let e = one.exp();
+    ///
+    /// // ln(e) - 1 == 0
+    /// let abs_difference = (e.ln() - 1.0).abs();
+    ///
+    /// assert!(abs_difference <= f16::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn exp(self) -> f16 {
+        unsafe { intrinsics::expf16(self) }
+    }
+
+    /// Returns `2^(self)`.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let f = 2.0f16;
+    ///
+    /// // 2^2 - 4 == 0
+    /// let abs_difference = (f.exp2() - 4.0).abs();
+    ///
+    /// assert!(abs_difference <= f16::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn exp2(self) -> f16 {
+        unsafe { intrinsics::exp2f16(self) }
+    }
+
+    /// Returns the natural logarithm of the number.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let one = 1.0f16;
+    /// // e^1
+    /// let e = one.exp();
+    ///
+    /// // ln(e) - 1 == 0
+    /// let abs_difference = (e.ln() - 1.0).abs();
+    ///
+    /// assert!(abs_difference <= f16::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn ln(self) -> f16 {
+        unsafe { intrinsics::logf16(self) }
+    }
+
+    /// Returns the logarithm of the number with respect to an arbitrary base.
+    ///
+    /// The result might not be correctly rounded owing to implementation details;
+    /// `self.log2()` can produce more accurate results for base 2, and
+    /// `self.log10()` can produce more accurate results for base 10.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let five = 5.0f16;
+    ///
+    /// // log5(5) - 1 == 0
+    /// let abs_difference = (five.log(5.0) - 1.0).abs();
+    ///
+    /// assert!(abs_difference <= f16::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn log(self, base: f16) -> f16 {
+        self.ln() / base.ln()
+    }
+
+    /// Returns the base 2 logarithm of the number.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let two = 2.0f16;
+    ///
+    /// // log2(2) - 1 == 0
+    /// let abs_difference = (two.log2() - 1.0).abs();
+    ///
+    /// assert!(abs_difference <= f16::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn log2(self) -> f16 {
+        unsafe { intrinsics::log2f16(self) }
+    }
+
+    /// Returns the base 10 logarithm of the number.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let ten = 10.0f16;
+    ///
+    /// // log10(10) - 1 == 0
+    /// let abs_difference = (ten.log10() - 1.0).abs();
+    ///
+    /// assert!(abs_difference <= f16::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn log10(self) -> f16 {
+        unsafe { intrinsics::log10f16(self) }
+    }
+
+    /// Returns the cube root of a number.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `cbrtf` from libc on Unix
+    /// and Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let x = 8.0f16;
+    ///
+    /// // x^(1/3) - 2 == 0
+    /// let abs_difference = (x.cbrt() - 2.0).abs();
+    ///
+    /// assert!(abs_difference <= f16::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn cbrt(self) -> f16 {
+        (unsafe { cmath::cbrtf(self as f32) }) as f16
+    }
+
+    /// Compute the distance between the origin and a point (`x`, `y`) on the
+    /// Euclidean plane. Equivalently, compute the length of the hypotenuse of a
+    /// right-angle triangle with other sides having length `x.abs()` and
+    /// `y.abs()`.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `hypotf` from libc on Unix
+    /// and Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let x = 2.0f16;
+    /// let y = 3.0f16;
+    ///
+    /// // sqrt(x^2 + y^2)
+    /// let abs_difference = (x.hypot(y) - (x.powi(2) + y.powi(2)).sqrt()).abs();
+    ///
+    /// assert!(abs_difference <= f16::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn hypot(self, other: f16) -> f16 {
+        (unsafe { cmath::hypotf(self as f32, other as f32) }) as f16
+    }
+
+    /// Computes the sine of a number (in radians).
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let x = std::f16::consts::FRAC_PI_2;
+    ///
+    /// let abs_difference = (x.sin() - 1.0).abs();
+    ///
+    /// assert!(abs_difference <= f16::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn sin(self) -> f16 {
+        unsafe { intrinsics::sinf16(self) }
+    }
+
+    /// Computes the cosine of a number (in radians).
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let x = 2.0 * std::f16::consts::PI;
+    ///
+    /// let abs_difference = (x.cos() - 1.0).abs();
+    ///
+    /// assert!(abs_difference <= f16::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn cos(self) -> f16 {
+        unsafe { intrinsics::cosf16(self) }
+    }
+
+    /// Computes the tangent of a number (in radians).
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `tanf` from libc on Unix and
+    /// Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let x = std::f16::consts::FRAC_PI_4;
+    /// let abs_difference = (x.tan() - 1.0).abs();
+    ///
+    /// assert!(abs_difference <= f16::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn tan(self) -> f16 {
+        (unsafe { cmath::tanf(self as f32) }) as f16
+    }
+
+    /// Computes the arcsine of a number. Return value is in radians in
+    /// the range [-pi/2, pi/2] or NaN if the number is outside the range
+    /// [-1, 1].
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `asinf` from libc on Unix
+    /// and Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let f = std::f16::consts::FRAC_PI_2;
+    ///
+    /// // asin(sin(pi/2))
+    /// let abs_difference = (f.sin().asin() - std::f16::consts::FRAC_PI_2).abs();
+    ///
+    /// assert!(abs_difference <= f16::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[doc(alias = "arcsin")]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn asin(self) -> f16 {
+        (unsafe { cmath::asinf(self as f32) }) as f16
+    }
+
+    /// Computes the arccosine of a number. Return value is in radians in
+    /// the range [0, pi] or NaN if the number is outside the range
+    /// [-1, 1].
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `acosf` from libc on Unix
+    /// and Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let f = std::f16::consts::FRAC_PI_4;
+    ///
+    /// // acos(cos(pi/4))
+    /// let abs_difference = (f.cos().acos() - std::f16::consts::FRAC_PI_4).abs();
+    ///
+    /// assert!(abs_difference <= f16::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[doc(alias = "arccos")]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn acos(self) -> f16 {
+        (unsafe { cmath::acosf(self as f32) }) as f16
+    }
+
+    /// Computes the arctangent of a number. Return value is in radians in the
+    /// range [-pi/2, pi/2];
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `atanf` from libc on Unix
+    /// and Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let f = 1.0f16;
+    ///
+    /// // atan(tan(1))
+    /// let abs_difference = (f.tan().atan() - 1.0).abs();
+    ///
+    /// assert!(abs_difference <= f16::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[doc(alias = "arctan")]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn atan(self) -> f16 {
+        (unsafe { cmath::atanf(self as f32) }) as f16
+    }
+
+    /// Computes the four quadrant arctangent of `self` (`y`) and `other` (`x`) in radians.
+    ///
+    /// * `x = 0`, `y = 0`: `0`
+    /// * `x >= 0`: `arctan(y/x)` -> `[-pi/2, pi/2]`
+    /// * `y >= 0`: `arctan(y/x) + pi` -> `(pi/2, pi]`
+    /// * `y < 0`: `arctan(y/x) - pi` -> `(-pi, -pi/2)`
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `atan2f` from libc on Unix
+    /// and Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// // Positive angles measured counter-clockwise
+    /// // from positive x axis
+    /// // -pi/4 radians (45 deg clockwise)
+    /// let x1 = 3.0f16;
+    /// let y1 = -3.0f16;
+    ///
+    /// // 3pi/4 radians (135 deg counter-clockwise)
+    /// let x2 = -3.0f16;
+    /// let y2 = 3.0f16;
+    ///
+    /// let abs_difference_1 = (y1.atan2(x1) - (-std::f16::consts::FRAC_PI_4)).abs();
+    /// let abs_difference_2 = (y2.atan2(x2) - (3.0 * std::f16::consts::FRAC_PI_4)).abs();
+    ///
+    /// assert!(abs_difference_1 <= f16::EPSILON);
+    /// assert!(abs_difference_2 <= f16::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn atan2(self, other: f16) -> f16 {
+        (unsafe { cmath::atan2f(self as f32, other as f32) }) as f16
+    }
+
+    /// Simultaneously computes the sine and cosine of the number, `x`. Returns
+    /// `(sin(x), cos(x))`.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `(f16::sin(x),
+    /// f16::cos(x))`. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let x = std::f16::consts::FRAC_PI_4;
+    /// let f = x.sin_cos();
+    ///
+    /// let abs_difference_0 = (f.0 - x.sin()).abs();
+    /// let abs_difference_1 = (f.1 - x.cos()).abs();
+    ///
+    /// assert!(abs_difference_0 <= f16::EPSILON);
+    /// assert!(abs_difference_1 <= f16::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[doc(alias = "sincos")]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    pub fn sin_cos(self) -> (f16, f16) {
+        (self.sin(), self.cos())
+    }
+
+    /// Returns `e^(self) - 1` in a way that is accurate even if the
+    /// number is close to zero.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `expm1f` from libc on Unix
+    /// and Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let x = 1e-4_f16;
+    ///
+    /// // for very small x, e^x is approximately 1 + x + x^2 / 2
+    /// let approx = x + x * x / 2.0;
+    /// let abs_difference = (x.exp_m1() - approx).abs();
+    ///
+    /// assert!(abs_difference < 1e-4);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn exp_m1(self) -> f16 {
+        (unsafe { cmath::expm1f(self as f32) }) as f16
+    }
+
+    /// Returns `ln(1+n)` (natural logarithm) more accurately than if
+    /// the operations were performed separately.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `log1pf` from libc on Unix
+    /// and Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let x = 1e-4_f16;
+    ///
+    /// // for very small x, ln(1 + x) is approximately x - x^2 / 2
+    /// let approx = x - x * x / 2.0;
+    /// let abs_difference = (x.ln_1p() - approx).abs();
+    ///
+    /// assert!(abs_difference < 1e-4);
+    /// # }
+    /// ```
+    #[inline]
+    #[doc(alias = "log1p")]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn ln_1p(self) -> f16 {
+        (unsafe { cmath::log1pf(self as f32) }) as f16
+    }
+
+    /// Hyperbolic sine function.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `sinhf` from libc on Unix
+    /// and Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let e = std::f16::consts::E;
+    /// let x = 1.0f16;
+    ///
+    /// let f = x.sinh();
+    /// // Solving sinh() at 1 gives `(e^2-1)/(2e)`
+    /// let g = ((e * e) - 1.0) / (2.0 * e);
+    /// let abs_difference = (f - g).abs();
+    ///
+    /// assert!(abs_difference <= f16::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn sinh(self) -> f16 {
+        (unsafe { cmath::sinhf(self as f32) }) as f16
+    }
+
+    /// Hyperbolic cosine function.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `coshf` from libc on Unix
+    /// and Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let e = std::f16::consts::E;
+    /// let x = 1.0f16;
+    /// let f = x.cosh();
+    /// // Solving cosh() at 1 gives this result
+    /// let g = ((e * e) + 1.0) / (2.0 * e);
+    /// let abs_difference = (f - g).abs();
+    ///
+    /// // Same result
+    /// assert!(abs_difference <= f16::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn cosh(self) -> f16 {
+        (unsafe { cmath::coshf(self as f32) }) as f16
+    }
+
+    /// Hyperbolic tangent function.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `tanhf` from libc on Unix
+    /// and Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let e = std::f16::consts::E;
+    /// let x = 1.0f16;
+    ///
+    /// let f = x.tanh();
+    /// // Solving tanh() at 1 gives `(1 - e^(-2))/(1 + e^(-2))`
+    /// let g = (1.0 - e.powi(-2)) / (1.0 + e.powi(-2));
+    /// let abs_difference = (f - g).abs();
+    ///
+    /// assert!(abs_difference <= f16::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn tanh(self) -> f16 {
+        (unsafe { cmath::tanhf(self as f32) }) as f16
+    }
+
+    /// Inverse hyperbolic sine function.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let x = 1.0f16;
+    /// let f = x.sinh().asinh();
+    ///
+    /// let abs_difference = (f - x).abs();
+    ///
+    /// assert!(abs_difference <= f16::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[doc(alias = "arcsinh")]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn asinh(self) -> f16 {
+        let ax = self.abs();
+        let ix = 1.0 / ax;
+        (ax + (ax / (Self::hypot(1.0, ix) + ix))).ln_1p().copysign(self)
+    }
+
+    /// Inverse hyperbolic cosine function.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let x = 1.0f16;
+    /// let f = x.cosh().acosh();
+    ///
+    /// let abs_difference = (f - x).abs();
+    ///
+    /// assert!(abs_difference <= f16::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[doc(alias = "arccosh")]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn acosh(self) -> f16 {
+        if self < 1.0 {
+            Self::NAN
+        } else {
+            (self + ((self - 1.0).sqrt() * (self + 1.0).sqrt())).ln()
+        }
+    }
+
+    /// Inverse hyperbolic tangent function.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let e = std::f16::consts::E;
+    /// let f = e.tanh().atanh();
+    ///
+    /// let abs_difference = (f - e).abs();
+    ///
+    /// assert!(abs_difference <= 0.01);
+    /// # }
+    /// ```
+    #[inline]
+    #[doc(alias = "arctanh")]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn atanh(self) -> f16 {
+        0.5 * ((2.0 * self) / (1.0 - self)).ln_1p()
+    }
+
+    /// Gamma function.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `tgammaf` from libc on Unix
+    /// and Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// #![feature(float_gamma)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let x = 5.0f16;
+    ///
+    /// let abs_difference = (x.gamma() - 24.0).abs();
+    ///
+    /// assert!(abs_difference <= f16::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn gamma(self) -> f16 {
+        (unsafe { cmath::tgammaf(self as f32) }) as f16
+    }
+
+    /// Natural logarithm of the absolute value of the gamma function
+    ///
+    /// The integer part of the tuple indicates the sign of the gamma function.
+    ///
+    /// # Unspecified precision
+    ///
+    /// The precision of this function is non-deterministic. This means it varies by platform,
+    /// Rust version, and can even differ within the same execution from one invocation to the next.
+    ///
+    /// This function currently corresponds to the `lgamma_r` from libc on Unix
+    /// and Windows. Note that this might change in the future.
+    ///
+    /// # Examples
+    ///
+    /// ```
+    /// #![feature(f16)]
+    /// #![feature(float_gamma)]
+    /// # #[cfg(reliable_f16_math)] {
+    ///
+    /// let x = 2.0f16;
+    ///
+    /// let abs_difference = (x.ln_gamma().0 - 0.0).abs();
+    ///
+    /// assert!(abs_difference <= f16::EPSILON);
+    /// # }
+    /// ```
+    #[inline]
+    #[rustc_allow_incoherent_impl]
+    #[unstable(feature = "f16", issue = "116909")]
+    #[must_use = "method returns a new number and does not mutate the original value"]
+    pub fn ln_gamma(self) -> (f16, i32) {
+        let mut signgamp: i32 = 0;
+        let x = (unsafe { cmath::lgammaf_r(self as f32, &mut signgamp) }) as f16;
+        (x, signgamp)
+    }
 }
diff --git a/library/std/src/f16/tests.rs b/library/std/src/f16/tests.rs
index f73bdf68e82..f0ef807dac1 100644
--- a/library/std/src/f16/tests.rs
+++ b/library/std/src/f16/tests.rs
@@ -4,11 +4,21 @@
 use crate::f16::consts;
 use crate::num::{FpCategory as Fp, *};
 
-// We run out of precision pretty quickly with f16
-// const F16_APPROX_L1: f16 = 0.001;
-const F16_APPROX_L2: f16 = 0.01;
-// const F16_APPROX_L3: f16 = 0.1;
-const F16_APPROX_L4: f16 = 0.5;
+/// Tolerance for results on the order of 10.0e-2;
+#[cfg(reliable_f16_math)]
+const TOL_N2: f16 = 0.0001;
+
+/// Tolerance for results on the order of 10.0e+0
+#[cfg(reliable_f16_math)]
+const TOL_0: f16 = 0.01;
+
+/// Tolerance for results on the order of 10.0e+2
+#[cfg(reliable_f16_math)]
+const TOL_P2: f16 = 0.5;
+
+/// Tolerance for results on the order of 10.0e+4
+#[cfg(reliable_f16_math)]
+const TOL_P4: f16 = 10.0;
 
 /// Smallest number
 const TINY_BITS: u16 = 0x1;
@@ -197,9 +207,100 @@ fn test_classify() {
     assert_eq!(1e-5f16.classify(), Fp::Subnormal);
 }
 
-// FIXME(f16_f128): add missing math functions when available
+#[test]
+#[cfg(reliable_f16_math)]
+fn test_floor() {
+    assert_approx_eq!(1.0f16.floor(), 1.0f16, TOL_0);
+    assert_approx_eq!(1.3f16.floor(), 1.0f16, TOL_0);
+    assert_approx_eq!(1.5f16.floor(), 1.0f16, TOL_0);
+    assert_approx_eq!(1.7f16.floor(), 1.0f16, TOL_0);
+    assert_approx_eq!(0.0f16.floor(), 0.0f16, TOL_0);
+    assert_approx_eq!((-0.0f16).floor(), -0.0f16, TOL_0);
+    assert_approx_eq!((-1.0f16).floor(), -1.0f16, TOL_0);
+    assert_approx_eq!((-1.3f16).floor(), -2.0f16, TOL_0);
+    assert_approx_eq!((-1.5f16).floor(), -2.0f16, TOL_0);
+    assert_approx_eq!((-1.7f16).floor(), -2.0f16, TOL_0);
+}
+
+#[test]
+#[cfg(reliable_f16_math)]
+fn test_ceil() {
+    assert_approx_eq!(1.0f16.ceil(), 1.0f16, TOL_0);
+    assert_approx_eq!(1.3f16.ceil(), 2.0f16, TOL_0);
+    assert_approx_eq!(1.5f16.ceil(), 2.0f16, TOL_0);
+    assert_approx_eq!(1.7f16.ceil(), 2.0f16, TOL_0);
+    assert_approx_eq!(0.0f16.ceil(), 0.0f16, TOL_0);
+    assert_approx_eq!((-0.0f16).ceil(), -0.0f16, TOL_0);
+    assert_approx_eq!((-1.0f16).ceil(), -1.0f16, TOL_0);
+    assert_approx_eq!((-1.3f16).ceil(), -1.0f16, TOL_0);
+    assert_approx_eq!((-1.5f16).ceil(), -1.0f16, TOL_0);
+    assert_approx_eq!((-1.7f16).ceil(), -1.0f16, TOL_0);
+}
+
+#[test]
+#[cfg(reliable_f16_math)]
+fn test_round() {
+    assert_approx_eq!(2.5f16.round(), 3.0f16, TOL_0);
+    assert_approx_eq!(1.0f16.round(), 1.0f16, TOL_0);
+    assert_approx_eq!(1.3f16.round(), 1.0f16, TOL_0);
+    assert_approx_eq!(1.5f16.round(), 2.0f16, TOL_0);
+    assert_approx_eq!(1.7f16.round(), 2.0f16, TOL_0);
+    assert_approx_eq!(0.0f16.round(), 0.0f16, TOL_0);
+    assert_approx_eq!((-0.0f16).round(), -0.0f16, TOL_0);
+    assert_approx_eq!((-1.0f16).round(), -1.0f16, TOL_0);
+    assert_approx_eq!((-1.3f16).round(), -1.0f16, TOL_0);
+    assert_approx_eq!((-1.5f16).round(), -2.0f16, TOL_0);
+    assert_approx_eq!((-1.7f16).round(), -2.0f16, TOL_0);
+}
+
+#[test]
+#[cfg(reliable_f16_math)]
+fn test_round_ties_even() {
+    assert_approx_eq!(2.5f16.round_ties_even(), 2.0f16, TOL_0);
+    assert_approx_eq!(1.0f16.round_ties_even(), 1.0f16, TOL_0);
+    assert_approx_eq!(1.3f16.round_ties_even(), 1.0f16, TOL_0);
+    assert_approx_eq!(1.5f16.round_ties_even(), 2.0f16, TOL_0);
+    assert_approx_eq!(1.7f16.round_ties_even(), 2.0f16, TOL_0);
+    assert_approx_eq!(0.0f16.round_ties_even(), 0.0f16, TOL_0);
+    assert_approx_eq!((-0.0f16).round_ties_even(), -0.0f16, TOL_0);
+    assert_approx_eq!((-1.0f16).round_ties_even(), -1.0f16, TOL_0);
+    assert_approx_eq!((-1.3f16).round_ties_even(), -1.0f16, TOL_0);
+    assert_approx_eq!((-1.5f16).round_ties_even(), -2.0f16, TOL_0);
+    assert_approx_eq!((-1.7f16).round_ties_even(), -2.0f16, TOL_0);
+}
 
 #[test]
+#[cfg(reliable_f16_math)]
+fn test_trunc() {
+    assert_approx_eq!(1.0f16.trunc(), 1.0f16, TOL_0);
+    assert_approx_eq!(1.3f16.trunc(), 1.0f16, TOL_0);
+    assert_approx_eq!(1.5f16.trunc(), 1.0f16, TOL_0);
+    assert_approx_eq!(1.7f16.trunc(), 1.0f16, TOL_0);
+    assert_approx_eq!(0.0f16.trunc(), 0.0f16, TOL_0);
+    assert_approx_eq!((-0.0f16).trunc(), -0.0f16, TOL_0);
+    assert_approx_eq!((-1.0f16).trunc(), -1.0f16, TOL_0);
+    assert_approx_eq!((-1.3f16).trunc(), -1.0f16, TOL_0);
+    assert_approx_eq!((-1.5f16).trunc(), -1.0f16, TOL_0);
+    assert_approx_eq!((-1.7f16).trunc(), -1.0f16, TOL_0);
+}
+
+#[test]
+#[cfg(reliable_f16_math)]
+fn test_fract() {
+    assert_approx_eq!(1.0f16.fract(), 0.0f16, TOL_0);
+    assert_approx_eq!(1.3f16.fract(), 0.3f16, TOL_0);
+    assert_approx_eq!(1.5f16.fract(), 0.5f16, TOL_0);
+    assert_approx_eq!(1.7f16.fract(), 0.7f16, TOL_0);
+    assert_approx_eq!(0.0f16.fract(), 0.0f16, TOL_0);
+    assert_approx_eq!((-0.0f16).fract(), -0.0f16, TOL_0);
+    assert_approx_eq!((-1.0f16).fract(), -0.0f16, TOL_0);
+    assert_approx_eq!((-1.3f16).fract(), -0.3f16, TOL_0);
+    assert_approx_eq!((-1.5f16).fract(), -0.5f16, TOL_0);
+    assert_approx_eq!((-1.7f16).fract(), -0.7f16, TOL_0);
+}
+
+#[test]
+#[cfg(reliable_f16_math)]
 fn test_abs() {
     assert_eq!(f16::INFINITY.abs(), f16::INFINITY);
     assert_eq!(1f16.abs(), 1f16);
@@ -299,6 +400,24 @@ fn test_next_down() {
 }
 
 #[test]
+#[cfg(reliable_f16_math)]
+fn test_mul_add() {
+    let nan: f16 = f16::NAN;
+    let inf: f16 = f16::INFINITY;
+    let neg_inf: f16 = f16::NEG_INFINITY;
+    assert_approx_eq!(12.3f16.mul_add(4.5, 6.7), 62.05, TOL_P2);
+    assert_approx_eq!((-12.3f16).mul_add(-4.5, -6.7), 48.65, TOL_P2);
+    assert_approx_eq!(0.0f16.mul_add(8.9, 1.2), 1.2, TOL_0);
+    assert_approx_eq!(3.4f16.mul_add(-0.0, 5.6), 5.6, TOL_0);
+    assert!(nan.mul_add(7.8, 9.0).is_nan());
+    assert_eq!(inf.mul_add(7.8, 9.0), inf);
+    assert_eq!(neg_inf.mul_add(7.8, 9.0), neg_inf);
+    assert_eq!(8.9f16.mul_add(inf, 3.2), inf);
+    assert_eq!((-3.2f16).mul_add(2.4, neg_inf), neg_inf);
+}
+
+#[test]
+#[cfg(reliable_f16_math)]
 fn test_recip() {
     let nan: f16 = f16::NAN;
     let inf: f16 = f16::INFINITY;
@@ -307,20 +426,166 @@ fn test_recip() {
     assert_eq!(2.0f16.recip(), 0.5);
     assert_eq!((-0.4f16).recip(), -2.5);
     assert_eq!(0.0f16.recip(), inf);
+    assert_approx_eq!(f16::MAX.recip(), 1.526624e-5f16, 1e-4);
     assert!(nan.recip().is_nan());
     assert_eq!(inf.recip(), 0.0);
     assert_eq!(neg_inf.recip(), 0.0);
 }
 
 #[test]
+#[cfg(reliable_f16_math)]
+fn test_powi() {
+    // FIXME(llvm19): LLVM misoptimizes `powi.f16`
+    // <https://github.com/llvm/llvm-project/issues/98665>
+    // let nan: f16 = f16::NAN;
+    // let inf: f16 = f16::INFINITY;
+    // let neg_inf: f16 = f16::NEG_INFINITY;
+    // assert_eq!(1.0f16.powi(1), 1.0);
+    // assert_approx_eq!((-3.1f16).powi(2), 9.61, TOL_0);
+    // assert_approx_eq!(5.9f16.powi(-2), 0.028727, TOL_N2);
+    // assert_eq!(8.3f16.powi(0), 1.0);
+    // assert!(nan.powi(2).is_nan());
+    // assert_eq!(inf.powi(3), inf);
+    // assert_eq!(neg_inf.powi(2), inf);
+}
+
+#[test]
+#[cfg(reliable_f16_math)]
+fn test_powf() {
+    let nan: f16 = f16::NAN;
+    let inf: f16 = f16::INFINITY;
+    let neg_inf: f16 = f16::NEG_INFINITY;
+    assert_eq!(1.0f16.powf(1.0), 1.0);
+    assert_approx_eq!(3.4f16.powf(4.5), 246.408183, TOL_P2);
+    assert_approx_eq!(2.7f16.powf(-3.2), 0.041652, TOL_N2);
+    assert_approx_eq!((-3.1f16).powf(2.0), 9.61, TOL_P2);
+    assert_approx_eq!(5.9f16.powf(-2.0), 0.028727, TOL_N2);
+    assert_eq!(8.3f16.powf(0.0), 1.0);
+    assert!(nan.powf(2.0).is_nan());
+    assert_eq!(inf.powf(2.0), inf);
+    assert_eq!(neg_inf.powf(3.0), neg_inf);
+}
+
+#[test]
+#[cfg(reliable_f16_math)]
+fn test_sqrt_domain() {
+    assert!(f16::NAN.sqrt().is_nan());
+    assert!(f16::NEG_INFINITY.sqrt().is_nan());
+    assert!((-1.0f16).sqrt().is_nan());
+    assert_eq!((-0.0f16).sqrt(), -0.0);
+    assert_eq!(0.0f16.sqrt(), 0.0);
+    assert_eq!(1.0f16.sqrt(), 1.0);
+    assert_eq!(f16::INFINITY.sqrt(), f16::INFINITY);
+}
+
+#[test]
+#[cfg(reliable_f16_math)]
+fn test_exp() {
+    assert_eq!(1.0, 0.0f16.exp());
+    assert_approx_eq!(2.718282, 1.0f16.exp(), TOL_0);
+    assert_approx_eq!(148.413159, 5.0f16.exp(), TOL_0);
+
+    let inf: f16 = f16::INFINITY;
+    let neg_inf: f16 = f16::NEG_INFINITY;
+    let nan: f16 = f16::NAN;
+    assert_eq!(inf, inf.exp());
+    assert_eq!(0.0, neg_inf.exp());
+    assert!(nan.exp().is_nan());
+}
+
+#[test]
+#[cfg(reliable_f16_math)]
+fn test_exp2() {
+    assert_eq!(32.0, 5.0f16.exp2());
+    assert_eq!(1.0, 0.0f16.exp2());
+
+    let inf: f16 = f16::INFINITY;
+    let neg_inf: f16 = f16::NEG_INFINITY;
+    let nan: f16 = f16::NAN;
+    assert_eq!(inf, inf.exp2());
+    assert_eq!(0.0, neg_inf.exp2());
+    assert!(nan.exp2().is_nan());
+}
+
+#[test]
+#[cfg(reliable_f16_math)]
+fn test_ln() {
+    let nan: f16 = f16::NAN;
+    let inf: f16 = f16::INFINITY;
+    let neg_inf: f16 = f16::NEG_INFINITY;
+    assert_approx_eq!(1.0f16.exp().ln(), 1.0, TOL_0);
+    assert!(nan.ln().is_nan());
+    assert_eq!(inf.ln(), inf);
+    assert!(neg_inf.ln().is_nan());
+    assert!((-2.3f16).ln().is_nan());
+    assert_eq!((-0.0f16).ln(), neg_inf);
+    assert_eq!(0.0f16.ln(), neg_inf);
+    assert_approx_eq!(4.0f16.ln(), 1.386294, TOL_0);
+}
+
+#[test]
+#[cfg(reliable_f16_math)]
+fn test_log() {
+    let nan: f16 = f16::NAN;
+    let inf: f16 = f16::INFINITY;
+    let neg_inf: f16 = f16::NEG_INFINITY;
+    assert_eq!(10.0f16.log(10.0), 1.0);
+    assert_approx_eq!(2.3f16.log(3.5), 0.664858, TOL_0);
+    assert_eq!(1.0f16.exp().log(1.0f16.exp()), 1.0);
+    assert!(1.0f16.log(1.0).is_nan());
+    assert!(1.0f16.log(-13.9).is_nan());
+    assert!(nan.log(2.3).is_nan());
+    assert_eq!(inf.log(10.0), inf);
+    assert!(neg_inf.log(8.8).is_nan());
+    assert!((-2.3f16).log(0.1).is_nan());
+    assert_eq!((-0.0f16).log(2.0), neg_inf);
+    assert_eq!(0.0f16.log(7.0), neg_inf);
+}
+
+#[test]
+#[cfg(reliable_f16_math)]
+fn test_log2() {
+    let nan: f16 = f16::NAN;
+    let inf: f16 = f16::INFINITY;
+    let neg_inf: f16 = f16::NEG_INFINITY;
+    assert_approx_eq!(10.0f16.log2(), 3.321928, TOL_0);
+    assert_approx_eq!(2.3f16.log2(), 1.201634, TOL_0);
+    assert_approx_eq!(1.0f16.exp().log2(), 1.442695, TOL_0);
+    assert!(nan.log2().is_nan());
+    assert_eq!(inf.log2(), inf);
+    assert!(neg_inf.log2().is_nan());
+    assert!((-2.3f16).log2().is_nan());
+    assert_eq!((-0.0f16).log2(), neg_inf);
+    assert_eq!(0.0f16.log2(), neg_inf);
+}
+
+#[test]
+#[cfg(reliable_f16_math)]
+fn test_log10() {
+    let nan: f16 = f16::NAN;
+    let inf: f16 = f16::INFINITY;
+    let neg_inf: f16 = f16::NEG_INFINITY;
+    assert_eq!(10.0f16.log10(), 1.0);
+    assert_approx_eq!(2.3f16.log10(), 0.361728, TOL_0);
+    assert_approx_eq!(1.0f16.exp().log10(), 0.434294, TOL_0);
+    assert_eq!(1.0f16.log10(), 0.0);
+    assert!(nan.log10().is_nan());
+    assert_eq!(inf.log10(), inf);
+    assert!(neg_inf.log10().is_nan());
+    assert!((-2.3f16).log10().is_nan());
+    assert_eq!((-0.0f16).log10(), neg_inf);
+    assert_eq!(0.0f16.log10(), neg_inf);
+}
+
+#[test]
 fn test_to_degrees() {
     let pi: f16 = consts::PI;
     let nan: f16 = f16::NAN;
     let inf: f16 = f16::INFINITY;
     let neg_inf: f16 = f16::NEG_INFINITY;
     assert_eq!(0.0f16.to_degrees(), 0.0);
-    assert_approx_eq!((-5.8f16).to_degrees(), -332.315521);
-    assert_approx_eq!(pi.to_degrees(), 180.0, F16_APPROX_L4);
+    assert_approx_eq!((-5.8f16).to_degrees(), -332.315521, TOL_P2);
+    assert_approx_eq!(pi.to_degrees(), 180.0, TOL_P2);
     assert!(nan.to_degrees().is_nan());
     assert_eq!(inf.to_degrees(), inf);
     assert_eq!(neg_inf.to_degrees(), neg_inf);
@@ -334,15 +599,113 @@ fn test_to_radians() {
     let inf: f16 = f16::INFINITY;
     let neg_inf: f16 = f16::NEG_INFINITY;
     assert_eq!(0.0f16.to_radians(), 0.0);
-    assert_approx_eq!(154.6f16.to_radians(), 2.698279);
-    assert_approx_eq!((-332.31f16).to_radians(), -5.799903);
-    assert_approx_eq!(180.0f16.to_radians(), pi, F16_APPROX_L2);
+    assert_approx_eq!(154.6f16.to_radians(), 2.698279, TOL_0);
+    assert_approx_eq!((-332.31f16).to_radians(), -5.799903, TOL_0);
+    assert_approx_eq!(180.0f16.to_radians(), pi, TOL_0);
     assert!(nan.to_radians().is_nan());
     assert_eq!(inf.to_radians(), inf);
     assert_eq!(neg_inf.to_radians(), neg_inf);
 }
 
 #[test]
+#[cfg(reliable_f16_math)]
+fn test_asinh() {
+    assert_eq!(0.0f16.asinh(), 0.0f16);
+    assert_eq!((-0.0f16).asinh(), -0.0f16);
+
+    let inf: f16 = f16::INFINITY;
+    let neg_inf: f16 = f16::NEG_INFINITY;
+    let nan: f16 = f16::NAN;
+    assert_eq!(inf.asinh(), inf);
+    assert_eq!(neg_inf.asinh(), neg_inf);
+    assert!(nan.asinh().is_nan());
+    assert!((-0.0f16).asinh().is_sign_negative());
+    // issue 63271
+    assert_approx_eq!(2.0f16.asinh(), 1.443635475178810342493276740273105f16, TOL_0);
+    assert_approx_eq!((-2.0f16).asinh(), -1.443635475178810342493276740273105f16, TOL_0);
+    // regression test for the catastrophic cancellation fixed in 72486
+    assert_approx_eq!((-200.0f16).asinh(), -5.991470797049389, TOL_0);
+
+    // test for low accuracy from issue 104548
+    assert_approx_eq!(10.0f16, 10.0f16.sinh().asinh(), TOL_0);
+    // mul needed for approximate comparison to be meaningful
+    assert_approx_eq!(1.0f16, 1e-3f16.sinh().asinh() * 1e3f16, TOL_0);
+}
+
+#[test]
+#[cfg(reliable_f16_math)]
+fn test_acosh() {
+    assert_eq!(1.0f16.acosh(), 0.0f16);
+    assert!(0.999f16.acosh().is_nan());
+
+    let inf: f16 = f16::INFINITY;
+    let neg_inf: f16 = f16::NEG_INFINITY;
+    let nan: f16 = f16::NAN;
+    assert_eq!(inf.acosh(), inf);
+    assert!(neg_inf.acosh().is_nan());
+    assert!(nan.acosh().is_nan());
+    assert_approx_eq!(2.0f16.acosh(), 1.31695789692481670862504634730796844f16, TOL_0);
+    assert_approx_eq!(3.0f16.acosh(), 1.76274717403908605046521864995958461f16, TOL_0);
+
+    // test for low accuracy from issue 104548
+    assert_approx_eq!(10.0f16, 10.0f16.cosh().acosh(), TOL_P2);
+}
+
+#[test]
+#[cfg(reliable_f16_math)]
+fn test_atanh() {
+    assert_eq!(0.0f16.atanh(), 0.0f16);
+    assert_eq!((-0.0f16).atanh(), -0.0f16);
+
+    let inf: f16 = f16::INFINITY;
+    let neg_inf: f16 = f16::NEG_INFINITY;
+    let nan: f16 = f16::NAN;
+    assert_eq!(1.0f16.atanh(), inf);
+    assert_eq!((-1.0f16).atanh(), neg_inf);
+    assert!(2f16.atanh().atanh().is_nan());
+    assert!((-2f16).atanh().atanh().is_nan());
+    assert!(inf.atanh().is_nan());
+    assert!(neg_inf.atanh().is_nan());
+    assert!(nan.atanh().is_nan());
+    assert_approx_eq!(0.5f16.atanh(), 0.54930614433405484569762261846126285f16, TOL_0);
+    assert_approx_eq!((-0.5f16).atanh(), -0.54930614433405484569762261846126285f16, TOL_0);
+}
+
+#[test]
+#[cfg(reliable_f16_math)]
+fn test_gamma() {
+    // precision can differ among platforms
+    assert_approx_eq!(1.0f16.gamma(), 1.0f16, TOL_0);
+    assert_approx_eq!(2.0f16.gamma(), 1.0f16, TOL_0);
+    assert_approx_eq!(3.0f16.gamma(), 2.0f16, TOL_0);
+    assert_approx_eq!(4.0f16.gamma(), 6.0f16, TOL_0);
+    assert_approx_eq!(5.0f16.gamma(), 24.0f16, TOL_0);
+    assert_approx_eq!(0.5f16.gamma(), consts::PI.sqrt(), TOL_0);
+    assert_approx_eq!((-0.5f16).gamma(), -2.0 * consts::PI.sqrt(), TOL_0);
+    assert_eq!(0.0f16.gamma(), f16::INFINITY);
+    assert_eq!((-0.0f16).gamma(), f16::NEG_INFINITY);
+    assert!((-1.0f16).gamma().is_nan());
+    assert!((-2.0f16).gamma().is_nan());
+    assert!(f16::NAN.gamma().is_nan());
+    assert!(f16::NEG_INFINITY.gamma().is_nan());
+    assert_eq!(f16::INFINITY.gamma(), f16::INFINITY);
+    assert_eq!(171.71f16.gamma(), f16::INFINITY);
+}
+
+#[test]
+#[cfg(reliable_f16_math)]
+fn test_ln_gamma() {
+    assert_approx_eq!(1.0f16.ln_gamma().0, 0.0f16, TOL_0);
+    assert_eq!(1.0f16.ln_gamma().1, 1);
+    assert_approx_eq!(2.0f16.ln_gamma().0, 0.0f16, TOL_0);
+    assert_eq!(2.0f16.ln_gamma().1, 1);
+    assert_approx_eq!(3.0f16.ln_gamma().0, 2.0f16.ln(), TOL_0);
+    assert_eq!(3.0f16.ln_gamma().1, 1);
+    assert_approx_eq!((-0.5f16).ln_gamma().0, (2.0 * consts::PI.sqrt()).ln(), TOL_0);
+    assert_eq!((-0.5f16).ln_gamma().1, -1);
+}
+
+#[test]
 fn test_real_consts() {
     // FIXME(f16_f128): add math tests when available
     use super::consts;
@@ -355,29 +718,34 @@ fn test_real_consts() {
     let frac_pi_8: f16 = consts::FRAC_PI_8;
     let frac_1_pi: f16 = consts::FRAC_1_PI;
     let frac_2_pi: f16 = consts::FRAC_2_PI;
-    // let frac_2_sqrtpi: f16 = consts::FRAC_2_SQRT_PI;
-    // let sqrt2: f16 = consts::SQRT_2;
-    // let frac_1_sqrt2: f16 = consts::FRAC_1_SQRT_2;
-    // let e: f16 = consts::E;
-    // let log2_e: f16 = consts::LOG2_E;
-    // let log10_e: f16 = consts::LOG10_E;
-    // let ln_2: f16 = consts::LN_2;
-    // let ln_10: f16 = consts::LN_10;
-
-    assert_approx_eq!(frac_pi_2, pi / 2f16);
-    assert_approx_eq!(frac_pi_3, pi / 3f16);
-    assert_approx_eq!(frac_pi_4, pi / 4f16);
-    assert_approx_eq!(frac_pi_6, pi / 6f16);
-    assert_approx_eq!(frac_pi_8, pi / 8f16);
-    assert_approx_eq!(frac_1_pi, 1f16 / pi);
-    assert_approx_eq!(frac_2_pi, 2f16 / pi);
-    // assert_approx_eq!(frac_2_sqrtpi, 2f16 / pi.sqrt());
-    // assert_approx_eq!(sqrt2, 2f16.sqrt());
-    // assert_approx_eq!(frac_1_sqrt2, 1f16 / 2f16.sqrt());
-    // assert_approx_eq!(log2_e, e.log2());
-    // assert_approx_eq!(log10_e, e.log10());
-    // assert_approx_eq!(ln_2, 2f16.ln());
-    // assert_approx_eq!(ln_10, 10f16.ln());
+
+    assert_approx_eq!(frac_pi_2, pi / 2f16, TOL_0);
+    assert_approx_eq!(frac_pi_3, pi / 3f16, TOL_0);
+    assert_approx_eq!(frac_pi_4, pi / 4f16, TOL_0);
+    assert_approx_eq!(frac_pi_6, pi / 6f16, TOL_0);
+    assert_approx_eq!(frac_pi_8, pi / 8f16, TOL_0);
+    assert_approx_eq!(frac_1_pi, 1f16 / pi, TOL_0);
+    assert_approx_eq!(frac_2_pi, 2f16 / pi, TOL_0);
+
+    #[cfg(reliable_f16_math)]
+    {
+        let frac_2_sqrtpi: f16 = consts::FRAC_2_SQRT_PI;
+        let sqrt2: f16 = consts::SQRT_2;
+        let frac_1_sqrt2: f16 = consts::FRAC_1_SQRT_2;
+        let e: f16 = consts::E;
+        let log2_e: f16 = consts::LOG2_E;
+        let log10_e: f16 = consts::LOG10_E;
+        let ln_2: f16 = consts::LN_2;
+        let ln_10: f16 = consts::LN_10;
+
+        assert_approx_eq!(frac_2_sqrtpi, 2f16 / pi.sqrt(), TOL_0);
+        assert_approx_eq!(sqrt2, 2f16.sqrt(), TOL_0);
+        assert_approx_eq!(frac_1_sqrt2, 1f16 / 2f16.sqrt(), TOL_0);
+        assert_approx_eq!(log2_e, e.log2(), TOL_0);
+        assert_approx_eq!(log10_e, e.log10(), TOL_0);
+        assert_approx_eq!(ln_2, 2f16.ln(), TOL_0);
+        assert_approx_eq!(ln_10, 10f16.ln(), TOL_0);
+    }
 }
 
 #[test]
@@ -386,10 +754,10 @@ fn test_float_bits_conv() {
     assert_eq!((12.5f16).to_bits(), 0x4a40);
     assert_eq!((1337f16).to_bits(), 0x6539);
     assert_eq!((-14.25f16).to_bits(), 0xcb20);
-    assert_approx_eq!(f16::from_bits(0x3c00), 1.0);
-    assert_approx_eq!(f16::from_bits(0x4a40), 12.5);
-    assert_approx_eq!(f16::from_bits(0x6539), 1337.0);
-    assert_approx_eq!(f16::from_bits(0xcb20), -14.25);
+    assert_approx_eq!(f16::from_bits(0x3c00), 1.0, TOL_0);
+    assert_approx_eq!(f16::from_bits(0x4a40), 12.5, TOL_0);
+    assert_approx_eq!(f16::from_bits(0x6539), 1337.0, TOL_P4);
+    assert_approx_eq!(f16::from_bits(0xcb20), -14.25, TOL_0);
 
     // Check that NaNs roundtrip their bits regardless of signaling-ness
     let masked_nan1 = f16::NAN.to_bits() ^ NAN_MASK1;
diff --git a/library/std/src/macros.rs b/library/std/src/macros.rs
index ba519afc62b..1b0d7f3dbf2 100644
--- a/library/std/src/macros.rs
+++ b/library/std/src/macros.rs
@@ -382,7 +382,7 @@ macro_rules! assert_approx_eq {
         let diff = (*a - *b).abs();
         assert!(
             diff < $lim,
-            "{a:?} is not approximately equal to {b:?} (threshold {lim:?}, actual {diff:?})",
+            "{a:?} is not approximately equal to {b:?} (threshold {lim:?}, difference {diff:?})",
             lim = $lim
         );
     }};
diff --git a/library/std/src/sys/cmath.rs b/library/std/src/sys/cmath.rs
index 99df503b82d..2997e908fa1 100644
--- a/library/std/src/sys/cmath.rs
+++ b/library/std/src/sys/cmath.rs
@@ -28,6 +28,21 @@ extern "C" {
     pub fn lgamma_r(n: f64, s: &mut i32) -> f64;
     pub fn lgammaf_r(n: f32, s: &mut i32) -> f32;
 
+    pub fn acosf128(n: f128) -> f128;
+    pub fn asinf128(n: f128) -> f128;
+    pub fn atanf128(n: f128) -> f128;
+    pub fn atan2f128(a: f128, b: f128) -> f128;
+    pub fn cbrtf128(n: f128) -> f128;
+    pub fn coshf128(n: f128) -> f128;
+    pub fn expm1f128(n: f128) -> f128;
+    pub fn hypotf128(x: f128, y: f128) -> f128;
+    pub fn log1pf128(n: f128) -> f128;
+    pub fn sinhf128(n: f128) -> f128;
+    pub fn tanf128(n: f128) -> f128;
+    pub fn tanhf128(n: f128) -> f128;
+    pub fn tgammaf128(n: f128) -> f128;
+    pub fn lgammaf128_r(n: f128, s: &mut i32) -> f128;
+
     cfg_if::cfg_if! {
     if #[cfg(not(all(target_os = "windows", target_env = "msvc", target_arch = "x86")))] {
         pub fn acosf(n: f32) -> f32;