1 | //! Utilities for Rust numbers. |
2 | |
3 | #![doc (hidden)] |
4 | |
5 | #[cfg (all(not(feature = "std" ), feature = "compact" ))] |
6 | use crate::libm::{powd, powf}; |
7 | #[cfg (not(feature = "compact" ))] |
8 | use crate::table::{SMALL_F32_POW10, SMALL_F64_POW10, SMALL_INT_POW10, SMALL_INT_POW5}; |
9 | #[cfg (not(feature = "compact" ))] |
10 | use core::hint; |
11 | use core::ops; |
12 | |
13 | /// Generic floating-point type, to be used in generic code for parsing. |
14 | /// |
15 | /// Although the trait is part of the public API, the trait provides methods |
16 | /// and constants that are effectively non-public: they may be removed |
17 | /// at any time without any breaking changes. |
18 | pub trait Float: |
19 | Sized |
20 | + Copy |
21 | + PartialEq |
22 | + PartialOrd |
23 | + Send |
24 | + Sync |
25 | + ops::Add<Output = Self> |
26 | + ops::AddAssign |
27 | + ops::Div<Output = Self> |
28 | + ops::DivAssign |
29 | + ops::Mul<Output = Self> |
30 | + ops::MulAssign |
31 | + ops::Rem<Output = Self> |
32 | + ops::RemAssign |
33 | + ops::Sub<Output = Self> |
34 | + ops::SubAssign |
35 | + ops::Neg<Output = Self> |
36 | { |
37 | /// Maximum number of digits that can contribute in the mantissa. |
38 | /// |
39 | /// We can exactly represent a float in radix `b` from radix 2 if |
40 | /// `b` is divisible by 2. This function calculates the exact number of |
41 | /// digits required to exactly represent that float. |
42 | /// |
43 | /// According to the "Handbook of Floating Point Arithmetic", |
44 | /// for IEEE754, with emin being the min exponent, p2 being the |
45 | /// precision, and b being the radix, the number of digits follows as: |
46 | /// |
47 | /// `−emin + p2 + ⌊(emin + 1) log(2, b) − log(1 − 2^(−p2), b)⌋` |
48 | /// |
49 | /// For f32, this follows as: |
50 | /// emin = -126 |
51 | /// p2 = 24 |
52 | /// |
53 | /// For f64, this follows as: |
54 | /// emin = -1022 |
55 | /// p2 = 53 |
56 | /// |
57 | /// In Python: |
58 | /// `-emin + p2 + math.floor((emin+1)*math.log(2, b) - math.log(1-2**(-p2), b))` |
59 | /// |
60 | /// This was used to calculate the maximum number of digits for [2, 36]. |
61 | const MAX_DIGITS: usize; |
62 | |
63 | // MASKS |
64 | |
65 | /// Bitmask for the sign bit. |
66 | const SIGN_MASK: u64; |
67 | /// Bitmask for the exponent, including the hidden bit. |
68 | const EXPONENT_MASK: u64; |
69 | /// Bitmask for the hidden bit in exponent, which is an implicit 1 in the fraction. |
70 | const HIDDEN_BIT_MASK: u64; |
71 | /// Bitmask for the mantissa (fraction), excluding the hidden bit. |
72 | const MANTISSA_MASK: u64; |
73 | |
74 | // PROPERTIES |
75 | |
76 | /// Size of the significand (mantissa) without hidden bit. |
77 | const MANTISSA_SIZE: i32; |
78 | /// Bias of the exponet |
79 | const EXPONENT_BIAS: i32; |
80 | /// Exponent portion of a denormal float. |
81 | const DENORMAL_EXPONENT: i32; |
82 | /// Maximum exponent value in float. |
83 | const MAX_EXPONENT: i32; |
84 | |
85 | // ROUNDING |
86 | |
87 | /// Mask to determine if a full-carry occurred (1 in bit above hidden bit). |
88 | const CARRY_MASK: u64; |
89 | |
90 | /// Bias for marking an invalid extended float. |
91 | // Value is `i16::MIN`, using hard-coded constants for older Rustc versions. |
92 | const INVALID_FP: i32 = -0x8000; |
93 | |
94 | // Maximum mantissa for the fast-path (`1 << 53` for f64). |
95 | const MAX_MANTISSA_FAST_PATH: u64 = 2_u64 << Self::MANTISSA_SIZE; |
96 | |
97 | // Largest exponent value `(1 << EXP_BITS) - 1`. |
98 | const INFINITE_POWER: i32 = Self::MAX_EXPONENT + Self::EXPONENT_BIAS; |
99 | |
100 | // Round-to-even only happens for negative values of q |
101 | // when q ≥ −4 in the 64-bit case and when q ≥ −17 in |
102 | // the 32-bitcase. |
103 | // |
104 | // When q ≥ 0,we have that 5^q ≤ 2m+1. In the 64-bit case,we |
105 | // have 5^q ≤ 2m+1 ≤ 2^54 or q ≤ 23. In the 32-bit case,we have |
106 | // 5^q ≤ 2m+1 ≤ 2^25 or q ≤ 10. |
107 | // |
108 | // When q < 0, we have w ≥ (2m+1)×5^−q. We must have that w < 2^64 |
109 | // so (2m+1)×5^−q < 2^64. We have that 2m+1 > 2^53 (64-bit case) |
110 | // or 2m+1 > 2^24 (32-bit case). Hence,we must have 2^53×5^−q < 2^64 |
111 | // (64-bit) and 2^24×5^−q < 2^64 (32-bit). Hence we have 5^−q < 2^11 |
112 | // or q ≥ −4 (64-bit case) and 5^−q < 2^40 or q ≥ −17 (32-bitcase). |
113 | // |
114 | // Thus we have that we only need to round ties to even when |
115 | // we have that q ∈ [−4,23](in the 64-bit case) or q∈[−17,10] |
116 | // (in the 32-bit case). In both cases,the power of five(5^|q|) |
117 | // fits in a 64-bit word. |
118 | const MIN_EXPONENT_ROUND_TO_EVEN: i32; |
119 | const MAX_EXPONENT_ROUND_TO_EVEN: i32; |
120 | |
121 | /// Minimum normal exponent value `-(1 << (EXPONENT_SIZE - 1)) + 1`. |
122 | const MINIMUM_EXPONENT: i32; |
123 | |
124 | /// Smallest decimal exponent for a non-zero value. |
125 | const SMALLEST_POWER_OF_TEN: i32; |
126 | |
127 | /// Largest decimal exponent for a non-infinite value. |
128 | const LARGEST_POWER_OF_TEN: i32; |
129 | |
130 | /// Minimum exponent that for a fast path case, or `-⌊(MANTISSA_SIZE+1)/log2(10)⌋` |
131 | const MIN_EXPONENT_FAST_PATH: i32; |
132 | |
133 | /// Maximum exponent that for a fast path case, or `⌊(MANTISSA_SIZE+1)/log2(5)⌋` |
134 | const MAX_EXPONENT_FAST_PATH: i32; |
135 | |
136 | /// Maximum exponent that can be represented for a disguised-fast path case. |
137 | /// This is `MAX_EXPONENT_FAST_PATH + ⌊(MANTISSA_SIZE+1)/log2(10)⌋` |
138 | const MAX_EXPONENT_DISGUISED_FAST_PATH: i32; |
139 | |
140 | /// Convert 64-bit integer to float. |
141 | fn from_u64(u: u64) -> Self; |
142 | |
143 | // Re-exported methods from std. |
144 | fn from_bits(u: u64) -> Self; |
145 | fn to_bits(self) -> u64; |
146 | |
147 | /// Get a small power-of-radix for fast-path multiplication. |
148 | /// |
149 | /// # Safety |
150 | /// |
151 | /// Safe as long as the exponent is smaller than the table size. |
152 | unsafe fn pow_fast_path(exponent: usize) -> Self; |
153 | |
154 | /// Get a small, integral power-of-radix for fast-path multiplication. |
155 | /// |
156 | /// # Safety |
157 | /// |
158 | /// Safe as long as the exponent is smaller than the table size. |
159 | #[inline (always)] |
160 | unsafe fn int_pow_fast_path(exponent: usize, radix: u32) -> u64 { |
161 | // SAFETY: safe as long as the exponent is smaller than the radix table. |
162 | #[cfg (not(feature = "compact" ))] |
163 | return match radix { |
164 | 5 => unsafe { *SMALL_INT_POW5.get_unchecked(exponent) }, |
165 | 10 => unsafe { *SMALL_INT_POW10.get_unchecked(exponent) }, |
166 | _ => unsafe { hint::unreachable_unchecked() }, |
167 | }; |
168 | |
169 | #[cfg (feature = "compact" )] |
170 | return (radix as u64).pow(exponent as u32); |
171 | } |
172 | |
173 | /// Returns true if the float is a denormal. |
174 | #[inline ] |
175 | fn is_denormal(self) -> bool { |
176 | self.to_bits() & Self::EXPONENT_MASK == 0 |
177 | } |
178 | |
179 | /// Get exponent component from the float. |
180 | #[inline ] |
181 | fn exponent(self) -> i32 { |
182 | if self.is_denormal() { |
183 | return Self::DENORMAL_EXPONENT; |
184 | } |
185 | |
186 | let bits = self.to_bits(); |
187 | let biased_e: i32 = ((bits & Self::EXPONENT_MASK) >> Self::MANTISSA_SIZE) as i32; |
188 | biased_e - Self::EXPONENT_BIAS |
189 | } |
190 | |
191 | /// Get mantissa (significand) component from float. |
192 | #[inline ] |
193 | fn mantissa(self) -> u64 { |
194 | let bits = self.to_bits(); |
195 | let s = bits & Self::MANTISSA_MASK; |
196 | if !self.is_denormal() { |
197 | s + Self::HIDDEN_BIT_MASK |
198 | } else { |
199 | s |
200 | } |
201 | } |
202 | } |
203 | |
204 | impl Float for f32 { |
205 | const MAX_DIGITS: usize = 114; |
206 | const SIGN_MASK: u64 = 0x80000000; |
207 | const EXPONENT_MASK: u64 = 0x7F800000; |
208 | const HIDDEN_BIT_MASK: u64 = 0x00800000; |
209 | const MANTISSA_MASK: u64 = 0x007FFFFF; |
210 | const MANTISSA_SIZE: i32 = 23; |
211 | const EXPONENT_BIAS: i32 = 127 + Self::MANTISSA_SIZE; |
212 | const DENORMAL_EXPONENT: i32 = 1 - Self::EXPONENT_BIAS; |
213 | const MAX_EXPONENT: i32 = 0xFF - Self::EXPONENT_BIAS; |
214 | const CARRY_MASK: u64 = 0x1000000; |
215 | const MIN_EXPONENT_ROUND_TO_EVEN: i32 = -17; |
216 | const MAX_EXPONENT_ROUND_TO_EVEN: i32 = 10; |
217 | const MINIMUM_EXPONENT: i32 = -127; |
218 | const SMALLEST_POWER_OF_TEN: i32 = -65; |
219 | const LARGEST_POWER_OF_TEN: i32 = 38; |
220 | const MIN_EXPONENT_FAST_PATH: i32 = -10; |
221 | const MAX_EXPONENT_FAST_PATH: i32 = 10; |
222 | const MAX_EXPONENT_DISGUISED_FAST_PATH: i32 = 17; |
223 | |
224 | #[inline (always)] |
225 | unsafe fn pow_fast_path(exponent: usize) -> Self { |
226 | // SAFETY: safe as long as the exponent is smaller than the radix table. |
227 | #[cfg (not(feature = "compact" ))] |
228 | return unsafe { *SMALL_F32_POW10.get_unchecked(exponent) }; |
229 | |
230 | #[cfg (feature = "compact" )] |
231 | return powf(10.0f32, exponent as f32); |
232 | } |
233 | |
234 | #[inline ] |
235 | fn from_u64(u: u64) -> f32 { |
236 | u as _ |
237 | } |
238 | |
239 | #[inline ] |
240 | fn from_bits(u: u64) -> f32 { |
241 | // Constant is `u32::MAX` for older Rustc versions. |
242 | debug_assert!(u <= 0xffff_ffff); |
243 | f32::from_bits(u as u32) |
244 | } |
245 | |
246 | #[inline ] |
247 | fn to_bits(self) -> u64 { |
248 | f32::to_bits(self) as u64 |
249 | } |
250 | } |
251 | |
252 | impl Float for f64 { |
253 | const MAX_DIGITS: usize = 769; |
254 | const SIGN_MASK: u64 = 0x8000000000000000; |
255 | const EXPONENT_MASK: u64 = 0x7FF0000000000000; |
256 | const HIDDEN_BIT_MASK: u64 = 0x0010000000000000; |
257 | const MANTISSA_MASK: u64 = 0x000FFFFFFFFFFFFF; |
258 | const MANTISSA_SIZE: i32 = 52; |
259 | const EXPONENT_BIAS: i32 = 1023 + Self::MANTISSA_SIZE; |
260 | const DENORMAL_EXPONENT: i32 = 1 - Self::EXPONENT_BIAS; |
261 | const MAX_EXPONENT: i32 = 0x7FF - Self::EXPONENT_BIAS; |
262 | const CARRY_MASK: u64 = 0x20000000000000; |
263 | const MIN_EXPONENT_ROUND_TO_EVEN: i32 = -4; |
264 | const MAX_EXPONENT_ROUND_TO_EVEN: i32 = 23; |
265 | const MINIMUM_EXPONENT: i32 = -1023; |
266 | const SMALLEST_POWER_OF_TEN: i32 = -342; |
267 | const LARGEST_POWER_OF_TEN: i32 = 308; |
268 | const MIN_EXPONENT_FAST_PATH: i32 = -22; |
269 | const MAX_EXPONENT_FAST_PATH: i32 = 22; |
270 | const MAX_EXPONENT_DISGUISED_FAST_PATH: i32 = 37; |
271 | |
272 | #[inline (always)] |
273 | unsafe fn pow_fast_path(exponent: usize) -> Self { |
274 | // SAFETY: safe as long as the exponent is smaller than the radix table. |
275 | #[cfg (not(feature = "compact" ))] |
276 | return unsafe { *SMALL_F64_POW10.get_unchecked(exponent) }; |
277 | |
278 | #[cfg (feature = "compact" )] |
279 | return powd(10.0f64, exponent as f64); |
280 | } |
281 | |
282 | #[inline ] |
283 | fn from_u64(u: u64) -> f64 { |
284 | u as _ |
285 | } |
286 | |
287 | #[inline ] |
288 | fn from_bits(u: u64) -> f64 { |
289 | f64::from_bits(u) |
290 | } |
291 | |
292 | #[inline ] |
293 | fn to_bits(self) -> u64 { |
294 | f64::to_bits(self) |
295 | } |
296 | } |
297 | |
298 | #[inline (always)] |
299 | #[cfg (all(feature = "std" , feature = "compact" ))] |
300 | pub fn powf(x: f32, y: f32) -> f32 { |
301 | x.powf(y) |
302 | } |
303 | |
304 | #[inline (always)] |
305 | #[cfg (all(feature = "std" , feature = "compact" ))] |
306 | pub fn powd(x: f64, y: f64) -> f64 { |
307 | x.powf(y) |
308 | } |
309 | |