float_cmp/
ulps.rs

1// Copyright 2014-2020 Optimal Computing (NZ) Ltd.
2// Licensed under the MIT license.  See LICENSE for details.
3
4#[cfg(feature="num_traits")]
5use num_traits::NumCast;
6
7/// A trait for floating point numbers which computes the number of representable
8/// values or ULPs (Units of Least Precision) that separate the two given values.
9#[cfg(feature="num_traits")]
10pub trait Ulps {
11    type U: Copy + NumCast;
12
13    /// The number of representable values or ULPs (Units of Least Precision) that
14    /// separate `self` and `other`.  The result `U` is an integral value, and will
15    /// be zero if `self` and `other` are exactly equal.
16    fn ulps(&self, other: &Self) -> <Self as Ulps>::U;
17
18    /// The next representable number above this one
19    fn next(&self) -> Self;
20
21    /// The previous representable number below this one
22    fn prev(&self) -> Self;
23}
24
25#[cfg(not(feature="num_traits"))]
26pub trait Ulps {
27    type U: Copy;
28
29    /// The number of representable values or ULPs (Units of Least Precision) that
30    /// separate `self` and `other`.  The result `U` is an integral value, and will
31    /// be zero if `self` and `other` are exactly equal.
32    fn ulps(&self, other: &Self) -> <Self as Ulps>::U;
33
34    /// The next representable number above this one
35    fn next(&self) -> Self;
36
37    /// The previous representable number below this one
38    fn prev(&self) -> Self;
39}
40
41impl Ulps for f32 {
42    type U = i32;
43
44    fn ulps(&self, other: &f32) -> i32 {
45
46        // IEEE754 defined floating point storage representation to
47        // maintain their order when their bit patterns are interpreted as
48        // integers.  This is a huge boon to the task at hand, as we can
49        // reinterpret them as integers to find out how many ULPs apart any
50        // two floats are
51
52        // Setup integer representations of the input
53        let ai32: i32 = self.to_bits() as i32;
54        let bi32: i32 = other.to_bits() as i32;
55
56        ai32.wrapping_sub(bi32)
57    }
58
59    fn next(&self) -> Self {
60        if self.is_infinite() && *self > 0.0 {
61            *self
62        } else if *self == -0.0 && self.is_sign_negative() {
63            0.0
64        } else {
65            let mut u = self.to_bits();
66            if *self >= 0.0 {
67                u += 1;
68            } else {
69                u -= 1;
70            }
71            f32::from_bits(u)
72        }
73    }
74
75    fn prev(&self) -> Self {
76        if self.is_infinite() && *self < 0.0 {
77            *self
78        } else if *self == 0.0 && self.is_sign_positive() {
79            -0.0
80        } else {
81            let mut u = self.to_bits();
82            if *self <= -0.0 {
83                u += 1;
84            } else {
85                u -= 1;
86            }
87            f32::from_bits(u)
88        }
89    }
90}
91
92#[test]
93fn f32_ulps_test1() {
94    let x: f32 = 1000000_f32;
95    let y: f32 = 1000000.1_f32;
96    assert!(x.ulps(&y) == -2);
97}
98
99#[test]
100fn f32_ulps_test2() {
101    let pzero: f32 = f32::from_bits(0x00000000_u32);
102    let nzero: f32 = f32::from_bits(0x80000000_u32);
103    assert!(pzero.ulps(&nzero) == -2147483648);
104}
105#[test]
106fn f32_ulps_test3() {
107    let pinf: f32 = f32::from_bits(0x7f800000_u32);
108    let ninf: f32 = f32::from_bits(0xff800000_u32);
109    assert!(pinf.ulps(&ninf) == -2147483648);
110}
111
112#[test]
113fn f32_ulps_test4() {
114    let x: f32 = f32::from_bits(0x63a7f026_u32);
115    let y: f32 = f32::from_bits(0x63a7f023_u32);
116    assert!(x.ulps(&y) == 3);
117}
118
119#[test]
120fn f32_ulps_test5() {
121    let x: f32 = 2.0;
122    let ulps: i32 = x.to_bits() as i32;
123    let x2: f32 = <f32>::from_bits(ulps as u32);
124    assert_eq!(x, x2);
125}
126
127#[test]
128fn f32_ulps_test6() {
129    let negzero: f32 = -0.;
130    let zero: f32 = 0.;
131    assert_eq!(negzero.next(), zero);
132    assert_eq!(zero.prev(), negzero);
133    assert!(negzero.prev() < 0.0);
134    assert!(zero.next() > 0.0);
135}
136
137impl Ulps for f64 {
138    type U = i64;
139
140    fn ulps(&self, other: &f64) -> i64 {
141
142        // IEEE754 defined floating point storage representation to
143        // maintain their order when their bit patterns are interpreted as
144        // integers.  This is a huge boon to the task at hand, as we can
145        // reinterpret them as integers to find out how many ULPs apart any
146        // two floats are
147
148        // Setup integer representations of the input
149        let ai64: i64 = self.to_bits() as i64;
150        let bi64: i64 = other.to_bits() as i64;
151
152        ai64.wrapping_sub(bi64)
153    }
154
155    fn next(&self) -> Self {
156        if self.is_infinite() && *self > 0.0 {
157            *self
158        } else if *self == -0.0 && self.is_sign_negative() {
159            0.0
160        } else {
161            let mut u = self.to_bits();
162            if *self >= 0.0 {
163                u += 1;
164            } else {
165                u -= 1;
166            }
167            f64::from_bits(u)
168        }
169    }
170
171    fn prev(&self) -> Self {
172        if self.is_infinite() && *self < 0.0 {
173            *self
174        } else if *self == 0.0 && self.is_sign_positive() {
175            -0.0
176        } else {
177            let mut u = self.to_bits();
178            if *self <= -0.0 {
179                u += 1;
180            } else {
181                u -= 1;
182            }
183            f64::from_bits(u)
184        }
185    }
186}
187
188#[test]
189fn f64_ulps_test1() {
190    let x: f64 = 1000000_f64;
191    let y: f64 = 1000000.00000001_f64;
192    assert!(x.ulps(&y) == -86);
193}
194
195#[test]
196fn f64_ulps_test2() {
197    let pzero: f64 = f64::from_bits(0x0000000000000000_u64);
198    let nzero: f64 = f64::from_bits(0x8000000000000000_u64);
199    assert!(pzero.ulps(&nzero) == -9223372036854775808i64);
200}
201#[test]
202fn f64_ulps_test3() {
203    let pinf: f64 = f64::from_bits(0x7f80000000000000_u64);
204    let ninf: f64 = f64::from_bits(0xff80000000000000_u64);
205    assert!(pinf.ulps(&ninf) == -9223372036854775808i64);
206}
207
208#[test]
209fn f64_ulps_test4() {
210    let x: f64 = f64::from_bits(0xd017f6cc63a7f026_u64);
211    let y: f64 = f64::from_bits(0xd017f6cc63a7f023_u64);
212    assert!(x.ulps(&y) == 3);
213}
214
215#[test]
216fn f64_ulps_test5() {
217    let x: f64 = 2.0;
218    let ulps: i64 = x.to_bits() as i64;
219    let x2: f64 = <f64>::from_bits(ulps as u64);
220    assert_eq!(x, x2);
221}
222
223#[test]
224fn f64_ulps_test6() {
225    let negzero: f64 = -0.;
226    let zero: f64 = 0.;
227    assert_eq!(negzero.next(), zero);
228    assert_eq!(zero.prev(), negzero);
229    assert!(negzero.prev() < 0.0);
230    assert!(zero.next() > 0.0);
231}
232