1 // Copyright 2014-2020 Optimal Computing (NZ) Ltd.
2 // Licensed under the MIT license. See LICENSE for details.
3
4 #[cfg(feature="num_traits")]
5 use 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")]
10 pub 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.
ulps(&self, other: &Self) -> <Self as Ulps>::U16 fn ulps(&self, other: &Self) -> <Self as Ulps>::U;
17
18 /// The next representable number above this one
next(&self) -> Self19 fn next(&self) -> Self;
20
21 /// The previous representable number below this one
prev(&self) -> Self22 fn prev(&self) -> Self;
23 }
24
25 #[cfg(not(feature="num_traits"))]
26 pub 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.
ulps(&self, other: &Self) -> <Self as Ulps>::U32 fn ulps(&self, other: &Self) -> <Self as Ulps>::U;
33
34 /// The next representable number above this one
next(&self) -> Self35 fn next(&self) -> Self;
36
37 /// The previous representable number below this one
prev(&self) -> Self38 fn prev(&self) -> Self;
39 }
40
41 impl Ulps for f32 {
42 type U = i32;
43
ulps(&self, other: &f32) -> i3244 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
next(&self) -> Self59 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
prev(&self) -> Self75 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]
f32_ulps_test1()93 fn 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]
f32_ulps_test2()100 fn 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]
f32_ulps_test3()106 fn 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]
f32_ulps_test4()113 fn 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]
f32_ulps_test5()120 fn 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]
f32_ulps_test6()128 fn 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
137 impl Ulps for f64 {
138 type U = i64;
139
ulps(&self, other: &f64) -> i64140 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
next(&self) -> Self155 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
prev(&self) -> Self171 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]
f64_ulps_test1()189 fn 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]
f64_ulps_test2()196 fn 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]
f64_ulps_test3()202 fn 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]
f64_ulps_test4()209 fn 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]
f64_ulps_test5()216 fn 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]
f64_ulps_test6()224 fn 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
233