?  search:
summary | log | graphiclog1 | graphiclog2 | commit | commitdiff | tree | refs | edit | fork
blame | history | raw | HEAD
hppa: Cleanup libm-test-ulps.
[glibc.git] / soft-fp / op-common.h
bloba9d94d7f1648629319dcae64ba03d5dcf9128f9d
1 /* Software floating-point emulation. Common operations.
2 Copyright (C) 1997-2013 Free Software Foundation, Inc.
3 This file is part of the GNU C Library.
4 Contributed by Richard Henderson (rth@cygnus.com),
5 Jakub Jelinek (jj@ultra.linux.cz),
6 David S. Miller (davem@redhat.com) and
7 Peter Maydell (pmaydell@chiark.greenend.org.uk).
8
9 The GNU C Library is free software; you can redistribute it and/or
10 modify it under the terms of the GNU Lesser General Public
11 License as published by the Free Software Foundation; either
12 version 2.1 of the License, or (at your option) any later version.
13
14 In addition to the permissions in the GNU Lesser General Public
15 License, the Free Software Foundation gives you unlimited
16 permission to link the compiled version of this file into
17 combinations with other programs, and to distribute those
18 combinations without any restriction coming from the use of this
19 file. (The Lesser General Public License restrictions do apply in
20 other respects; for example, they cover modification of the file,
21 and distribution when not linked into a combine executable.)
22
23 The GNU C Library is distributed in the hope that it will be useful,
24 but WITHOUT ANY WARRANTY; without even the implied warranty of
25 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
26 Lesser General Public License for more details.
27
28 You should have received a copy of the GNU Lesser General Public
29 License along with the GNU C Library; if not, see
30 <http://www.gnu.org/licenses/>. */
31
32 #define _FP_DECL(wc, X) \
33 _FP_I_TYPE X##_c __attribute__((unused)); \
34 _FP_I_TYPE X##_s __attribute__((unused)); \
35 _FP_I_TYPE X##_e; \
36 _FP_FRAC_DECL_##wc(X)
37
38 /*
39 * Finish truely unpacking a native fp value by classifying the kind
40 * of fp value and normalizing both the exponent and the fraction.
41 */
42
43 #define _FP_UNPACK_CANONICAL(fs, wc, X) \
44 do { \
45 switch (X##_e) \
46 { \
47 default: \
48 _FP_FRAC_HIGH_RAW_##fs(X) |= _FP_IMPLBIT_##fs; \
49 _FP_FRAC_SLL_##wc(X, _FP_WORKBITS); \
50 X##_e -= _FP_EXPBIAS_##fs; \
51 X##_c = FP_CLS_NORMAL; \
52 break; \
53 \
54 case 0: \
55 if (_FP_FRAC_ZEROP_##wc(X)) \
56 X##_c = FP_CLS_ZERO; \
57 else \
58 { \
59 /* a denormalized number */ \
60 _FP_I_TYPE _shift; \
61 _FP_FRAC_CLZ_##wc(_shift, X); \
62 _shift -= _FP_FRACXBITS_##fs; \
63 _FP_FRAC_SLL_##wc(X, (_shift+_FP_WORKBITS)); \
64 X##_e -= _FP_EXPBIAS_##fs - 1 + _shift; \
65 X##_c = FP_CLS_NORMAL; \
66 FP_SET_EXCEPTION(FP_EX_DENORM); \
67 } \
68 break; \
69 \
70 case _FP_EXPMAX_##fs: \
71 if (_FP_FRAC_ZEROP_##wc(X)) \
72 X##_c = FP_CLS_INF; \
73 else \
74 { \
75 X##_c = FP_CLS_NAN; \
76 /* Check for signaling NaN */ \
77 if (!(_FP_FRAC_HIGH_RAW_##fs(X) & _FP_QNANBIT_##fs)) \
78 FP_SET_EXCEPTION(FP_EX_INVALID); \
79 } \
80 break; \
81 } \
82 } while (0)
83
84 /* Finish unpacking an fp value in semi-raw mode: the mantissa is
85 shifted by _FP_WORKBITS but the implicit MSB is not inserted and
86 other classification is not done. */
87 #define _FP_UNPACK_SEMIRAW(fs, wc, X) _FP_FRAC_SLL_##wc(X, _FP_WORKBITS)
88
89 /* A semi-raw value has overflowed to infinity. Adjust the mantissa
90 and exponent appropriately. */
91 #define _FP_OVERFLOW_SEMIRAW(fs, wc, X) \
92 do { \
93 if (FP_ROUNDMODE == FP_RND_NEAREST \
94 || (FP_ROUNDMODE == FP_RND_PINF && !X##_s) \
95 || (FP_ROUNDMODE == FP_RND_MINF && X##_s)) \
96 { \
97 X##_e = _FP_EXPMAX_##fs; \
98 _FP_FRAC_SET_##wc(X, _FP_ZEROFRAC_##wc); \
99 } \
100 else \
101 { \
102 X##_e = _FP_EXPMAX_##fs - 1; \
103 _FP_FRAC_SET_##wc(X, _FP_MAXFRAC_##wc); \
104 } \
105 FP_SET_EXCEPTION(FP_EX_INEXACT); \
106 FP_SET_EXCEPTION(FP_EX_OVERFLOW); \
107 } while (0)
108
109 /* Check for a semi-raw value being a signaling NaN and raise the
110 invalid exception if so. */
111 #define _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, X) \
112 do { \
113 if (X##_e == _FP_EXPMAX_##fs \
114 && !_FP_FRAC_ZEROP_##wc(X) \
115 && !(_FP_FRAC_HIGH_##fs(X) & _FP_QNANBIT_SH_##fs)) \
116 FP_SET_EXCEPTION(FP_EX_INVALID); \
117 } while (0)
118
119 /* Choose a NaN result from an operation on two semi-raw NaN
120 values. */
121 #define _FP_CHOOSENAN_SEMIRAW(fs, wc, R, X, Y, OP) \
122 do { \
123 /* _FP_CHOOSENAN expects raw values, so shift as required. */ \
124 _FP_FRAC_SRL_##wc(X, _FP_WORKBITS); \
125 _FP_FRAC_SRL_##wc(Y, _FP_WORKBITS); \
126 _FP_CHOOSENAN(fs, wc, R, X, Y, OP); \
127 _FP_FRAC_SLL_##wc(R, _FP_WORKBITS); \
128 } while (0)
129
130 /* Test whether a biased exponent is normal (not zero or maximum). */
131 #define _FP_EXP_NORMAL(fs, wc, X) (((X##_e + 1) & _FP_EXPMAX_##fs) > 1)
132
133 /* Prepare to pack an fp value in semi-raw mode: the mantissa is
134 rounded and shifted right, with the rounding possibly increasing
135 the exponent (including changing a finite value to infinity). */
136 #define _FP_PACK_SEMIRAW(fs, wc, X) \
137 do { \
138 _FP_ROUND(wc, X); \
139 if (X##_e == 0 && !_FP_FRAC_ZEROP_##wc(X)) \
140 { \
141 if ((FP_CUR_EXCEPTIONS & FP_EX_INEXACT) \
142 || (FP_TRAPPING_EXCEPTIONS & FP_EX_UNDERFLOW)) \
143 FP_SET_EXCEPTION(FP_EX_UNDERFLOW); \
144 } \
145 if (_FP_FRAC_HIGH_##fs(X) \
146 & (_FP_OVERFLOW_##fs >> 1)) \
147 { \
148 _FP_FRAC_HIGH_##fs(X) &= ~(_FP_OVERFLOW_##fs >> 1); \
149 X##_e++; \
150 if (X##_e == _FP_EXPMAX_##fs) \
151 _FP_OVERFLOW_SEMIRAW(fs, wc, X); \
152 } \
153 _FP_FRAC_SRL_##wc(X, _FP_WORKBITS); \
154 if (X##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc(X)) \
155 { \
156 if (!_FP_KEEPNANFRACP) \
157 { \
158 _FP_FRAC_SET_##wc(X, _FP_NANFRAC_##fs); \
159 X##_s = _FP_NANSIGN_##fs; \
160 } \
161 else \
162 _FP_FRAC_HIGH_RAW_##fs(X) |= _FP_QNANBIT_##fs; \
163 } \
164 } while (0)
165
166 /*
167 * Before packing the bits back into the native fp result, take care
168 * of such mundane things as rounding and overflow. Also, for some
169 * kinds of fp values, the original parts may not have been fully
170 * extracted -- but that is ok, we can regenerate them now.
171 */
172
173 #define _FP_PACK_CANONICAL(fs, wc, X) \
174 do { \
175 switch (X##_c) \
176 { \
177 case FP_CLS_NORMAL: \
178 X##_e += _FP_EXPBIAS_##fs; \
179 if (X##_e > 0) \
180 { \
181 _FP_ROUND(wc, X); \
182 if (_FP_FRAC_OVERP_##wc(fs, X)) \
183 { \
184 _FP_FRAC_CLEAR_OVERP_##wc(fs, X); \
185 X##_e++; \
186 } \
187 _FP_FRAC_SRL_##wc(X, _FP_WORKBITS); \
188 if (X##_e >= _FP_EXPMAX_##fs) \
189 { \
190 /* overflow */ \
191 switch (FP_ROUNDMODE) \
192 { \
193 case FP_RND_NEAREST: \
194 X##_c = FP_CLS_INF; \
195 break; \
196 case FP_RND_PINF: \
197 if (!X##_s) X##_c = FP_CLS_INF; \
198 break; \
199 case FP_RND_MINF: \
200 if (X##_s) X##_c = FP_CLS_INF; \
201 break; \
202 } \
203 if (X##_c == FP_CLS_INF) \
204 { \
205 /* Overflow to infinity */ \
206 X##_e = _FP_EXPMAX_##fs; \
207 _FP_FRAC_SET_##wc(X, _FP_ZEROFRAC_##wc); \
208 } \
209 else \
210 { \
211 /* Overflow to maximum normal */ \
212 X##_e = _FP_EXPMAX_##fs - 1; \
213 _FP_FRAC_SET_##wc(X, _FP_MAXFRAC_##wc); \
214 } \
215 FP_SET_EXCEPTION(FP_EX_OVERFLOW); \
216 FP_SET_EXCEPTION(FP_EX_INEXACT); \
217 } \
218 } \
219 else \
220 { \
221 /* we've got a denormalized number */ \
222 X##_e = -X##_e + 1; \
223 if (X##_e <= _FP_WFRACBITS_##fs) \
224 { \
225 _FP_FRAC_SRS_##wc(X, X##_e, _FP_WFRACBITS_##fs); \
226 _FP_ROUND(wc, X); \
227 if (_FP_FRAC_HIGH_##fs(X) \
228 & (_FP_OVERFLOW_##fs >> 1)) \
229 { \
230 X##_e = 1; \
231 _FP_FRAC_SET_##wc(X, _FP_ZEROFRAC_##wc); \
232 FP_SET_EXCEPTION(FP_EX_INEXACT); \
233 } \
234 else \
235 { \
236 X##_e = 0; \
237 _FP_FRAC_SRL_##wc(X, _FP_WORKBITS); \
238 } \
239 if ((FP_CUR_EXCEPTIONS & FP_EX_INEXACT) \
240 || (FP_TRAPPING_EXCEPTIONS & FP_EX_UNDERFLOW)) \
241 FP_SET_EXCEPTION(FP_EX_UNDERFLOW); \
242 } \
243 else \
244 { \
245 /* underflow to zero */ \
246 X##_e = 0; \
247 if (!_FP_FRAC_ZEROP_##wc(X)) \
248 { \
249 _FP_FRAC_SET_##wc(X, _FP_MINFRAC_##wc); \
250 _FP_ROUND(wc, X); \
251 _FP_FRAC_LOW_##wc(X) >>= (_FP_WORKBITS); \
252 } \
253 FP_SET_EXCEPTION(FP_EX_UNDERFLOW); \
254 } \
255 } \
256 break; \
257 \
258 case FP_CLS_ZERO: \
259 X##_e = 0; \
260 _FP_FRAC_SET_##wc(X, _FP_ZEROFRAC_##wc); \
261 break; \
262 \
263 case FP_CLS_INF: \
264 X##_e = _FP_EXPMAX_##fs; \
265 _FP_FRAC_SET_##wc(X, _FP_ZEROFRAC_##wc); \
266 break; \
267 \
268 case FP_CLS_NAN: \
269 X##_e = _FP_EXPMAX_##fs; \
270 if (!_FP_KEEPNANFRACP) \
271 { \
272 _FP_FRAC_SET_##wc(X, _FP_NANFRAC_##fs); \
273 X##_s = _FP_NANSIGN_##fs; \
274 } \
275 else \
276 _FP_FRAC_HIGH_RAW_##fs(X) |= _FP_QNANBIT_##fs; \
277 break; \
278 } \
279 } while (0)
280
281 /* This one accepts raw argument and not cooked, returns
282 * 1 if X is a signaling NaN.
283 */
284 #define _FP_ISSIGNAN(fs, wc, X) \
285 ({ \
286 int __ret = 0; \
287 if (X##_e == _FP_EXPMAX_##fs) \
288 { \
289 if (!_FP_FRAC_ZEROP_##wc(X) \
290 && !(_FP_FRAC_HIGH_RAW_##fs(X) & _FP_QNANBIT_##fs)) \
291 __ret = 1; \
292 } \
293 __ret; \
294 })
295
296
297
298
299
300 /* Addition on semi-raw values. */
301 #define _FP_ADD_INTERNAL(fs, wc, R, X, Y, OP) \
302 do { \
303 if (X##_s == Y##_s) \
304 { \
305 /* Addition. */ \
306 R##_s = X##_s; \
307 int ediff = X##_e - Y##_e; \
308 if (ediff > 0) \
309 { \
310 R##_e = X##_e; \
311 if (Y##_e == 0) \
312 { \
313 /* Y is zero or denormalized. */ \
314 if (_FP_FRAC_ZEROP_##wc(Y)) \
315 { \
316 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, X); \
317 _FP_FRAC_COPY_##wc(R, X); \
318 goto add_done; \
319 } \
320 else \
321 { \
322 FP_SET_EXCEPTION(FP_EX_DENORM); \
323 ediff--; \
324 if (ediff == 0) \
325 { \
326 _FP_FRAC_ADD_##wc(R, X, Y); \
327 goto add3; \
328 } \
329 if (X##_e == _FP_EXPMAX_##fs) \
330 { \
331 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, X); \
332 _FP_FRAC_COPY_##wc(R, X); \
333 goto add_done; \
334 } \
335 goto add1; \
336 } \
337 } \
338 else if (X##_e == _FP_EXPMAX_##fs) \
339 { \
340 /* X is NaN or Inf, Y is normal. */ \
341 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, X); \
342 _FP_FRAC_COPY_##wc(R, X); \
343 goto add_done; \
344 } \
345 \
346 /* Insert implicit MSB of Y. */ \
347 _FP_FRAC_HIGH_##fs(Y) |= _FP_IMPLBIT_SH_##fs; \
348 \
349 add1: \
350 /* Shift the mantissa of Y to the right EDIFF steps; \
351 remember to account later for the implicit MSB of X. */ \
352 if (ediff <= _FP_WFRACBITS_##fs) \
353 _FP_FRAC_SRS_##wc(Y, ediff, _FP_WFRACBITS_##fs); \
354 else if (!_FP_FRAC_ZEROP_##wc(Y)) \
355 _FP_FRAC_SET_##wc(Y, _FP_MINFRAC_##wc); \
356 _FP_FRAC_ADD_##wc(R, X, Y); \
357 } \
358 else if (ediff < 0) \
359 { \
360 ediff = -ediff; \
361 R##_e = Y##_e; \
362 if (X##_e == 0) \
363 { \
364 /* X is zero or denormalized. */ \
365 if (_FP_FRAC_ZEROP_##wc(X)) \
366 { \
367 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, Y); \
368 _FP_FRAC_COPY_##wc(R, Y); \
369 goto add_done; \
370 } \
371 else \
372 { \
373 FP_SET_EXCEPTION(FP_EX_DENORM); \
374 ediff--; \
375 if (ediff == 0) \
376 { \
377 _FP_FRAC_ADD_##wc(R, Y, X); \
378 goto add3; \
379 } \
380 if (Y##_e == _FP_EXPMAX_##fs) \
381 { \
382 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, Y); \
383 _FP_FRAC_COPY_##wc(R, Y); \
384 goto add_done; \
385 } \
386 goto add2; \
387 } \
388 } \
389 else if (Y##_e == _FP_EXPMAX_##fs) \
390 { \
391 /* Y is NaN or Inf, X is normal. */ \
392 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, Y); \
393 _FP_FRAC_COPY_##wc(R, Y); \
394 goto add_done; \
395 } \
396 \
397 /* Insert implicit MSB of X. */ \
398 _FP_FRAC_HIGH_##fs(X) |= _FP_IMPLBIT_SH_##fs; \
399 \
400 add2: \
401 /* Shift the mantissa of X to the right EDIFF steps; \
402 remember to account later for the implicit MSB of Y. */ \
403 if (ediff <= _FP_WFRACBITS_##fs) \
404 _FP_FRAC_SRS_##wc(X, ediff, _FP_WFRACBITS_##fs); \
405 else if (!_FP_FRAC_ZEROP_##wc(X)) \
406 _FP_FRAC_SET_##wc(X, _FP_MINFRAC_##wc); \
407 _FP_FRAC_ADD_##wc(R, Y, X); \
408 } \
409 else \
410 { \
411 /* ediff == 0. */ \
412 if (!_FP_EXP_NORMAL(fs, wc, X)) \
413 { \
414 if (X##_e == 0) \
415 { \
416 /* X and Y are zero or denormalized. */ \
417 R##_e = 0; \
418 if (_FP_FRAC_ZEROP_##wc(X)) \
419 { \
420 if (!_FP_FRAC_ZEROP_##wc(Y)) \
421 FP_SET_EXCEPTION(FP_EX_DENORM); \
422 _FP_FRAC_COPY_##wc(R, Y); \
423 goto add_done; \
424 } \
425 else if (_FP_FRAC_ZEROP_##wc(Y)) \
426 { \
427 FP_SET_EXCEPTION(FP_EX_DENORM); \
428 _FP_FRAC_COPY_##wc(R, X); \
429 goto add_done; \
430 } \
431 else \
432 { \
433 FP_SET_EXCEPTION(FP_EX_DENORM); \
434 _FP_FRAC_ADD_##wc(R, X, Y); \
435 if (_FP_FRAC_HIGH_##fs(R) & _FP_IMPLBIT_SH_##fs) \
436 { \
437 /* Normalized result. */ \
438 _FP_FRAC_HIGH_##fs(R) \
439 &= ~(_FP_W_TYPE)_FP_IMPLBIT_SH_##fs; \
440 R##_e = 1; \
441 } \
442 goto add_done; \
443 } \
444 } \
445 else \
446 { \
447 /* X and Y are NaN or Inf. */ \
448 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, X); \
449 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, Y); \
450 R##_e = _FP_EXPMAX_##fs; \
451 if (_FP_FRAC_ZEROP_##wc(X)) \
452 _FP_FRAC_COPY_##wc(R, Y); \
453 else if (_FP_FRAC_ZEROP_##wc(Y)) \
454 _FP_FRAC_COPY_##wc(R, X); \
455 else \
456 _FP_CHOOSENAN_SEMIRAW(fs, wc, R, X, Y, OP); \
457 goto add_done; \
458 } \
459 } \
460 /* The exponents of X and Y, both normal, are equal. The \
461 implicit MSBs will always add to increase the \
462 exponent. */ \
463 _FP_FRAC_ADD_##wc(R, X, Y); \
464 R##_e = X##_e + 1; \
465 _FP_FRAC_SRS_##wc(R, 1, _FP_WFRACBITS_##fs); \
466 if (R##_e == _FP_EXPMAX_##fs) \
467 /* Overflow to infinity (depending on rounding mode). */ \
468 _FP_OVERFLOW_SEMIRAW(fs, wc, R); \
469 goto add_done; \
470 } \
471 add3: \
472 if (_FP_FRAC_HIGH_##fs(R) & _FP_IMPLBIT_SH_##fs) \
473 { \
474 /* Overflow. */ \
475 _FP_FRAC_HIGH_##fs(R) &= ~(_FP_W_TYPE)_FP_IMPLBIT_SH_##fs; \
476 R##_e++; \
477 _FP_FRAC_SRS_##wc(R, 1, _FP_WFRACBITS_##fs); \
478 if (R##_e == _FP_EXPMAX_##fs) \
479 /* Overflow to infinity (depending on rounding mode). */ \
480 _FP_OVERFLOW_SEMIRAW(fs, wc, R); \
481 } \
482 add_done: ; \
483 } \
484 else \
485 { \
486 /* Subtraction. */ \
487 int ediff = X##_e - Y##_e; \
488 if (ediff > 0) \
489 { \
490 R##_e = X##_e; \
491 R##_s = X##_s; \
492 if (Y##_e == 0) \
493 { \
494 /* Y is zero or denormalized. */ \
495 if (_FP_FRAC_ZEROP_##wc(Y)) \
496 { \
497 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, X); \
498 _FP_FRAC_COPY_##wc(R, X); \
499 goto sub_done; \
500 } \
501 else \
502 { \
503 FP_SET_EXCEPTION(FP_EX_DENORM); \
504 ediff--; \
505 if (ediff == 0) \
506 { \
507 _FP_FRAC_SUB_##wc(R, X, Y); \
508 goto sub3; \
509 } \
510 if (X##_e == _FP_EXPMAX_##fs) \
511 { \
512 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, X); \
513 _FP_FRAC_COPY_##wc(R, X); \
514 goto sub_done; \
515 } \
516 goto sub1; \
517 } \
518 } \
519 else if (X##_e == _FP_EXPMAX_##fs) \
520 { \
521 /* X is NaN or Inf, Y is normal. */ \
522 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, X); \
523 _FP_FRAC_COPY_##wc(R, X); \
524 goto sub_done; \
525 } \
526 \
527 /* Insert implicit MSB of Y. */ \
528 _FP_FRAC_HIGH_##fs(Y) |= _FP_IMPLBIT_SH_##fs; \
529 \
530 sub1: \
531 /* Shift the mantissa of Y to the right EDIFF steps; \
532 remember to account later for the implicit MSB of X. */ \
533 if (ediff <= _FP_WFRACBITS_##fs) \
534 _FP_FRAC_SRS_##wc(Y, ediff, _FP_WFRACBITS_##fs); \
535 else if (!_FP_FRAC_ZEROP_##wc(Y)) \
536 _FP_FRAC_SET_##wc(Y, _FP_MINFRAC_##wc); \
537 _FP_FRAC_SUB_##wc(R, X, Y); \
538 } \
539 else if (ediff < 0) \
540 { \
541 ediff = -ediff; \
542 R##_e = Y##_e; \
543 R##_s = Y##_s; \
544 if (X##_e == 0) \
545 { \
546 /* X is zero or denormalized. */ \
547 if (_FP_FRAC_ZEROP_##wc(X)) \
548 { \
549 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, Y); \
550 _FP_FRAC_COPY_##wc(R, Y); \
551 goto sub_done; \
552 } \
553 else \
554 { \
555 FP_SET_EXCEPTION(FP_EX_DENORM); \
556 ediff--; \
557 if (ediff == 0) \
558 { \
559 _FP_FRAC_SUB_##wc(R, Y, X); \
560 goto sub3; \
561 } \
562 if (Y##_e == _FP_EXPMAX_##fs) \
563 { \
564 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, Y); \
565 _FP_FRAC_COPY_##wc(R, Y); \
566 goto sub_done; \
567 } \
568 goto sub2; \
569 } \
570 } \
571 else if (Y##_e == _FP_EXPMAX_##fs) \
572 { \
573 /* Y is NaN or Inf, X is normal. */ \
574 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, Y); \
575 _FP_FRAC_COPY_##wc(R, Y); \
576 goto sub_done; \
577 } \
578 \
579 /* Insert implicit MSB of X. */ \
580 _FP_FRAC_HIGH_##fs(X) |= _FP_IMPLBIT_SH_##fs; \
581 \
582 sub2: \
583 /* Shift the mantissa of X to the right EDIFF steps; \
584 remember to account later for the implicit MSB of Y. */ \
585 if (ediff <= _FP_WFRACBITS_##fs) \
586 _FP_FRAC_SRS_##wc(X, ediff, _FP_WFRACBITS_##fs); \
587 else if (!_FP_FRAC_ZEROP_##wc(X)) \
588 _FP_FRAC_SET_##wc(X, _FP_MINFRAC_##wc); \
589 _FP_FRAC_SUB_##wc(R, Y, X); \
590 } \
591 else \
592 { \
593 /* ediff == 0. */ \
594 if (!_FP_EXP_NORMAL(fs, wc, X)) \
595 { \
596 if (X##_e == 0) \
597 { \
598 /* X and Y are zero or denormalized. */ \
599 R##_e = 0; \
600 if (_FP_FRAC_ZEROP_##wc(X)) \
601 { \
602 _FP_FRAC_COPY_##wc(R, Y); \
603 if (_FP_FRAC_ZEROP_##wc(Y)) \
604 R##_s = (FP_ROUNDMODE == FP_RND_MINF); \
605 else \
606 { \
607 FP_SET_EXCEPTION(FP_EX_DENORM); \
608 R##_s = Y##_s; \
609 } \
610 goto sub_done; \
611 } \
612 else if (_FP_FRAC_ZEROP_##wc(Y)) \
613 { \
614 FP_SET_EXCEPTION(FP_EX_DENORM); \
615 _FP_FRAC_COPY_##wc(R, X); \
616 R##_s = X##_s; \
617 goto sub_done; \
618 } \
619 else \
620 { \
621 FP_SET_EXCEPTION(FP_EX_DENORM); \
622 _FP_FRAC_SUB_##wc(R, X, Y); \
623 R##_s = X##_s; \
624 if (_FP_FRAC_HIGH_##fs(R) & _FP_IMPLBIT_SH_##fs) \
625 { \
626 /* |X| < |Y|, negate result. */ \
627 _FP_FRAC_SUB_##wc(R, Y, X); \
628 R##_s = Y##_s; \
629 } \
630 else if (_FP_FRAC_ZEROP_##wc(R)) \
631 R##_s = (FP_ROUNDMODE == FP_RND_MINF); \
632 goto sub_done; \
633 } \
634 } \
635 else \
636 { \
637 /* X and Y are NaN or Inf, of opposite signs. */ \
638 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, X); \
639 _FP_CHECK_SIGNAN_SEMIRAW(fs, wc, Y); \
640 R##_e = _FP_EXPMAX_##fs; \
641 if (_FP_FRAC_ZEROP_##wc(X)) \
642 { \
643 if (_FP_FRAC_ZEROP_##wc(Y)) \
644 { \
645 /* Inf - Inf. */ \
646 R##_s = _FP_NANSIGN_##fs; \
647 _FP_FRAC_SET_##wc(R, _FP_NANFRAC_##fs); \
648 _FP_FRAC_SLL_##wc(R, _FP_WORKBITS); \
649 FP_SET_EXCEPTION(FP_EX_INVALID); \
650 } \
651 else \
652 { \
653 /* Inf - NaN. */ \
654 R##_s = Y##_s; \
655 _FP_FRAC_COPY_##wc(R, Y); \
656 } \
657 } \
658 else \
659 { \
660 if (_FP_FRAC_ZEROP_##wc(Y)) \
661 { \
662 /* NaN - Inf. */ \
663 R##_s = X##_s; \
664 _FP_FRAC_COPY_##wc(R, X); \
665 } \
666 else \
667 { \
668 /* NaN - NaN. */ \
669 _FP_CHOOSENAN_SEMIRAW(fs, wc, R, X, Y, OP); \
670 } \
671 } \
672 goto sub_done; \
673 } \
674 } \
675 /* The exponents of X and Y, both normal, are equal. The \
676 implicit MSBs cancel. */ \
677 R##_e = X##_e; \
678 _FP_FRAC_SUB_##wc(R, X, Y); \
679 R##_s = X##_s; \
680 if (_FP_FRAC_HIGH_##fs(R) & _FP_IMPLBIT_SH_##fs) \
681 { \
682 /* |X| < |Y|, negate result. */ \
683 _FP_FRAC_SUB_##wc(R, Y, X); \
684 R##_s = Y##_s; \
685 } \
686 else if (_FP_FRAC_ZEROP_##wc(R)) \
687 { \
688 R##_e = 0; \
689 R##_s = (FP_ROUNDMODE == FP_RND_MINF); \
690 goto sub_done; \
691 } \
692 goto norm; \
693 } \
694 sub3: \
695 if (_FP_FRAC_HIGH_##fs(R) & _FP_IMPLBIT_SH_##fs) \
696 { \
697 int diff; \
698 /* Carry into most significant bit of larger one of X and Y, \
699 canceling it; renormalize. */ \
700 _FP_FRAC_HIGH_##fs(R) &= _FP_IMPLBIT_SH_##fs - 1; \
701 norm: \
702 _FP_FRAC_CLZ_##wc(diff, R); \
703 diff -= _FP_WFRACXBITS_##fs; \
704 _FP_FRAC_SLL_##wc(R, diff); \
705 if (R##_e <= diff) \
706 { \
707 /* R is denormalized. */ \
708 diff = diff - R##_e + 1; \
709 _FP_FRAC_SRS_##wc(R, diff, _FP_WFRACBITS_##fs); \
710 R##_e = 0; \
711 } \
712 else \
713 { \
714 R##_e -= diff; \
715 _FP_FRAC_HIGH_##fs(R) &= ~(_FP_W_TYPE)_FP_IMPLBIT_SH_##fs; \
716 } \
717 } \
718 sub_done: ; \
719 } \
720 } while (0)
721
722 #define _FP_ADD(fs, wc, R, X, Y) _FP_ADD_INTERNAL(fs, wc, R, X, Y, '+')
723 #define _FP_SUB(fs, wc, R, X, Y) \
724 do { \
725 if (!(Y##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc(Y))) Y##_s ^= 1; \
726 _FP_ADD_INTERNAL(fs, wc, R, X, Y, '-'); \
727 } while (0)
728
729
730 /*
731 * Main negation routine. FIXME -- when we care about setting exception
732 * bits reliably, this will not do. We should examine all of the fp classes.
733 */
734
735 #define _FP_NEG(fs, wc, R, X) \
736 do { \
737 _FP_FRAC_COPY_##wc(R, X); \
738 R##_c = X##_c; \
739 R##_e = X##_e; \
740 R##_s = 1 ^ X##_s; \
741 } while (0)
742
743
744 /*
745 * Main multiplication routine. The input values should be cooked.
746 */
747
748 #define _FP_MUL(fs, wc, R, X, Y) \
749 do { \
750 R##_s = X##_s ^ Y##_s; \
751 R##_e = X##_e + Y##_e + 1; \
752 switch (_FP_CLS_COMBINE(X##_c, Y##_c)) \
753 { \
754 case _FP_CLS_COMBINE(FP_CLS_NORMAL,FP_CLS_NORMAL): \
755 R##_c = FP_CLS_NORMAL; \
756 \
757 _FP_MUL_MEAT_##fs(R,X,Y); \
758 \
759 if (_FP_FRAC_OVERP_##wc(fs, R)) \
760 _FP_FRAC_SRS_##wc(R, 1, _FP_WFRACBITS_##fs); \
761 else \
762 R##_e--; \
763 break; \
764 \
765 case _FP_CLS_COMBINE(FP_CLS_NAN,FP_CLS_NAN): \
766 _FP_CHOOSENAN(fs, wc, R, X, Y, '*'); \
767 break; \
768 \
769 case _FP_CLS_COMBINE(FP_CLS_NAN,FP_CLS_NORMAL): \
770 case _FP_CLS_COMBINE(FP_CLS_NAN,FP_CLS_INF): \
771 case _FP_CLS_COMBINE(FP_CLS_NAN,FP_CLS_ZERO): \
772 R##_s = X##_s; \
773 \
774 case _FP_CLS_COMBINE(FP_CLS_INF,FP_CLS_INF): \
775 case _FP_CLS_COMBINE(FP_CLS_INF,FP_CLS_NORMAL): \
776 case _FP_CLS_COMBINE(FP_CLS_ZERO,FP_CLS_NORMAL): \
777 case _FP_CLS_COMBINE(FP_CLS_ZERO,FP_CLS_ZERO): \
778 _FP_FRAC_COPY_##wc(R, X); \
779 R##_c = X##_c; \
780 break; \
781 \
782 case _FP_CLS_COMBINE(FP_CLS_NORMAL,FP_CLS_NAN): \
783 case _FP_CLS_COMBINE(FP_CLS_INF,FP_CLS_NAN): \
784 case _FP_CLS_COMBINE(FP_CLS_ZERO,FP_CLS_NAN): \
785 R##_s = Y##_s; \
786 \
787 case _FP_CLS_COMBINE(FP_CLS_NORMAL,FP_CLS_INF): \
788 case _FP_CLS_COMBINE(FP_CLS_NORMAL,FP_CLS_ZERO): \
789 _FP_FRAC_COPY_##wc(R, Y); \
790 R##_c = Y##_c; \
791 break; \
792 \
793 case _FP_CLS_COMBINE(FP_CLS_INF,FP_CLS_ZERO): \
794 case _FP_CLS_COMBINE(FP_CLS_ZERO,FP_CLS_INF): \
795 R##_s = _FP_NANSIGN_##fs; \
796 R##_c = FP_CLS_NAN; \
797 _FP_FRAC_SET_##wc(R, _FP_NANFRAC_##fs); \
798 FP_SET_EXCEPTION(FP_EX_INVALID); \
799 break; \
800 \
801 default: \
802 abort(); \
803 } \
804 } while (0)
805
806
807 /*
808 * Main division routine. The input values should be cooked.
809 */
810
811 #define _FP_DIV(fs, wc, R, X, Y) \
812 do { \
813 R##_s = X##_s ^ Y##_s; \
814 R##_e = X##_e - Y##_e; \
815 switch (_FP_CLS_COMBINE(X##_c, Y##_c)) \
816 { \
817 case _FP_CLS_COMBINE(FP_CLS_NORMAL,FP_CLS_NORMAL): \
818 R##_c = FP_CLS_NORMAL; \
819 \
820 _FP_DIV_MEAT_##fs(R,X,Y); \
821 break; \
822 \
823 case _FP_CLS_COMBINE(FP_CLS_NAN,FP_CLS_NAN): \
824 _FP_CHOOSENAN(fs, wc, R, X, Y, '/'); \
825 break; \
826 \
827 case _FP_CLS_COMBINE(FP_CLS_NAN,FP_CLS_NORMAL): \
828 case _FP_CLS_COMBINE(FP_CLS_NAN,FP_CLS_INF): \
829 case _FP_CLS_COMBINE(FP_CLS_NAN,FP_CLS_ZERO): \
830 R##_s = X##_s; \
831 _FP_FRAC_COPY_##wc(R, X); \
832 R##_c = X##_c; \
833 break; \
834 \
835 case _FP_CLS_COMBINE(FP_CLS_NORMAL,FP_CLS_NAN): \
836 case _FP_CLS_COMBINE(FP_CLS_INF,FP_CLS_NAN): \
837 case _FP_CLS_COMBINE(FP_CLS_ZERO,FP_CLS_NAN): \
838 R##_s = Y##_s; \
839 _FP_FRAC_COPY_##wc(R, Y); \
840 R##_c = Y##_c; \
841 break; \
842 \
843 case _FP_CLS_COMBINE(FP_CLS_NORMAL,FP_CLS_INF): \
844 case _FP_CLS_COMBINE(FP_CLS_ZERO,FP_CLS_INF): \
845 case _FP_CLS_COMBINE(FP_CLS_ZERO,FP_CLS_NORMAL): \
846 R##_c = FP_CLS_ZERO; \
847 break; \
848 \
849 case _FP_CLS_COMBINE(FP_CLS_NORMAL,FP_CLS_ZERO): \
850 FP_SET_EXCEPTION(FP_EX_DIVZERO); \
851 case _FP_CLS_COMBINE(FP_CLS_INF,FP_CLS_ZERO): \
852 case _FP_CLS_COMBINE(FP_CLS_INF,FP_CLS_NORMAL): \
853 R##_c = FP_CLS_INF; \
854 break; \
855 \
856 case _FP_CLS_COMBINE(FP_CLS_INF,FP_CLS_INF): \
857 case _FP_CLS_COMBINE(FP_CLS_ZERO,FP_CLS_ZERO): \
858 R##_s = _FP_NANSIGN_##fs; \
859 R##_c = FP_CLS_NAN; \
860 _FP_FRAC_SET_##wc(R, _FP_NANFRAC_##fs); \
861 FP_SET_EXCEPTION(FP_EX_INVALID); \
862 break; \
863 \
864 default: \
865 abort(); \
866 } \
867 } while (0)
868
869
870 /*
871 * Main differential comparison routine. The inputs should be raw not
872 * cooked. The return is -1,0,1 for normal values, 2 otherwise.
873 */
874
875 #define _FP_CMP(fs, wc, ret, X, Y, un) \
876 do { \
877 /* NANs are unordered */ \
878 if ((X##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc(X)) \
879 || (Y##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc(Y))) \
880 { \
881 ret = un; \
882 } \
883 else \
884 { \
885 int __is_zero_x; \
886 int __is_zero_y; \
887 \
888 __is_zero_x = (!X##_e && _FP_FRAC_ZEROP_##wc(X)) ? 1 : 0; \
889 __is_zero_y = (!Y##_e && _FP_FRAC_ZEROP_##wc(Y)) ? 1 : 0; \
890 \
891 if (__is_zero_x && __is_zero_y) \
892 ret = 0; \
893 else if (__is_zero_x) \
894 ret = Y##_s ? 1 : -1; \
895 else if (__is_zero_y) \
896 ret = X##_s ? -1 : 1; \
897 else if (X##_s != Y##_s) \
898 ret = X##_s ? -1 : 1; \
899 else if (X##_e > Y##_e) \
900 ret = X##_s ? -1 : 1; \
901 else if (X##_e < Y##_e) \
902 ret = X##_s ? 1 : -1; \
903 else if (_FP_FRAC_GT_##wc(X, Y)) \
904 ret = X##_s ? -1 : 1; \
905 else if (_FP_FRAC_GT_##wc(Y, X)) \
906 ret = X##_s ? 1 : -1; \
907 else \
908 ret = 0; \
909 } \
910 } while (0)
911
912
913 /* Simplification for strict equality. */
914
915 #define _FP_CMP_EQ(fs, wc, ret, X, Y) \
916 do { \
917 /* NANs are unordered */ \
918 if ((X##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc(X)) \
919 || (Y##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc(Y))) \
920 { \
921 ret = 1; \
922 } \
923 else \
924 { \
925 ret = !(X##_e == Y##_e \
926 && _FP_FRAC_EQ_##wc(X, Y) \
927 && (X##_s == Y##_s || (!X##_e && _FP_FRAC_ZEROP_##wc(X)))); \
928 } \
929 } while (0)
930
931 /* Version to test unordered. */
932
933 #define _FP_CMP_UNORD(fs, wc, ret, X, Y) \
934 do { \
935 ret = ((X##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc(X)) \
936 || (Y##_e == _FP_EXPMAX_##fs && !_FP_FRAC_ZEROP_##wc(Y))); \
937 } while (0)
938
939 /*
940 * Main square root routine. The input value should be cooked.
941 */
942
943 #define _FP_SQRT(fs, wc, R, X) \
944 do { \
945 _FP_FRAC_DECL_##wc(T); _FP_FRAC_DECL_##wc(S); \
946 _FP_W_TYPE q; \
947 switch (X##_c) \
948 { \
949 case FP_CLS_NAN: \
950 _FP_FRAC_COPY_##wc(R, X); \
951 R##_s = X##_s; \
952 R##_c = FP_CLS_NAN; \
953 break; \
954 case FP_CLS_INF: \
955 if (X##_s) \
956 { \
957 R##_s = _FP_NANSIGN_##fs; \
958 R##_c = FP_CLS_NAN; /* NAN */ \
959 _FP_FRAC_SET_##wc(R, _FP_NANFRAC_##fs); \
960 FP_SET_EXCEPTION(FP_EX_INVALID); \
961 } \
962 else \
963 { \
964 R##_s = 0; \
965 R##_c = FP_CLS_INF; /* sqrt(+inf) = +inf */ \
966 } \
967 break; \
968 case FP_CLS_ZERO: \
969 R##_s = X##_s; \
970 R##_c = FP_CLS_ZERO; /* sqrt(+-0) = +-0 */ \
971 break; \
972 case FP_CLS_NORMAL: \
973 R##_s = 0; \
974 if (X##_s) \
975 { \
976 R##_c = FP_CLS_NAN; /* NAN */ \
977 R##_s = _FP_NANSIGN_##fs; \
978 _FP_FRAC_SET_##wc(R, _FP_NANFRAC_##fs); \
979 FP_SET_EXCEPTION(FP_EX_INVALID); \
980 break; \
981 } \
982 R##_c = FP_CLS_NORMAL; \
983 if (X##_e & 1) \
984 _FP_FRAC_SLL_##wc(X, 1); \
985 R##_e = X##_e >> 1; \
986 _FP_FRAC_SET_##wc(S, _FP_ZEROFRAC_##wc); \
987 _FP_FRAC_SET_##wc(R, _FP_ZEROFRAC_##wc); \
988 q = _FP_OVERFLOW_##fs >> 1; \
989 _FP_SQRT_MEAT_##wc(R, S, T, X, q); \
990 } \
991 } while (0)
992
993 /*
994 * Convert from FP to integer. Input is raw.
995 */
996
997 /* RSIGNED can have following values:
998 * 0: the number is required to be 0..(2^rsize)-1, if not, NV is set plus
999 * the result is either 0 or (2^rsize)-1 depending on the sign in such
1000 * case.
1001 * 1: the number is required to be -(2^(rsize-1))..(2^(rsize-1))-1, if not,
1002 * NV is set plus the result is either -(2^(rsize-1)) or (2^(rsize-1))-1
1003 * depending on the sign in such case.
1004 * -1: the number is required to be -(2^(rsize-1))..(2^rsize)-1, if not, NV is
1005 * set plus the result is either -(2^(rsize-1)) or (2^(rsize-1))-1
1006 * depending on the sign in such case.
1007 */
1008 #define _FP_TO_INT(fs, wc, r, X, rsize, rsigned) \
1009 do { \
1010 if (X##_e < _FP_EXPBIAS_##fs) \
1011 { \
1012 r = 0; \
1013 if (X##_e == 0) \
1014 { \
1015 if (!_FP_FRAC_ZEROP_##wc(X)) \
1016 { \
1017 FP_SET_EXCEPTION(FP_EX_INEXACT); \
1018 FP_SET_EXCEPTION(FP_EX_DENORM); \
1019 } \
1020 } \
1021 else \
1022 FP_SET_EXCEPTION(FP_EX_INEXACT); \
1023 } \
1024 else if (X##_e >= _FP_EXPBIAS_##fs + rsize - (rsigned > 0 || X##_s) \
1025 || (!rsigned && X##_s)) \
1026 { \
1027 /* Overflow or converting to the most negative integer. */ \
1028 if (rsigned) \
1029 { \
1030 r = 1; \
1031 r <<= rsize - 1; \
1032 r -= 1 - X##_s; \
1033 } else { \
1034 r = 0; \
1035 if (X##_s) \
1036 r = ~r; \
1037 } \
1038 \
1039 if (rsigned && X##_s && X##_e == _FP_EXPBIAS_##fs + rsize - 1) \
1040 { \
1041 /* Possibly converting to most negative integer; check the \
1042 mantissa. */ \
1043 int inexact = 0; \
1044 (void)((_FP_FRACBITS_##fs > rsize) \
1045 ? ({ _FP_FRAC_SRST_##wc(X, inexact, \
1046 _FP_FRACBITS_##fs - rsize, \
1047 _FP_FRACBITS_##fs); 0; }) \
1048 : 0); \
1049 if (!_FP_FRAC_ZEROP_##wc(X)) \
1050 FP_SET_EXCEPTION(FP_EX_INVALID); \
1051 else if (inexact) \
1052 FP_SET_EXCEPTION(FP_EX_INEXACT); \
1053 } \
1054 else \
1055 FP_SET_EXCEPTION(FP_EX_INVALID); \
1056 } \
1057 else \
1058 { \
1059 _FP_FRAC_HIGH_RAW_##fs(X) |= _FP_IMPLBIT_##fs; \
1060 if (X##_e >= _FP_EXPBIAS_##fs + _FP_FRACBITS_##fs - 1) \
1061 { \
1062 _FP_FRAC_ASSEMBLE_##wc(r, X, rsize); \
1063 r <<= X##_e - _FP_EXPBIAS_##fs - _FP_FRACBITS_##fs + 1; \
1064 } \
1065 else \
1066 { \
1067 int inexact; \
1068 _FP_FRAC_SRST_##wc(X, inexact, \
1069 (_FP_FRACBITS_##fs + _FP_EXPBIAS_##fs - 1 \
1070 - X##_e), \
1071 _FP_FRACBITS_##fs); \
1072 if (inexact) \
1073 FP_SET_EXCEPTION(FP_EX_INEXACT); \
1074 _FP_FRAC_ASSEMBLE_##wc(r, X, rsize); \
1075 } \
1076 if (rsigned && X##_s) \
1077 r = -r; \
1078 } \
1079 } while (0)
1080
1081 /* Convert integer to fp. Output is raw. RTYPE is unsigned even if
1082 input is signed. */
1083 #define _FP_FROM_INT(fs, wc, X, r, rsize, rtype) \
1084 do { \
1085 if (r) \
1086 { \
1087 rtype ur_; \
1088 \
1089 if ((X##_s = (r < 0))) \
1090 r = -(rtype)r; \
1091 \
1092 ur_ = (rtype) r; \
1093 (void)((rsize <= _FP_W_TYPE_SIZE) \
1094 ? ({ \
1095 int lz_; \
1096 __FP_CLZ(lz_, (_FP_W_TYPE)ur_); \
1097 X##_e = _FP_EXPBIAS_##fs + _FP_W_TYPE_SIZE - 1 - lz_; \
1098 }) \
1099 : ((rsize <= 2 * _FP_W_TYPE_SIZE) \
1100 ? ({ \
1101 int lz_; \
1102 __FP_CLZ_2(lz_, (_FP_W_TYPE)(ur_ >> _FP_W_TYPE_SIZE), \
1103 (_FP_W_TYPE)ur_); \
1104 X##_e = (_FP_EXPBIAS_##fs + 2 * _FP_W_TYPE_SIZE - 1 \
1105 - lz_); \
1106 }) \
1107 : (abort(), 0))); \
1108 \
1109 if (rsize - 1 + _FP_EXPBIAS_##fs >= _FP_EXPMAX_##fs \
1110 && X##_e >= _FP_EXPMAX_##fs) \
1111 { \
1112 /* Exponent too big; overflow to infinity. (May also \
1113 happen after rounding below.) */ \
1114 _FP_OVERFLOW_SEMIRAW(fs, wc, X); \
1115 goto pack_semiraw; \
1116 } \
1117 \
1118 if (rsize <= _FP_FRACBITS_##fs \
1119 || X##_e < _FP_EXPBIAS_##fs + _FP_FRACBITS_##fs) \
1120 { \
1121 /* Exactly representable; shift left. */ \
1122 _FP_FRAC_DISASSEMBLE_##wc(X, ur_, rsize); \
1123 _FP_FRAC_SLL_##wc(X, (_FP_EXPBIAS_##fs \
1124 + _FP_FRACBITS_##fs - 1 - X##_e)); \
1125 } \
1126 else \
1127 { \
1128 /* More bits in integer than in floating type; need to \
1129 round. */ \
1130 if (_FP_EXPBIAS_##fs + _FP_WFRACBITS_##fs - 1 < X##_e) \
1131 ur_ = ((ur_ >> (X##_e - _FP_EXPBIAS_##fs \
1132 - _FP_WFRACBITS_##fs + 1)) \
1133 | ((ur_ << (rsize - (X##_e - _FP_EXPBIAS_##fs \
1134 - _FP_WFRACBITS_##fs + 1))) \
1135 != 0)); \
1136 _FP_FRAC_DISASSEMBLE_##wc(X, ur_, rsize); \
1137 if ((_FP_EXPBIAS_##fs + _FP_WFRACBITS_##fs - 1 - X##_e) > 0) \
1138 _FP_FRAC_SLL_##wc(X, (_FP_EXPBIAS_##fs \
1139 + _FP_WFRACBITS_##fs - 1 - X##_e)); \
1140 _FP_FRAC_HIGH_##fs(X) &= ~(_FP_W_TYPE)_FP_IMPLBIT_SH_##fs; \
1141 pack_semiraw: \
1142 _FP_PACK_SEMIRAW(fs, wc, X); \
1143 } \
1144 } \
1145 else \
1146 { \
1147 X##_s = 0; \
1148 X##_e = 0; \
1149 _FP_FRAC_SET_##wc(X, _FP_ZEROFRAC_##wc); \
1150 } \
1151 } while (0)
1152
1153
1154 /* Extend from a narrower floating-point format to a wider one. Input
1155 and output are raw. */
1156 #define FP_EXTEND(dfs,sfs,dwc,swc,D,S) \
1157 do { \
1158 if (_FP_FRACBITS_##dfs < _FP_FRACBITS_##sfs \
1159 || (_FP_EXPMAX_##dfs - _FP_EXPBIAS_##dfs \
1160 < _FP_EXPMAX_##sfs - _FP_EXPBIAS_##sfs) \
1161 || (_FP_EXPBIAS_##dfs < _FP_EXPBIAS_##sfs + _FP_FRACBITS_##sfs - 1 \
1162 && _FP_EXPBIAS_##dfs != _FP_EXPBIAS_##sfs)) \
1163 abort(); \
1164 D##_s = S##_s; \
1165 _FP_FRAC_COPY_##dwc##_##swc(D, S); \
1166 if (_FP_EXP_NORMAL(sfs, swc, S)) \
1167 { \
1168 D##_e = S##_e + _FP_EXPBIAS_##dfs - _FP_EXPBIAS_##sfs; \
1169 _FP_FRAC_SLL_##dwc(D, (_FP_FRACBITS_##dfs - _FP_FRACBITS_##sfs)); \
1170 } \
1171 else \
1172 { \
1173 if (S##_e == 0) \
1174 { \
1175 if (_FP_FRAC_ZEROP_##swc(S)) \
1176 D##_e = 0; \
1177 else if (_FP_EXPBIAS_##dfs \
1178 < _FP_EXPBIAS_##sfs + _FP_FRACBITS_##sfs - 1) \
1179 { \
1180 FP_SET_EXCEPTION(FP_EX_DENORM); \
1181 _FP_FRAC_SLL_##dwc(D, (_FP_FRACBITS_##dfs \
1182 - _FP_FRACBITS_##sfs)); \
1183 D##_e = 0; \
1184 } \
1185 else \
1186 { \
1187 int _lz; \
1188 FP_SET_EXCEPTION(FP_EX_DENORM); \
1189 _FP_FRAC_CLZ_##swc(_lz, S); \
1190 _FP_FRAC_SLL_##dwc(D, \
1191 _lz + _FP_FRACBITS_##dfs \
1192 - _FP_FRACTBITS_##sfs); \
1193 D##_e = (_FP_EXPBIAS_##dfs - _FP_EXPBIAS_##sfs + 1 \
1194 + _FP_FRACXBITS_##sfs - _lz); \
1195 } \
1196 } \
1197 else \
1198 { \
1199 D##_e = _FP_EXPMAX_##dfs; \
1200 if (!_FP_FRAC_ZEROP_##swc(S)) \
1201 { \
1202 if (!(_FP_FRAC_HIGH_RAW_##sfs(S) & _FP_QNANBIT_##sfs)) \
1203 FP_SET_EXCEPTION(FP_EX_INVALID); \
1204 _FP_FRAC_SLL_##dwc(D, (_FP_FRACBITS_##dfs \
1205 - _FP_FRACBITS_##sfs)); \
1206 } \
1207 } \
1208 } \
1209 } while (0)
1210
1211 /* Truncate from a wider floating-point format to a narrower one.
1212 Input and output are semi-raw. */
1213 #define FP_TRUNC(dfs,sfs,dwc,swc,D,S) \
1214 do { \
1215 if (_FP_FRACBITS_##sfs < _FP_FRACBITS_##dfs \
1216 || (_FP_EXPBIAS_##sfs < _FP_EXPBIAS_##dfs + _FP_FRACBITS_##dfs - 1 \
1217 && _FP_EXPBIAS_##sfs != _FP_EXPBIAS_##dfs)) \
1218 abort(); \
1219 D##_s = S##_s; \
1220 if (_FP_EXP_NORMAL(sfs, swc, S)) \
1221 { \
1222 D##_e = S##_e + _FP_EXPBIAS_##dfs - _FP_EXPBIAS_##sfs; \
1223 if (D##_e >= _FP_EXPMAX_##dfs) \
1224 _FP_OVERFLOW_SEMIRAW(dfs, dwc, D); \
1225 else \
1226 { \
1227 if (D##_e <= 0) \
1228 { \
1229 if (D##_e < 1 - _FP_FRACBITS_##dfs) \
1230 { \
1231 _FP_FRAC_SET_##swc(S, _FP_ZEROFRAC_##swc); \
1232 _FP_FRAC_LOW_##swc(S) |= 1; \
1233 } \
1234 else \
1235 { \
1236 _FP_FRAC_HIGH_##sfs(S) |= _FP_IMPLBIT_SH_##sfs; \
1237 _FP_FRAC_SRS_##swc(S, (_FP_WFRACBITS_##sfs \
1238 - _FP_WFRACBITS_##dfs + 1 - D##_e), \
1239 _FP_WFRACBITS_##sfs); \
1240 } \
1241 D##_e = 0; \
1242 } \
1243 else \
1244 _FP_FRAC_SRS_##swc(S, (_FP_WFRACBITS_##sfs \
1245 - _FP_WFRACBITS_##dfs), \
1246 _FP_WFRACBITS_##sfs); \
1247 _FP_FRAC_COPY_##dwc##_##swc(D, S); \
1248 } \
1249 } \
1250 else \
1251 { \
1252 if (S##_e == 0) \
1253 { \
1254 D##_e = 0; \
1255 if (_FP_FRAC_ZEROP_##swc(S)) \
1256 _FP_FRAC_SET_##dwc(D, _FP_ZEROFRAC_##dwc); \
1257 else \
1258 { \
1259 FP_SET_EXCEPTION(FP_EX_DENORM); \
1260 if (_FP_EXPBIAS_##sfs \
1261 < _FP_EXPBIAS_##dfs + _FP_FRACBITS_##dfs - 1) \
1262 { \
1263 _FP_FRAC_SRS_##swc(S, (_FP_WFRACBITS_##sfs \
1264 - _FP_WFRACBITS_##dfs), \
1265 _FP_WFRACBITS_##sfs); \
1266 _FP_FRAC_COPY_##dwc##_##swc(D, S); \
1267 } \
1268 else \
1269 { \
1270 _FP_FRAC_SET_##dwc(D, _FP_ZEROFRAC_##dwc); \
1271 _FP_FRAC_LOW_##dwc(D) |= 1; \
1272 } \
1273 } \
1274 } \
1275 else \
1276 { \
1277 D##_e = _FP_EXPMAX_##dfs; \
1278 if (_FP_FRAC_ZEROP_##swc(S)) \
1279 _FP_FRAC_SET_##dwc(D, _FP_ZEROFRAC_##dwc); \
1280 else \
1281 { \
1282 _FP_CHECK_SIGNAN_SEMIRAW(sfs, swc, S); \
1283 _FP_FRAC_SRL_##swc(S, (_FP_WFRACBITS_##sfs \
1284 - _FP_WFRACBITS_##dfs)); \
1285 _FP_FRAC_COPY_##dwc##_##swc(D, S); \
1286 /* Semi-raw NaN must have all workbits cleared. */ \
1287 _FP_FRAC_LOW_##dwc(D) \
1288 &= ~(_FP_W_TYPE) ((1 << _FP_WORKBITS) - 1); \
1289 _FP_FRAC_HIGH_##dfs(D) |= _FP_QNANBIT_SH_##dfs; \
1290 } \
1291 } \
1292 } \
1293 } while (0)
1294
1295 /*
1296 * Helper primitives.
1297 */
1298
1299 /* Count leading zeros in a word. */
1300
1301 #ifndef __FP_CLZ
1302 /* GCC 3.4 and later provide the builtins for us. */
1303 #define __FP_CLZ(r, x) \
1304 do { \
1305 if (sizeof (_FP_W_TYPE) == sizeof (unsigned int)) \
1306 r = __builtin_clz (x); \
1307 else if (sizeof (_FP_W_TYPE) == sizeof (unsigned long)) \
1308 r = __builtin_clzl (x); \
1309 else if (sizeof (_FP_W_TYPE) == sizeof (unsigned long long)) \
1310 r = __builtin_clzll (x); \
1311 else \
1312 abort (); \
1313 } while (0)
1314 #endif /* ndef __FP_CLZ */
1315
1316 #define _FP_DIV_HELP_imm(q, r, n, d) \
1317 do { \
1318 q = n / d, r = n % d; \
1319 } while (0)
1320
1321
1322 /* A restoring bit-by-bit division primitive. */
1323
1324 #define _FP_DIV_MEAT_N_loop(fs, wc, R, X, Y) \
1325 do { \
1326 int count = _FP_WFRACBITS_##fs; \
1327 _FP_FRAC_DECL_##wc (u); \
1328 _FP_FRAC_DECL_##wc (v); \
1329 _FP_FRAC_COPY_##wc (u, X); \
1330 _FP_FRAC_COPY_##wc (v, Y); \
1331 _FP_FRAC_SET_##wc (R, _FP_ZEROFRAC_##wc); \
1332 /* Normalize U and V. */ \
1333 _FP_FRAC_SLL_##wc (u, _FP_WFRACXBITS_##fs); \
1334 _FP_FRAC_SLL_##wc (v, _FP_WFRACXBITS_##fs); \
1335 /* First round. Since the operands are normalized, either the \
1336 first or second bit will be set in the fraction. Produce a \
1337 normalized result by checking which and adjusting the loop \
1338 count and exponent accordingly. */ \
1339 if (_FP_FRAC_GE_1 (u, v)) \
1340 { \
1341 _FP_FRAC_SUB_##wc (u, u, v); \
1342 _FP_FRAC_LOW_##wc (R) |= 1; \
1343 count--; \
1344 } \
1345 else \
1346 R##_e--; \
1347 /* Subsequent rounds. */ \
1348 do { \
1349 int msb = (_FP_WS_TYPE) _FP_FRAC_HIGH_##wc (u) < 0; \
1350 _FP_FRAC_SLL_##wc (u, 1); \
1351 _FP_FRAC_SLL_##wc (R, 1); \
1352 if (msb || _FP_FRAC_GE_1 (u, v)) \
1353 { \
1354 _FP_FRAC_SUB_##wc (u, u, v); \
1355 _FP_FRAC_LOW_##wc (R) |= 1; \
1356 } \
1357 } while (--count > 0); \
1358 /* If there's anything left in U, the result is inexact. */ \
1359 _FP_FRAC_LOW_##wc (R) |= !_FP_FRAC_ZEROP_##wc (u); \
1360 } while (0)
1361
1362 #define _FP_DIV_MEAT_1_loop(fs, R, X, Y) _FP_DIV_MEAT_N_loop (fs, 1, R, X, Y)
1363 #define _FP_DIV_MEAT_2_loop(fs, R, X, Y) _FP_DIV_MEAT_N_loop (fs, 2, R, X, Y)
1364 #define _FP_DIV_MEAT_4_loop(fs, R, X, Y) _FP_DIV_MEAT_N_loop (fs, 4, R, X, Y)
sources.redhat.com git mirror of glibc CVS