1 // Translated from C to Rust. The original C code can be found at
2 // https://github.com/ulfjack/ryu and carries the following license:
3 //
4 // Copyright 2018 Ulf Adams
5 //
6 // The contents of this file may be used under the terms of the Apache License,
7 // Version 2.0.
8 //
9 //    (See accompanying file LICENSE-Apache or copy at
10 //     http://www.apache.org/licenses/LICENSE-2.0)
11 //
12 // Alternatively, the contents of this file may be used under the terms of
13 // the Boost Software License, Version 1.0.
14 //    (See accompanying file LICENSE-Boost or copy at
15 //     https://www.boost.org/LICENSE_1_0.txt)
16 //
17 // Unless required by applicable law or agreed to in writing, this software
18 // is distributed on an "AS IS" BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
19 // KIND, either express or implied.
20 
21 use crate::common::*;
22 use crate::f2s_intrinsics::*;
23 
24 pub const FLOAT_MANTISSA_BITS: u32 = 23;
25 pub const FLOAT_EXPONENT_BITS: u32 = 8;
26 const FLOAT_BIAS: i32 = 127;
27 pub use crate::f2s_intrinsics::{FLOAT_POW5_BITCOUNT, FLOAT_POW5_INV_BITCOUNT};
28 
29 // A floating decimal representing m * 10^e.
30 pub struct FloatingDecimal32 {
31     pub mantissa: u32,
32     // Decimal exponent's range is -45 to 38
33     // inclusive, and can fit in i16 if needed.
34     pub exponent: i32,
35 }
36 
37 #[cfg_attr(feature = "no-panic", inline)]
f2d(ieee_mantissa: u32, ieee_exponent: u32) -> FloatingDecimal3238 pub fn f2d(ieee_mantissa: u32, ieee_exponent: u32) -> FloatingDecimal32 {
39     let (e2, m2) = if ieee_exponent == 0 {
40         (
41             // We subtract 2 so that the bounds computation has 2 additional bits.
42             1 - FLOAT_BIAS - FLOAT_MANTISSA_BITS as i32 - 2,
43             ieee_mantissa,
44         )
45     } else {
46         (
47             ieee_exponent as i32 - FLOAT_BIAS - FLOAT_MANTISSA_BITS as i32 - 2,
48             (1u32 << FLOAT_MANTISSA_BITS) | ieee_mantissa,
49         )
50     };
51     let even = (m2 & 1) == 0;
52     let accept_bounds = even;
53 
54     // Step 2: Determine the interval of valid decimal representations.
55     let mv = 4 * m2;
56     let mp = 4 * m2 + 2;
57     // Implicit bool -> int conversion. True is 1, false is 0.
58     let mm_shift = (ieee_mantissa != 0 || ieee_exponent <= 1) as u32;
59     let mm = 4 * m2 - 1 - mm_shift;
60 
61     // Step 3: Convert to a decimal power base using 64-bit arithmetic.
62     let mut vr: u32;
63     let mut vp: u32;
64     let mut vm: u32;
65     let e10: i32;
66     let mut vm_is_trailing_zeros = false;
67     let mut vr_is_trailing_zeros = false;
68     let mut last_removed_digit = 0u8;
69     if e2 >= 0 {
70         let q = log10_pow2(e2);
71         e10 = q as i32;
72         let k = FLOAT_POW5_INV_BITCOUNT + pow5bits(q as i32) - 1;
73         let i = -e2 + q as i32 + k;
74         vr = mul_pow5_inv_div_pow2(mv, q, i);
75         vp = mul_pow5_inv_div_pow2(mp, q, i);
76         vm = mul_pow5_inv_div_pow2(mm, q, i);
77         if q != 0 && (vp - 1) / 10 <= vm / 10 {
78             // We need to know one removed digit even if we are not going to loop below. We could use
79             // q = X - 1 above, except that would require 33 bits for the result, and we've found that
80             // 32-bit arithmetic is faster even on 64-bit machines.
81             let l = FLOAT_POW5_INV_BITCOUNT + pow5bits(q as i32 - 1) - 1;
82             last_removed_digit =
83                 (mul_pow5_inv_div_pow2(mv, q - 1, -e2 + q as i32 - 1 + l) % 10) as u8;
84         }
85         if q <= 9 {
86             // The largest power of 5 that fits in 24 bits is 5^10, but q <= 9 seems to be safe as well.
87             // Only one of mp, mv, and mm can be a multiple of 5, if any.
88             if mv % 5 == 0 {
89                 vr_is_trailing_zeros = multiple_of_power_of_5_32(mv, q);
90             } else if accept_bounds {
91                 vm_is_trailing_zeros = multiple_of_power_of_5_32(mm, q);
92             } else {
93                 vp -= multiple_of_power_of_5_32(mp, q) as u32;
94             }
95         }
96     } else {
97         let q = log10_pow5(-e2);
98         e10 = q as i32 + e2;
99         let i = -e2 - q as i32;
100         let k = pow5bits(i) - FLOAT_POW5_BITCOUNT;
101         let mut j = q as i32 - k;
102         vr = mul_pow5_div_pow2(mv, i as u32, j);
103         vp = mul_pow5_div_pow2(mp, i as u32, j);
104         vm = mul_pow5_div_pow2(mm, i as u32, j);
105         if q != 0 && (vp - 1) / 10 <= vm / 10 {
106             j = q as i32 - 1 - (pow5bits(i + 1) - FLOAT_POW5_BITCOUNT);
107             last_removed_digit = (mul_pow5_div_pow2(mv, (i + 1) as u32, j) % 10) as u8;
108         }
109         if q <= 1 {
110             // {vr,vp,vm} is trailing zeros if {mv,mp,mm} has at least q trailing 0 bits.
111             // mv = 4 * m2, so it always has at least two trailing 0 bits.
112             vr_is_trailing_zeros = true;
113             if accept_bounds {
114                 // mm = mv - 1 - mm_shift, so it has 1 trailing 0 bit iff mm_shift == 1.
115                 vm_is_trailing_zeros = mm_shift == 1;
116             } else {
117                 // mp = mv + 2, so it always has at least one trailing 0 bit.
118                 vp -= 1;
119             }
120         } else if q < 31 {
121             // TODO(ulfjack): Use a tighter bound here.
122             vr_is_trailing_zeros = multiple_of_power_of_2_32(mv, q - 1);
123         }
124     }
125 
126     // Step 4: Find the shortest decimal representation in the interval of valid representations.
127     let mut removed = 0i32;
128     let output = if vm_is_trailing_zeros || vr_is_trailing_zeros {
129         // General case, which happens rarely (~4.0%).
130         while vp / 10 > vm / 10 {
131             vm_is_trailing_zeros &= vm - (vm / 10) * 10 == 0;
132             vr_is_trailing_zeros &= last_removed_digit == 0;
133             last_removed_digit = (vr % 10) as u8;
134             vr /= 10;
135             vp /= 10;
136             vm /= 10;
137             removed += 1;
138         }
139         if vm_is_trailing_zeros {
140             while vm % 10 == 0 {
141                 vr_is_trailing_zeros &= last_removed_digit == 0;
142                 last_removed_digit = (vr % 10) as u8;
143                 vr /= 10;
144                 vp /= 10;
145                 vm /= 10;
146                 removed += 1;
147             }
148         }
149         if vr_is_trailing_zeros && last_removed_digit == 5 && vr % 2 == 0 {
150             // Round even if the exact number is .....50..0.
151             last_removed_digit = 4;
152         }
153         // We need to take vr + 1 if vr is outside bounds or we need to round up.
154         vr + ((vr == vm && (!accept_bounds || !vm_is_trailing_zeros)) || last_removed_digit >= 5)
155             as u32
156     } else {
157         // Specialized for the common case (~96.0%). Percentages below are relative to this.
158         // Loop iterations below (approximately):
159         // 0: 13.6%, 1: 70.7%, 2: 14.1%, 3: 1.39%, 4: 0.14%, 5+: 0.01%
160         while vp / 10 > vm / 10 {
161             last_removed_digit = (vr % 10) as u8;
162             vr /= 10;
163             vp /= 10;
164             vm /= 10;
165             removed += 1;
166         }
167         // We need to take vr + 1 if vr is outside bounds or we need to round up.
168         vr + (vr == vm || last_removed_digit >= 5) as u32
169     };
170     let exp = e10 + removed;
171 
172     FloatingDecimal32 {
173         exponent: exp,
174         mantissa: output,
175     }
176 }
177