const_soft_float/soft_f64/
mul.rs1use crate::soft_f64::{u64_widen_mul, SoftF64};
2
3type F = SoftF64;
4
5type FInt = u64;
6
7const fn widen_mul(a: FInt, b: FInt) -> (FInt, FInt) {
8 u64_widen_mul(a, b)
9}
10
11pub(crate) const fn mul(a: F, b: F) -> F {
12 let one: FInt = 1;
13 let zero: FInt = 0;
14
15 let bits = F::BITS;
16 let significand_bits = F::SIGNIFICAND_BITS;
17 let max_exponent = F::EXPONENT_MAX;
18
19 let exponent_bias = F::EXPONENT_BIAS;
20
21 let implicit_bit = F::IMPLICIT_BIT;
22 let significand_mask = F::SIGNIFICAND_MASK;
23 let sign_bit = F::SIGN_MASK as FInt;
24 let abs_mask = sign_bit - one;
25 let exponent_mask = F::EXPONENT_MASK;
26 let inf_rep = exponent_mask;
27 let quiet_bit = implicit_bit >> 1;
28 let qnan_rep = exponent_mask | quiet_bit;
29 let exponent_bits = F::EXPONENT_BITS;
30
31 let a_rep = a.repr();
32 let b_rep = b.repr();
33
34 let a_exponent = (a_rep >> significand_bits) & max_exponent as FInt;
35 let b_exponent = (b_rep >> significand_bits) & max_exponent as FInt;
36 let product_sign = (a_rep ^ b_rep) & sign_bit;
37
38 let mut a_significand = a_rep & significand_mask;
39 let mut b_significand = b_rep & significand_mask;
40 let mut scale = 0;
41
42 if a_exponent.wrapping_sub(one) >= (max_exponent - 1) as FInt
44 || b_exponent.wrapping_sub(one) >= (max_exponent - 1) as FInt
45 {
46 let a_abs = a_rep & abs_mask;
47 let b_abs = b_rep & abs_mask;
48
49 if a_abs > inf_rep {
51 return F::from_repr(a_rep | quiet_bit);
52 }
53 if b_abs > inf_rep {
55 return F::from_repr(b_rep | quiet_bit);
56 }
57
58 if a_abs == inf_rep {
59 if b_abs != zero {
60 return F::from_repr(a_abs | product_sign);
62 } else {
63 return F::from_repr(qnan_rep);
65 }
66 }
67
68 if b_abs == inf_rep {
69 if a_abs != zero {
70 return F::from_repr(b_abs | product_sign);
72 } else {
73 return F::from_repr(qnan_rep);
75 }
76 }
77
78 if a_abs == zero {
80 return F::from_repr(product_sign);
81 }
82
83 if b_abs == zero {
85 return F::from_repr(product_sign);
86 }
87
88 if a_abs < implicit_bit {
92 let (exponent, significand) = F::normalize(a_significand);
93 scale += exponent;
94 a_significand = significand;
95 }
96
97 if b_abs < implicit_bit {
98 let (exponent, significand) = F::normalize(b_significand);
99 scale += exponent;
100 b_significand = significand;
101 }
102 }
103
104 a_significand |= implicit_bit;
108 b_significand |= implicit_bit;
109
110 let (mut product_low, mut product_high) =
116 widen_mul(a_significand, b_significand << exponent_bits);
117
118 let a_exponent_i32: i32 = a_exponent as _;
119 let b_exponent_i32: i32 = b_exponent as _;
120 let mut product_exponent: i32 = a_exponent_i32
121 .wrapping_add(b_exponent_i32)
122 .wrapping_add(scale)
123 .wrapping_sub(exponent_bias as i32);
124
125 if (product_high & implicit_bit) != zero {
127 product_exponent = product_exponent.wrapping_add(1);
128 } else {
129 product_high = (product_high << 1) | (product_low >> (bits - 1));
130 product_low <<= 1;
131 }
132
133 if product_exponent >= max_exponent as i32 {
135 return F::from_repr(inf_rep | product_sign);
136 }
137
138 if product_exponent <= 0 {
139 let shift = one.wrapping_sub(product_exponent as FInt) as u32;
146 if shift >= bits {
147 return F::from_repr(product_sign);
148 }
149
150 if shift < bits {
153 let sticky = product_low << (bits - shift);
154 product_low = product_high << (bits - shift) | product_low >> shift | sticky;
155 product_high >>= shift;
156 } else if shift < (2 * bits) {
157 let sticky = product_high << (2 * bits - shift) | product_low;
158 product_low = product_high >> (shift - bits) | sticky;
159 product_high = zero;
160 } else {
161 product_high = zero;
162 }
163 } else {
164 product_high &= significand_mask;
166 product_high |= (product_exponent as FInt) << significand_bits;
167 }
168
169 product_high |= product_sign;
171
172 if product_low > sign_bit {
176 product_high += one;
177 }
178
179 if product_low == sign_bit {
180 product_high += product_high & one;
181 }
182
183 F::from_repr(product_high)
184}
185
186#[cfg(test)]
187mod test {
188 use crate::soft_f64::SoftF64;
189
190 #[test]
191 fn sanity_check() {
192 assert_eq!(SoftF64(2.0).mul(SoftF64(2.0)).0, 4.0)
193 }
194}