@@ -6,73 +6,160 @@ use super::Float;
66
77/// Conversions from integers to floats.
88///
9- /// These are hand-optimized bit twiddling code,
10- /// which unfortunately isn't the easiest kind of code to read.
9+ /// The algorithm is explained here: <https://blog.m-ou.se/floats/>. It roughly does the following:
10+ /// - Calculate a base mantissa by shifting the integer into mantissa position. This gives us a
11+ /// mantissa _with the implicit bit set_!
12+ /// - Figure out if rounding needs to occour by classifying truncated bits. Some patterns are used
13+ /// to simplify this. Adjust the mantissa if needed.
14+ /// - Calculate the exponent based on the base-2 logarithm of `i` (leading zeros) and subtract one.
15+ /// - Shift the exponent and add the mantissa to create the final representation. Subtracting one
16+ /// from the exponent (above) accounts for the explicit bit being set in the mantissa.
1117///
12- /// The algorithm is explained here: <https://blog.m-ou.se/floats/>
18+ /// # Terminology
19+ ///
20+ /// - `i`: the original integer
21+ /// - `i_m`: the integer, shifted fully left (no leading zeros)
22+ /// - `n`: number of leading zeroes
23+ /// - `e`: the resulting exponent. Usually 1 is subtracted to offset the mantissa implicit bit.
24+ /// - `m`: the resulting mantissa
25+ /// - `m_base`: the mantissa before adjusting for truncated bits. Implicit bit is usually set.
26+ /// - `adj`: the bits that will be truncated, possibly compressed in some way.
1327mod int_to_float {
28+ use super :: * ;
29+
30+ /// Calculate the exponent from the number of leading zeros.
31+ ///
32+ /// Usually 1 is subtracted from this function's result, so that a mantissa with the implicit
33+ /// bit set can be added back later.
34+ fn exp < I : Int , F : Float < Int : CastFrom < u32 > > > ( n : u32 ) -> F :: Int {
35+ F :: Int :: cast_from ( F :: EXPONENT_BIAS - 1 + I :: BITS - n)
36+ }
37+
38+ /// Adjust a mantissa with dropped bits to perform correct rounding.
39+ ///
40+ /// The dropped bits should be exactly the bits that get truncated (left-aligned),
41+ fn m_adj < F : Float > ( m_base : F :: Int , dropped_bits : F :: Int ) -> F :: Int {
42+ // Branchlessly extract a `1` if rounding up should happen, 0 otherwise
43+ // This accounts for rounding to even.
44+ let adj = ( dropped_bits - ( dropped_bits >> ( F :: BITS - 1 ) & !m_base) ) >> ( F :: BITS - 1 ) ;
45+
46+ // Add one when we need to round up. Break ties to even.
47+ m_base + adj
48+ }
49+
50+ /// Shift the exponent to its position and add the mantissa.
51+ ///
52+ /// If the mantissa has the implicit bit set, the exponent should be one less than its actual
53+ /// value to cancel it out.
54+ fn repr < F : Float > ( e : F :: Int , m : F :: Int ) -> F :: Int {
55+ // + rather than | so the mantissa can overflow into the exponent
56+ ( e << F :: SIGNIFICAND_BITS ) + m
57+ }
58+
59+ /// Shift distance for a left-aligned integer to a smaller float.
60+ fn shift_f_lt_i < I : Int , F : Float > ( ) -> u32 {
61+ ( I :: BITS - F :: BITS ) + F :: EXPONENT_BITS
62+ }
63+
64+ /// Shift distance for an integer with `n` leading zeros to a smaller float.
65+ fn shift_f_gt_i < I : Int , F : Float > ( n : u32 ) -> u32 {
66+ F :: SIGNIFICAND_BITS - I :: BITS + 1 + n
67+ }
68+
69+ /// Perform a signed operation as unsigned, then add the sign back.
70+ pub fn signed < I , F , Conv > ( i : I , conv : Conv ) -> F
71+ where
72+ F : Float ,
73+ I : Int ,
74+ F :: Int : CastFrom < I > ,
75+ Conv : Fn ( I :: UnsignedInt ) -> F :: Int ,
76+ {
77+ let sign_bit = F :: Int :: cast_from ( i >> ( I :: BITS - 1 ) ) << ( F :: BITS - 1 ) ;
78+ F :: from_repr ( conv ( i. unsigned_abs ( ) ) | sign_bit)
79+ }
80+
1481 pub fn u32_to_f32_bits ( i : u32 ) -> u32 {
1582 if i == 0 {
1683 return 0 ;
1784 }
1885 let n = i. leading_zeros ( ) ;
19- let a = ( i << n) >> 8 ; // Significant bits, with bit 24 still in tact.
20- let b = ( i << n) << 24 ; // Insignificant bits, only relevant for rounding.
21- let m = a + ( ( b - ( b >> 31 & !a) ) >> 31 ) ; // Add one when we need to round up. Break ties to even.
22- let e = 157 - n; // Exponent plus 127, minus one.
23- ( e << 23 ) + m // + not |, so the mantissa can overflow into the exponent.
86+ // Mantissa with implicit bit set
87+ let m_base = ( i << n) >> f32:: EXPONENT_BITS ;
88+ // Bits that will be dropped
89+ let adj = ( i << n) << ( f32:: SIGNIFICAND_BITS + 1 ) ;
90+ let m = m_adj :: < f32 > ( m_base, adj) ;
91+ let e = exp :: < u32 , f32 > ( n) - 1 ;
92+ repr :: < f32 > ( e, m)
2493 }
2594
2695 pub fn u32_to_f64_bits ( i : u32 ) -> u64 {
2796 if i == 0 {
2897 return 0 ;
2998 }
3099 let n = i. leading_zeros ( ) ;
31- let m = ( i as u64 ) << ( 21 + n) ; // Significant bits, with bit 53 still in tact.
32- let e = 1053 - n as u64 ; // Exponent plus 1023, minus one.
33- ( e << 52 ) + m // Bit 53 of m will overflow into e.
100+ // Mantissa with implicit bit set
101+ let m = ( i as u64 ) << shift_f_gt_i :: < u32 , f64 > ( n) ;
102+ let e = exp :: < u32 , f64 > ( n) - 1 ;
103+ repr :: < f64 > ( e, m)
34104 }
35105
36106 pub fn u64_to_f32_bits ( i : u64 ) -> u32 {
37107 let n = i. leading_zeros ( ) ;
38- let y = i. wrapping_shl ( n) ;
39- let a = ( y >> 40 ) as u32 ; // Significant bits, with bit 24 still in tact.
40- let b = ( y >> 8 | y & 0xFFFF ) as u32 ; // Insignificant bits, only relevant for rounding.
41- let m = a + ( ( b - ( b >> 31 & !a) ) >> 31 ) ; // Add one when we need to round up. Break ties to even.
42- let e = if i == 0 { 0 } else { 189 - n } ; // Exponent plus 127, minus one, except for zero.
43- ( e << 23 ) + m // + not |, so the mantissa can overflow into the exponent.
108+ let i_m = i. wrapping_shl ( n) ;
109+ // Mantissa with implicit bit set
110+ let m_base: u32 = ( i_m >> shift_f_lt_i :: < u64 , f32 > ( ) ) as u32 ;
111+ // The entire lower half of `i` will be truncated (masked portion), plus the
112+ // next `EXPONENT_BITS` bits.
113+ let adj = ( i_m >> f32:: EXPONENT_BITS | i_m & 0xFFFF ) as u32 ;
114+ let m = m_adj :: < f32 > ( m_base, adj) ;
115+ let e = if i == 0 { 0 } else { exp :: < u64 , f32 > ( n) - 1 } ;
116+ repr :: < f32 > ( e, m)
44117 }
45118
46119 pub fn u64_to_f64_bits ( i : u64 ) -> u64 {
47120 if i == 0 {
48121 return 0 ;
49122 }
50123 let n = i. leading_zeros ( ) ;
51- let a = ( i << n) >> 11 ; // Significant bits, with bit 53 still in tact.
52- let b = ( i << n) << 53 ; // Insignificant bits, only relevant for rounding.
53- let m = a + ( ( b - ( b >> 63 & !a) ) >> 63 ) ; // Add one when we need to round up. Break ties to even.
54- let e = 1085 - n as u64 ; // Exponent plus 1023, minus one.
55- ( e << 52 ) + m // + not |, so the mantissa can overflow into the exponent.
124+ // Mantissa with implicit bit set
125+ let m_base = ( i << n) >> f64:: EXPONENT_BITS ;
126+ let adj = ( i << n) << ( f64:: SIGNIFICAND_BITS + 1 ) ;
127+ let m = m_adj :: < f64 > ( m_base, adj) ;
128+ let e = exp :: < u64 , f64 > ( n) - 1 ;
129+ repr :: < f64 > ( e, m)
56130 }
57131
58132 pub fn u128_to_f32_bits ( i : u128 ) -> u32 {
59133 let n = i. leading_zeros ( ) ;
60- let y = i. wrapping_shl ( n) ;
61- let a = ( y >> 104 ) as u32 ; // Significant bits, with bit 24 still in tact.
62- let b = ( y >> 72 ) as u32 | ( ( y << 32 ) >> 32 != 0 ) as u32 ; // Insignificant bits, only relevant for rounding.
63- let m = a + ( ( b - ( b >> 31 & !a) ) >> 31 ) ; // Add one when we need to round up. Break ties to even.
64- let e = if i == 0 { 0 } else { 253 - n } ; // Exponent plus 127, minus one, except for zero.
65- ( e << 23 ) + m // + not |, so the mantissa can overflow into the exponent.
134+ let i_m = i. wrapping_shl ( n) ; // Mantissa, shifted so the first bit is nonzero
135+ let m_base: u32 = ( i_m >> shift_f_lt_i :: < u128 , f32 > ( ) ) as u32 ;
136+
137+ // Within the upper `F::BITS`, everything except for the signifcand
138+ // gets truncated
139+ let d1: u32 = ( i_m >> ( u128:: BITS - f32:: BITS - f32:: SIGNIFICAND_BITS - 1 ) ) . cast ( ) ;
140+
141+ // The entire rest of `i_m` gets truncated. Zero the upper `F::BITS` then just
142+ // check if it is nonzero.
143+ let d2: u32 = ( i_m << f32:: BITS >> f32:: BITS != 0 ) . into ( ) ;
144+ let adj = d1 | d2;
145+
146+ // Mantissa with implicit bit set
147+ let m = m_adj :: < f32 > ( m_base, adj) ;
148+ let e = if i == 0 { 0 } else { exp :: < u128 , f32 > ( n) - 1 } ;
149+ repr :: < f32 > ( e, m)
66150 }
67151
68152 pub fn u128_to_f64_bits ( i : u128 ) -> u64 {
69153 let n = i. leading_zeros ( ) ;
70- let y = i. wrapping_shl ( n) ;
71- let a = ( y >> 75 ) as u64 ; // Significant bits, with bit 53 still in tact.
72- let b = ( y >> 11 | y & 0xFFFF_FFFF ) as u64 ; // Insignificant bits, only relevant for rounding.
73- let m = a + ( ( b - ( b >> 63 & !a) ) >> 63 ) ; // Add one when we need to round up. Break ties to even.
74- let e = if i == 0 { 0 } else { 1149 - n as u64 } ; // Exponent plus 1023, minus one, except for zero.
75- ( e << 52 ) + m // + not |, so the mantissa can overflow into the exponent.
154+ let i_m = i. wrapping_shl ( n) ;
155+ // Mantissa with implicit bit set
156+ let m_base: u64 = ( i_m >> shift_f_lt_i :: < u128 , f64 > ( ) ) as u64 ;
157+ // The entire lower half of `i` will be truncated (masked portion), plus the
158+ // next `EXPONENT_BITS` bits.
159+ let adj = ( i_m >> f64:: EXPONENT_BITS | i_m & 0xFFFF_FFFF ) as u64 ;
160+ let m = m_adj :: < f64 > ( m_base, adj) ;
161+ let e = if i == 0 { 0 } else { exp :: < u128 , f64 > ( n) - 1 } ;
162+ repr :: < f64 > ( e, m)
76163 }
77164}
78165
@@ -113,38 +200,32 @@ intrinsics! {
113200intrinsics ! {
114201 #[ arm_aeabi_alias = __aeabi_i2f]
115202 pub extern "C" fn __floatsisf( i: i32 ) -> f32 {
116- let sign_bit = ( ( i >> 31 ) as u32 ) << 31 ;
117- f32 :: from_bits( int_to_float:: u32_to_f32_bits( i. unsigned_abs( ) ) | sign_bit)
203+ int_to_float:: signed( i, int_to_float:: u32_to_f32_bits)
118204 }
119205
120206 #[ arm_aeabi_alias = __aeabi_i2d]
121207 pub extern "C" fn __floatsidf( i: i32 ) -> f64 {
122- let sign_bit = ( ( i >> 31 ) as u64 ) << 63 ;
123- f64 :: from_bits( int_to_float:: u32_to_f64_bits( i. unsigned_abs( ) ) | sign_bit)
208+ int_to_float:: signed( i, int_to_float:: u32_to_f64_bits)
124209 }
125210
126211 #[ arm_aeabi_alias = __aeabi_l2f]
127212 pub extern "C" fn __floatdisf( i: i64 ) -> f32 {
128- let sign_bit = ( ( i >> 63 ) as u32 ) << 31 ;
129- f32 :: from_bits( int_to_float:: u64_to_f32_bits( i. unsigned_abs( ) ) | sign_bit)
213+ int_to_float:: signed( i, int_to_float:: u64_to_f32_bits)
130214 }
131215
132216 #[ arm_aeabi_alias = __aeabi_l2d]
133217 pub extern "C" fn __floatdidf( i: i64 ) -> f64 {
134- let sign_bit = ( ( i >> 63 ) as u64 ) << 63 ;
135- f64 :: from_bits( int_to_float:: u64_to_f64_bits( i. unsigned_abs( ) ) | sign_bit)
218+ int_to_float:: signed( i, int_to_float:: u64_to_f64_bits)
136219 }
137220
138221 #[ cfg_attr( target_os = "uefi" , unadjusted_on_win64) ]
139222 pub extern "C" fn __floattisf( i: i128 ) -> f32 {
140- let sign_bit = ( ( i >> 127 ) as u32 ) << 31 ;
141- f32 :: from_bits( int_to_float:: u128_to_f32_bits( i. unsigned_abs( ) ) | sign_bit)
223+ int_to_float:: signed( i, int_to_float:: u128_to_f32_bits)
142224 }
143225
144226 #[ cfg_attr( target_os = "uefi" , unadjusted_on_win64) ]
145227 pub extern "C" fn __floattidf( i: i128 ) -> f64 {
146- let sign_bit = ( ( i >> 127 ) as u64 ) << 63 ;
147- f64 :: from_bits( int_to_float:: u128_to_f64_bits( i. unsigned_abs( ) ) | sign_bit)
228+ int_to_float:: signed( i, int_to_float:: u128_to_f64_bits)
148229 }
149230}
150231
0 commit comments