1// Copyright 2014-2020 Optimal Computing (NZ) Ltd.
2// Licensed under the MIT license. See LICENSE for details.
3
4use super::Ulps;
5#[cfg(feature = "num-traits")]
6#[allow(unused_imports)]
7use num_traits::float::FloatCore;
8
9/// ApproxEqUlps is a trait for approximate equality comparisons.
10/// The associated type Flt is a floating point type which implements Ulps, and is
11/// required so that this trait can be implemented for compound types (e.g. vectors),
12/// not just for the floats themselves.
13pub trait ApproxEqUlps {
14 type Flt: Ulps;
15
16 /// This method tests for `self` and `other` values to be approximately equal
17 /// within ULPs (Units of Least Precision) floating point representations.
18 /// Differing signs are always unequal with this method, and zeroes are only
19 /// equal to zeroes. Use approx_eq() from the ApproxEq trait if that is more
20 /// appropriate.
21 fn approx_eq_ulps(&self, other: &Self, ulps: <Self::Flt as Ulps>::U) -> bool;
22
23 /// This method tests for `self` and `other` values to be not approximately
24 /// equal within ULPs (Units of Least Precision) floating point representations.
25 /// Differing signs are always unequal with this method, and zeroes are only
26 /// equal to zeroes. Use approx_eq() from the ApproxEq trait if that is more
27 /// appropriate.
28 #[inline]
29 fn approx_ne_ulps(&self, other: &Self, ulps: <Self::Flt as Ulps>::U) -> bool {
30 !self.approx_eq_ulps(other, ulps)
31 }
32}
33
34impl ApproxEqUlps for f32 {
35 type Flt = f32;
36
37 fn approx_eq_ulps(&self, other: &f32, ulps: i32) -> bool {
38 // -0 and +0 are drastically far in ulps terms, so
39 // we need a special case for that.
40 if *self == *other {
41 return true;
42 }
43
44 // Handle differing signs as a special case, even if
45 // they are very close, most people consider them
46 // unequal.
47 if self.is_sign_positive() != other.is_sign_positive() {
48 return false;
49 }
50
51 let diff: i32 = self.ulps(other);
52 diff >= -ulps && diff <= ulps
53 }
54}
55
56#[test]
57fn f32_approx_eq_ulps_test1() {
58 let f: f32 = 0.1_f32;
59 let mut sum: f32 = 0.0_f32;
60 for _ in 0_isize..10_isize {
61 sum += f;
62 }
63 let product: f32 = f * 10.0_f32;
64 assert!(sum != product); // Should not be directly equal:
65 assert!(sum.approx_eq_ulps(&product, 1) == true); // But should be close
66 assert!(sum.approx_eq_ulps(&product, 0) == false);
67}
68#[test]
69fn f32_approx_eq_ulps_test2() {
70 let x: f32 = 1000000_f32;
71 let y: f32 = 1000000.1_f32;
72 assert!(x != y); // Should not be directly equal
73 assert!(x.approx_eq_ulps(&y, 2) == true);
74 assert!(x.approx_eq_ulps(&y, 1) == false);
75}
76#[test]
77fn f32_approx_eq_ulps_test_zeroes() {
78 let x: f32 = 0.0_f32;
79 let y: f32 = -0.0_f32;
80 assert!(x.approx_eq_ulps(&y, 0) == true);
81}
82
83impl ApproxEqUlps for f64 {
84 type Flt = f64;
85
86 fn approx_eq_ulps(&self, other: &f64, ulps: i64) -> bool {
87 // -0 and +0 are drastically far in ulps terms, so
88 // we need a special case for that.
89 if *self == *other {
90 return true;
91 }
92
93 // Handle differing signs as a special case, even if
94 // they are very close, most people consider them
95 // unequal.
96 if self.is_sign_positive() != other.is_sign_positive() {
97 return false;
98 }
99
100 let diff: i64 = self.ulps(other);
101 diff >= -ulps && diff <= ulps
102 }
103}
104
105#[test]
106fn f64_approx_eq_ulps_test1() {
107 let f: f64 = 0.1_f64;
108 let mut sum: f64 = 0.0_f64;
109 for _ in 0_isize..10_isize {
110 sum += f;
111 }
112 let product: f64 = f * 10.0_f64;
113 assert!(sum != product); // Should not be directly equal:
114 assert!(sum.approx_eq_ulps(&product, 1) == true); // But should be close
115 assert!(sum.approx_eq_ulps(&product, 0) == false);
116}
117#[test]
118fn f64_approx_eq_ulps_test2() {
119 let x: f64 = 1000000_f64;
120 let y: f64 = 1000000.0000000003_f64;
121 assert!(x != y); // Should not be directly equal
122 assert!(x.approx_eq_ulps(&y, 3) == true);
123 assert!(x.approx_eq_ulps(&y, 2) == false);
124}
125#[test]
126fn f64_approx_eq_ulps_test_zeroes() {
127 let x: f64 = 0.0_f64;
128 let y: f64 = -0.0_f64;
129 assert!(x.approx_eq_ulps(&y, 0) == true);
130}
131