1/* 2 * QEMU float support 3 * 4 * The code in this source file is derived from release 2a of the SoftFloat 5 * IEC/IEEE Floating-point Arithmetic Package. Those parts of the code (and 6 * some later contributions) are provided under that license, as detailed below. 7 * It has subsequently been modified by contributors to the QEMU Project, 8 * so some portions are provided under: 9 * the SoftFloat-2a license 10 * the BSD license 11 * GPL-v2-or-later 12 * 13 * Any future contributions to this file after December 1st 2014 will be 14 * taken to be licensed under the Softfloat-2a license unless specifically 15 * indicated otherwise. 16 */ 17 18static void partsN(return_nan)(FloatPartsN *a, float_status *s) 19{ 20 switch (a->cls) { 21 case float_class_snan: 22 float_raise(float_flag_invalid, s); 23 if (s->default_nan_mode) { 24 parts_default_nan(a, s); 25 } else { 26 parts_silence_nan(a, s); 27 } 28 break; 29 case float_class_qnan: 30 if (s->default_nan_mode) { 31 parts_default_nan(a, s); 32 } 33 break; 34 default: 35 g_assert_not_reached(); 36 } 37} 38 39static FloatPartsN *partsN(pick_nan)(FloatPartsN *a, FloatPartsN *b, 40 float_status *s) 41{ 42 if (is_snan(a->cls) || is_snan(b->cls)) { 43 float_raise(float_flag_invalid, s); 44 } 45 46 if (s->default_nan_mode) { 47 parts_default_nan(a, s); 48 } else { 49 int cmp = frac_cmp(a, b); 50 if (cmp == 0) { 51 cmp = a->sign < b->sign; 52 } 53 54 if (pickNaN(a->cls, b->cls, cmp > 0, s)) { 55 a = b; 56 } 57 if (is_snan(a->cls)) { 58 parts_silence_nan(a, s); 59 } 60 } 61 return a; 62} 63 64static FloatPartsN *partsN(pick_nan_muladd)(FloatPartsN *a, FloatPartsN *b, 65 FloatPartsN *c, float_status *s, 66 int ab_mask, int abc_mask) 67{ 68 int which; 69 70 if (unlikely(abc_mask & float_cmask_snan)) { 71 float_raise(float_flag_invalid, s); 72 } 73 74 which = pickNaNMulAdd(a->cls, b->cls, c->cls, 75 ab_mask == float_cmask_infzero, s); 76 77 if (s->default_nan_mode || which == 3) { 78 /* 79 * Note that this check is after pickNaNMulAdd so that function 80 * has an opportunity to set the Invalid flag for infzero. 81 */ 82 parts_default_nan(a, s); 83 return a; 84 } 85 86 switch (which) { 87 case 0: 88 break; 89 case 1: 90 a = b; 91 break; 92 case 2: 93 a = c; 94 break; 95 default: 96 g_assert_not_reached(); 97 } 98 if (is_snan(a->cls)) { 99 parts_silence_nan(a, s); 100 } 101 return a; 102} 103 104/* 105 * Canonicalize the FloatParts structure. Determine the class, 106 * unbias the exponent, and normalize the fraction. 107 */ 108static void partsN(canonicalize)(FloatPartsN *p, float_status *status, 109 const FloatFmt *fmt) 110{ 111 if (unlikely(p->exp == 0)) { 112 if (likely(frac_eqz(p))) { 113 p->cls = float_class_zero; 114 } else if (status->flush_inputs_to_zero) { 115 float_raise(float_flag_input_denormal, status); 116 p->cls = float_class_zero; 117 frac_clear(p); 118 } else { 119 int shift = frac_normalize(p); 120 p->cls = float_class_normal; 121 p->exp = fmt->frac_shift - fmt->exp_bias - shift + 1; 122 } 123 } else if (likely(p->exp < fmt->exp_max) || fmt->arm_althp) { 124 p->cls = float_class_normal; 125 p->exp -= fmt->exp_bias; 126 frac_shl(p, fmt->frac_shift); 127 p->frac_hi |= DECOMPOSED_IMPLICIT_BIT; 128 } else if (likely(frac_eqz(p))) { 129 p->cls = float_class_inf; 130 } else { 131 frac_shl(p, fmt->frac_shift); 132 p->cls = (parts_is_snan_frac(p->frac_hi, status) 133 ? float_class_snan : float_class_qnan); 134 } 135} 136 137/* 138 * Round and uncanonicalize a floating-point number by parts. There 139 * are FRAC_SHIFT bits that may require rounding at the bottom of the 140 * fraction; these bits will be removed. The exponent will be biased 141 * by EXP_BIAS and must be bounded by [EXP_MAX-1, 0]. 142 */ 143static void partsN(uncanon)(FloatPartsN *p, float_status *s, 144 const FloatFmt *fmt) 145{ 146 const int exp_max = fmt->exp_max; 147 const int frac_shift = fmt->frac_shift; 148 const uint64_t frac_lsb = fmt->frac_lsb; 149 const uint64_t frac_lsbm1 = fmt->frac_lsbm1; 150 const uint64_t round_mask = fmt->round_mask; 151 const uint64_t roundeven_mask = fmt->roundeven_mask; 152 uint64_t inc; 153 bool overflow_norm; 154 int exp, flags = 0; 155 156 if (unlikely(p->cls != float_class_normal)) { 157 switch (p->cls) { 158 case float_class_zero: 159 p->exp = 0; 160 frac_clear(p); 161 return; 162 case float_class_inf: 163 g_assert(!fmt->arm_althp); 164 p->exp = fmt->exp_max; 165 frac_clear(p); 166 return; 167 case float_class_qnan: 168 case float_class_snan: 169 g_assert(!fmt->arm_althp); 170 p->exp = fmt->exp_max; 171 frac_shr(p, fmt->frac_shift); 172 return; 173 default: 174 break; 175 } 176 g_assert_not_reached(); 177 } 178 179 switch (s->float_rounding_mode) { 180 case float_round_nearest_even: 181 overflow_norm = false; 182 inc = ((p->frac_lo & roundeven_mask) != frac_lsbm1 ? frac_lsbm1 : 0); 183 break; 184 case float_round_ties_away: 185 overflow_norm = false; 186 inc = frac_lsbm1; 187 break; 188 case float_round_to_zero: 189 overflow_norm = true; 190 inc = 0; 191 break; 192 case float_round_up: 193 inc = p->sign ? 0 : round_mask; 194 overflow_norm = p->sign; 195 break; 196 case float_round_down: 197 inc = p->sign ? round_mask : 0; 198 overflow_norm = !p->sign; 199 break; 200 case float_round_to_odd: 201 overflow_norm = true; 202 inc = p->frac_lo & frac_lsb ? 0 : round_mask; 203 break; 204 default: 205 g_assert_not_reached(); 206 } 207 208 exp = p->exp + fmt->exp_bias; 209 if (likely(exp > 0)) { 210 if (p->frac_lo & round_mask) { 211 flags |= float_flag_inexact; 212 if (frac_addi(p, p, inc)) { 213 frac_shr(p, 1); 214 p->frac_hi |= DECOMPOSED_IMPLICIT_BIT; 215 exp++; 216 } 217 } 218 frac_shr(p, frac_shift); 219 220 if (fmt->arm_althp) { 221 /* ARM Alt HP eschews Inf and NaN for a wider exponent. */ 222 if (unlikely(exp > exp_max)) { 223 /* Overflow. Return the maximum normal. */ 224 flags = float_flag_invalid; 225 exp = exp_max; 226 frac_allones(p); 227 } 228 } else if (unlikely(exp >= exp_max)) { 229 flags |= float_flag_overflow | float_flag_inexact; 230 if (overflow_norm) { 231 exp = exp_max - 1; 232 frac_allones(p); 233 } else { 234 p->cls = float_class_inf; 235 exp = exp_max; 236 frac_clear(p); 237 } 238 } 239 } else if (s->flush_to_zero) { 240 flags |= float_flag_output_denormal; 241 p->cls = float_class_zero; 242 exp = 0; 243 frac_clear(p); 244 } else { 245 bool is_tiny = s->tininess_before_rounding || exp < 0; 246 247 if (!is_tiny) { 248 FloatPartsN discard; 249 is_tiny = !frac_addi(&discard, p, inc); 250 } 251 252 frac_shrjam(p, 1 - exp); 253 254 if (p->frac_lo & round_mask) { 255 /* Need to recompute round-to-even/round-to-odd. */ 256 switch (s->float_rounding_mode) { 257 case float_round_nearest_even: 258 inc = ((p->frac_lo & roundeven_mask) != frac_lsbm1 259 ? frac_lsbm1 : 0); 260 break; 261 case float_round_to_odd: 262 inc = p->frac_lo & frac_lsb ? 0 : round_mask; 263 break; 264 default: 265 break; 266 } 267 flags |= float_flag_inexact; 268 frac_addi(p, p, inc); 269 } 270 271 exp = (p->frac_hi & DECOMPOSED_IMPLICIT_BIT) != 0; 272 frac_shr(p, frac_shift); 273 274 if (is_tiny && (flags & float_flag_inexact)) { 275 flags |= float_flag_underflow; 276 } 277 if (exp == 0 && frac_eqz(p)) { 278 p->cls = float_class_zero; 279 } 280 } 281 p->exp = exp; 282 float_raise(flags, s); 283} 284