1// Copyright 2014-2020 Optimal Computing (NZ) Ltd.
2// Licensed under the MIT license. See LICENSE for details.
3
4use core::{f32, f64};
5#[cfg(feature = "num-traits")]
6#[allow(unused_imports)]
7use num_traits::float::FloatCore;
8use super::Ulps;
9
10/// A trait for approximate equality comparisons.
11pub trait ApproxEq: Sized {
12 /// This type type defines a margin within which two values are to be
13 /// considered approximately equal. It must implement `Default` so that
14 /// `approx_eq()` can be called on unknown types.
15 type Margin: Copy + Default;
16
17 /// This method tests that the `self` and `other` values are equal within `margin`
18 /// of each other.
19 fn approx_eq<M: Into<Self::Margin>>(self, other: Self, margin: M) -> bool;
20
21 /// This method tests that the `self` and `other` values are not within `margin`
22 /// of each other.
23 fn approx_ne<M: Into<Self::Margin>>(self, other: Self, margin: M) -> bool {
24 !self.approx_eq(other, margin)
25 }
26}
27
28/// This type defines a margin within two `f32` values might be considered equal,
29/// and is intended as the associated type for the `ApproxEq` trait.
30///
31/// Two tests are used to determine approximate equality.
32///
33/// The first test considers two values approximately equal if they differ by <=
34/// `epsilon`. This will only succeed for very small numbers. Note that it may
35/// succeed even if the parameters are of differing signs, straddling zero.
36///
37/// The second test considers how many ULPs (units of least precision, units in
38/// the last place, which is the integer number of floating-point representations
39/// that the parameters are separated by) different the parameters are and considers
40/// them approximately equal if this is <= `ulps`. For large floating-point numbers,
41/// an ULP can be a rather large gap, but this kind of comparison is necessary
42/// because floating-point operations must round to the nearest representable value
43/// and so larger floating-point values accumulate larger errors.
44#[repr(C)]
45#[derive(Debug, Clone, Copy)]
46pub struct F32Margin {
47 pub epsilon: f32,
48 pub ulps: i32
49}
50impl Default for F32Margin {
51 #[inline]
52 fn default() -> F32Margin {
53 F32Margin {
54 epsilon: f32::EPSILON,
55 ulps: 4
56 }
57 }
58}
59impl F32Margin {
60 #[inline]
61 pub fn zero() -> F32Margin {
62 F32Margin {
63 epsilon: 0.0,
64 ulps: 0
65 }
66 }
67 pub fn epsilon(self, epsilon: f32) -> Self {
68 F32Margin {
69 epsilon,
70 ..self
71 }
72 }
73 pub fn ulps(self, ulps: i32) -> Self {
74 F32Margin {
75 ulps,
76 ..self
77 }
78 }
79}
80impl From<(f32, i32)> for F32Margin {
81 fn from(m: (f32, i32)) -> F32Margin {
82 F32Margin {
83 epsilon: m.0,
84 ulps: m.1
85 }
86 }
87}
88
89impl ApproxEq for f32 {
90 type Margin = F32Margin;
91
92 fn approx_eq<M: Into<Self::Margin>>(self, other: f32, margin: M) -> bool {
93 let margin = margin.into();
94
95 // Check for exact equality first. This is often true, and so we get the
96 // performance benefit of only doing one compare in most cases.
97 self==other ||
98
99 // Perform epsilon comparison next
100 ((self - other).abs() <= margin.epsilon) ||
101
102 {
103 // Perform ulps comparion last
104 let diff: i32 = self.ulps(&other);
105 saturating_abs_i32!(diff) <= margin.ulps
106 }
107 }
108}
109
110#[test]
111fn f32_approx_eq_test1() {
112 let f: f32 = 0.0_f32;
113 let g: f32 = -0.0000000000000005551115123125783_f32;
114 assert!(f != g); // Should not be directly equal
115 assert!(f.approx_eq(g, (f32::EPSILON, 0)) == true);
116}
117#[test]
118fn f32_approx_eq_test2() {
119 let f: f32 = 0.0_f32;
120 let g: f32 = -0.0_f32;
121 assert!(f.approx_eq(g, (f32::EPSILON, 0)) == true);
122}
123#[test]
124fn f32_approx_eq_test3() {
125 let f: f32 = 0.0_f32;
126 let g: f32 = 0.00000000000000001_f32;
127 assert!(f.approx_eq(g, (f32::EPSILON, 0)) == true);
128}
129#[test]
130fn f32_approx_eq_test4() {
131 let f: f32 = 0.00001_f32;
132 let g: f32 = 0.00000000000000001_f32;
133 assert!(f.approx_eq(g, (f32::EPSILON, 0)) == false);
134}
135#[test]
136fn f32_approx_eq_test5() {
137 let f: f32 = 0.1_f32;
138 let mut sum: f32 = 0.0_f32;
139 for _ in 0_isize..10_isize { sum += f; }
140 let product: f32 = f * 10.0_f32;
141 assert!(sum != product); // Should not be directly equal:
142 assert!(sum.approx_eq(product, (f32::EPSILON, 1)) == true);
143 assert!(sum.approx_eq(product, F32Margin::zero()) == false);
144}
145#[test]
146fn f32_approx_eq_test6() {
147 let x: f32 = 1000000_f32;
148 let y: f32 = 1000000.1_f32;
149 assert!(x != y); // Should not be directly equal
150 assert!(x.approx_eq(y, (0.0, 2)) == true); // 2 ulps does it
151 // epsilon method no good here:
152 assert!(x.approx_eq(y, (1000.0 * f32::EPSILON, 0)) == false);
153}
154
155/// This type defines a margin within two `f64` values might be considered equal,
156/// and is intended as the associated type for the `ApproxEq` trait.
157///
158/// Two tests are used to determine approximate equality.
159///
160/// The first test considers two values approximately equal if they differ by <=
161/// `epsilon`. This will only succeed for very small numbers. Note that it may
162/// succeed even if the parameters are of differing signs, straddling zero.
163///
164/// The second test considers how many ULPs (units of least precision, units in
165/// the last place, which is the integer number of floating-point representations
166/// that the parameters are separated by) different the parameters are and considers
167/// them approximately equal if this is <= `ulps`. For large floating-point numbers,
168/// an ULP can be a rather large gap, but this kind of comparison is necessary
169/// because floating-point operations must round to the nearest representable value
170/// and so larger floating-point values accumulate larger errors.
171#[derive(Debug, Clone, Copy)]
172pub struct F64Margin {
173 pub epsilon: f64,
174 pub ulps: i64
175}
176impl Default for F64Margin {
177 #[inline]
178 fn default() -> F64Margin {
179 F64Margin {
180 epsilon: f64::EPSILON,
181 ulps: 4
182 }
183 }
184}
185impl F64Margin {
186 #[inline]
187 pub fn zero() -> F64Margin {
188 F64Margin {
189 epsilon: 0.0,
190 ulps: 0
191 }
192 }
193 pub fn epsilon(self, epsilon: f64) -> Self {
194 F64Margin {
195 epsilon,
196 ..self
197 }
198 }
199 pub fn ulps(self, ulps: i64) -> Self {
200 F64Margin {
201 ulps,
202 ..self
203 }
204 }
205}
206impl From<(f64, i64)> for F64Margin {
207 fn from(m: (f64, i64)) -> F64Margin {
208 F64Margin {
209 epsilon: m.0,
210 ulps: m.1
211 }
212 }
213}
214
215impl ApproxEq for f64 {
216 type Margin = F64Margin;
217
218 fn approx_eq<M: Into<Self::Margin>>(self, other: f64, margin: M) -> bool {
219 let margin = margin.into();
220
221 // Check for exact equality first. This is often true, and so we get the
222 // performance benefit of only doing one compare in most cases.
223 self == other ||
224
225 // Perform epsilon comparison next
226 ((self - other).abs() <= margin.epsilon) ||
227
228 {
229 // Perform ulps comparion last
230 let diff: i64 = self.ulps(&other);
231 saturating_abs_i64!(diff) <= margin.ulps
232 }
233 }
234}
235
236#[test]
237fn f64_approx_eq_test1() {
238 let f: f64 = 0.0_f64;
239 let g: f64 = -0.0000000000000005551115123125783_f64;
240 assert!(f != g); // Should not be precisely equal.
241 assert!(f.approx_eq(g, (3.0 * f64::EPSILON, 0)) == true); // 3e is enough.
242 // ULPs test won't ever call these equal.
243}
244#[test]
245fn f64_approx_eq_test2() {
246 let f: f64 = 0.0_f64;
247 let g: f64 = -0.0_f64;
248 assert!(f.approx_eq(g, (f64::EPSILON, 0)) == true);
249}
250#[test]
251fn f64_approx_eq_test3() {
252 let f: f64 = 0.0_f64;
253 let g: f64 = 1e-17_f64;
254 assert!(f.approx_eq(g, (f64::EPSILON, 0)) == true);
255}
256#[test]
257fn f64_approx_eq_test4() {
258 let f: f64 = 0.00001_f64;
259 let g: f64 = 0.00000000000000001_f64;
260 assert!(f.approx_eq(g, (f64::EPSILON, 0)) == false);
261}
262#[test]
263fn f64_approx_eq_test5() {
264 let f: f64 = 0.1_f64;
265 let mut sum: f64 = 0.0_f64;
266 for _ in 0_isize..10_isize { sum += f; }
267 let product: f64 = f * 10.0_f64;
268 assert!(sum != product); // Should not be precisely equaly.
269 assert!(sum.approx_eq(product, (f64::EPSILON, 0)) == true);
270 assert!(sum.approx_eq(product, (0.0, 1)) == true);
271}
272#[test]
273fn f64_approx_eq_test6() {
274 let x: f64 = 1000000_f64;
275 let y: f64 = 1000000.0000000003_f64;
276 assert!(x != y); // Should not be precisely equal.
277 assert!(x.approx_eq(y, (0.0, 3)) == true);
278}
279#[test]
280fn f64_code_triggering_issue_20() {
281 assert_eq!((-25.0f64).approx_eq(25.0, (0.00390625, 1)), false);
282}
283
284impl<T> ApproxEq for &[T]
285where T: Copy + ApproxEq {
286 type Margin = <T as ApproxEq>::Margin;
287
288 fn approx_eq<M: Into<Self::Margin>>(self, other: Self, margin: M) -> bool {
289 let margin = margin.into();
290 if self.len() != other.len() { return false; }
291 self.iter().zip(other.iter()).all(|(a: &T,b: &T)| {
292 a.approx_eq(*b, margin)
293 })
294 }
295}
296
297#[test]
298fn test_slices() {
299 assert!( [1.33, 2.4, 2.5].approx_eq(&[1.33, 2.4, 2.5], (0.0, 0_i64)) );
300 assert!( ! [1.33, 2.4, 2.6].approx_eq(&[1.33, 2.4, 2.5], (0.0, 0_i64)) );
301 assert!( ! [1.33, 2.4].approx_eq(&[1.33, 2.4, 2.5], (0.0, 0_i64)) );
302 assert!( ! [1.33, 2.4, 2.5].approx_eq(&[1.33, 2.4], (0.0, 0_i64)) );
303}
304
305