Implement trait for `is_close` for types that have an `EPSILON`?

I was hoping to implement a generic function (or trait) to allow me to assert!(is_close(first_float, second_float)) for a few tests, specifically for types that have already defined EPSILON.

(I know there are crates like is_close - Rust and https://crates.io/crates/float_eq, but looking through their source they seem more complicated than I think should be possible for this simple case, and I wanted to see if I could do it.)

Is there a way for me to write a generic function or trait that requires that a type as an associated constant (in this case EPSILON)? If a trait, I'd also like to provide a default implementation that relies on this constant.

I was thinking about something like the following:

trait Closeable<T> {
    const EPSILON: T;
    fn is_close (lhs: T, rhs: T) -> bool where T: std::ops::Sub + PartialOrd {
        use std::cmp::Ordering::{Equal, Greater, Less};
        let ord = if let Some(ord) = lhs.partial_cmp(&rhs) {
            ord
        } else {
            panic!("Did I catch a NaN in there somewhere?")
        };
        match ord {
            Equal => true,
            Greater => (lhs - rhs) <= T::EPSILON,
            Less => (rhs - lhs) <= T::EPSILON,
        }
    }
}

impl Closeable<f32> for f32 {}

This fails to compile, complaining that EPSILON is not defined. Is there a way for me to constrain it to use the already defined EPSILON? Or is the proper way to manually const EPSILON = f32::EPSILON? for each case?

TIA for any help!

You can use the num-traits crate. It defines the FloatCore::epsilon() method.

If you want to code against an EPSILON associated constant it'll need to be part of a trait (either Closeable or num_traits::FloatCore or whatever) because Rust only lets you use traits to reason about a type (methods, associated types/constants, etc.).

This is in contrast with something like C++ where you can write the equivalent of a generic impl and if there is no T::EPSILON associated constant your impl will just be ignored (SFINAE).

trait Closeable {
    const EPSILON: Self;
    fn is_close(self, rhs: Self) -> bool
    where
        Self: std::ops::Sub<Output = Self> + PartialOrd + Sized,
    {
        use std::cmp::Ordering::{Equal, Greater, Less};
        let ord = if let Some(ord) = self.partial_cmp(&rhs) {
            ord
        } else {
            panic!("Did I catch a NaN in there somewhere?")
        };

        let diff = match ord {
            Equal => return true,
            Greater => self - rhs,
            Less => rhs - self,
        };

        if let Some(Greater) = diff.partial_cmp(&Self::EPSILON) {
            return false;
        }
        return true;
    }
}

impl Closeable for f32 {
    const EPSILON: Self = f32::EPSILON;
}

Okay, here's the first working version I could come up with. It's works when I call it as: assert!(content.duration.is_close(212.091));, (where content.duration is an f32).

The only thing I don't really care for is that Self::EPSILON isn't automatically defined, which was one of the original points of my question. Is there any way to have a default "implementation" for Self::EPSILON so I could just write:

impl Closeable for f32 {}

instead of

impl Closeable for f32 {
    const EPSILON: Self = f32::EPSILON;
}

?
const EPSILON: Self = f32::EPSILON;
}

The other crates you mention are complicated for a reason. It's often more robust to ask if two floats are within N units-of-least-precision. Comparing against epsilon will work ok-ish for small values but give unexpected answers for large values.
It might not matter for your use case but you should be aware of the minefield ahead.

5 Likes

As @lemmih suggests, this is almost certainly not the way you want to approximately compare floating point numbers. In particular, the above is equivalent to content.duration == 212.091, because no other f32 number is within f32::EPSILON of 212.091. The smallest possible nonzero difference for numbers of this size is 128.0 * f32::EPSILON.

3 Likes

Thanks for the replies!

I'm sure! Thought I assume part of the reason includes scope; I was only hoping to have something to quickly implement reasonably adequate floating-point comparison for f16/32/64, as I thought comparing them with == was generally a bad idea and a better way was using a margin of error (which is what I thought e.g. f32::EPSILON was for -- a recommended margin of error for f32 comparisons).

Huh, why is that? For large-but-still-f32 numbers? Thanks for the heads up!

Why 128? Something about bits and 2^7?

Floating point numbers aren't evenly spaced on the number line. As you move away from zero (and numbers become larger, in terms of absolute value), they become more and more sparse, which means that the distance between consecutive floating point numbers increases. This distance (also sometimes referred to as ULP, or units-in-the-last place) is (roughly) the size of the number times EPSILON. It's equal to EPSILON when the number is 1.0.

WikiPedia has a good chart here:

1 Like

That's not what it is. Unfortunately, some examples for f32 in the documentation seem to use it that way, which is unfortunate. They just got lucky. Some other examples there use a different margin (perhaps those that didn't get lucky with incorrectly using EPSILON).

f32::EPSILON is, roughly speaking, the smallest possible relative error. For large numbers that translates to larger absolute errors, for small numbers it translates to smaller absolute errors.

That's because f32 is a floating point number. It has a certain number of bits of precision (24 bits).

Imagine you have a base-10 floating point number with 5 digits of precision. So you might have something like 0.34123, or you might have something like 135.23. Both numbers have 5 digits of precision, but the latter one has only 2 digits after the decimal point, because it's larger.

A better way to check for approximate equality would be to check for relative error, and allow a larger multiple of f32::EPSILON. For instance, if you want to allow errors in the last 4 bits, you might check that your value is between x * (1 - 16.0 * f32::EPSILON) and x * (1 + 16.0 * f32::EPSILON).

Numbers between 1.0 and 2.0 have 1 bit before the binary point, and 23 bits of precision after the binary point. f32::EPSILON is 2^-23 to signify that.

Numbers between 128.0 and 256.0 have 8 bits before the binary point, and only 16 bits of precision after the binary point. Hence, their absolute accuracy is 128 times worse than for numbers between 1.0 and 2.0.

For instance, your number 212.091 is represented as 11010100.0001011101001100 in binary in an f32.

2 Likes

Wow, thank you both!

To follow up, here's what I eventually came up with for my "quick and dirty" version.

If anyone ends up here after a search, please read the helpful and informative comments above before any copy and pasting, as this is a problematic implementation that would benefit from multiplying the error margin by a variable degree.

https://play.rust-lang.org/?version=stable&mode=debug&edition=2018&gist=bde7878cb6f48eab5ed4a380e605730f

trait Close {
    const EPSILON: Self;

    fn is_close(&self, other: Self) -> bool
    where
        Self: std::ops::Sub<Output = Self> + PartialOrd + Sized + Copy,
    {
        use std::cmp::Ordering::{Equal, Greater, Less};
        let ord = if let Some(ord) = self.partial_cmp(&other) {
            ord
        } else {
            panic!("Did I catch a NaN in there somewhere?")
        };

        let diff = match ord {
            Equal => return true,
            Greater => *self - other,
            Less => other - *self,
        };

        if let Some(Greater) = diff.partial_cmp(&Self::EPSILON) {
            return false;
        }
        true
    }
}

macro_rules! close_enough {
    ($($ty:ty),*) => {
        $(
        impl Close for $ty {
            const EPSILON: $ty = <$ty>::EPSILON;
        }
        )*
    }
}

close_enough!(f32, f64);

#[test]
#[allow(clippy::float_cmp, clippy::eq_op)]
fn test_floats_are_weird() {
    assert_eq!(0.1_f32 + 0.2_f32, 0.3_f32);
    assert_ne!(0.1_f64 + 0.2_f64, 0.3_f64);
}

#[test]
fn test_is_close() {
    assert!((0.1_f32 + 0.2_f32).is_close(0.3_f32));
    assert!((0.1_f64 + 0.2_f64).is_close(0.3_f64));
}

Also, try as I might, I couldn't find any way to restrict the trait to types that defined T::EPSILON, and I couldn't find any way to give it a variable default definition. My numerous failed attempts included variations of the below that gave compiler errors about T::EPSILON being undefined, or Self::EPSILON running into recursion issues.

trait Close where Self::EPSILON // this is invalid syntax

trait Close { const EPSILON: Self = Self::EPSILON } // `<Self as Close>::EPSILON` is recursive

trait Close<T=Self> { const EPSILON: Self = <Self as T>::EPSILON } // no associated constant `EPSILON` in `T`

As you can see, I gave up and just made a macro.

Here's what I use for float compares:

/// Returns true if the given float values are nearly equal, taking into
/// account relative error and infinites.
#[allow(clippy::float_cmp)]
#[inline]
pub fn nearly_equal<T>(a: T, b: T) -> bool
    where T: Into<f64>
{
    use std::f64;
    let (a, b) = (a.into(), b.into());
    let abs_a = a.abs();
    let abs_b = b.abs();
    let diff = (a - b).abs();

    if a == b { // Shortcut, handles infinities.
        true
    } else if a == 0.0 || b == 0.0 || diff < f64::MIN_POSITIVE {
        // a or b is zero or both are extremely close to it
        // relative error is less meaningful here
        diff < (f64::EPSILON * f64::MIN_POSITIVE)
    } else { // Use relative error.
        (diff / f64::min(abs_a + abs_b, f64::MAX)) < f64::EPSILON
    }
}

A small nitpick: this is equivalent to diff == 0.0.

Your code allows up to 2 or 3 ulp difference in typical cases, but has inconsistent behavior for tiny numbers. For example, if:

    let a = 5.010420900022432e-293;
    let b = 5.010420900022433e-293;
    let c = 5.010420900022434e-293;

Then according to your criteria nearly_equal(a, c) is true, but nearly_equal(a, b) is not.

1 Like

TL;DR: comparing floats in a general way is hard. In specific cases, you should have an idea of roughly what scale is "close enough," and likely have a bound on number scale.

My personal helper crate of choice is approx. It also links to three articles I think any developer interacting with floating point should read. You don't need to fully digest every bit of minutia (especially from the Oracle one), but you should be familiar with all of it.

Comparing Floating Point Numbers, 2012 Edition | Random ASCII – tech blog of Bruce Dawson (wordpress.com)

The Floating-Point Guide - Comparison (floating-point-gui.de)

What Every Computer Scientist Should Know About Floating-Point Arithmetic (oracle.com)

2 Likes

Thanks for all the insightful comments, examples, and links. I've read through several of them and have learned a lot!

1 Like

This topic was automatically closed 90 days after the last reply. We invite you to open a new topic if you have further questions or comments.