1 /* Copyright (C) 2007-2024 Free Software Foundation, Inc.
3 This file is part of GCC.
5 GCC is free software; you can redistribute it and/or modify it under
6 the terms of the GNU General Public License as published by the Free
7 Software Foundation; either version 3, or (at your option) any later
10 GCC is distributed in the hope that it will be useful, but WITHOUT ANY
11 WARRANTY; without even the implied warranty of MERCHANTABILITY or
12 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
15 Under Section 7 of GPL version 3, you are granted additional
16 permissions described in the GCC Runtime Library Exception, version
17 3.1, as published by the Free Software Foundation.
19 You should have received a copy of the GNU General Public License and
20 a copy of the GCC Runtime Library Exception along with this program;
21 see the files COPYING3 and COPYING.RUNTIME respectively. If not, see
22 <http://www.gnu.org/licenses/>. */
24 #include "bid_internal.h"
26 /*****************************************************************************
28 * BID128 non-computational functions:
31 * - bid128_isSubnormal
35 * - bid128_isSignaling
36 * - bid128_isCanonical
44 * - bid128_totalOrderMag
45 * - bid128_sameQuantum
47 ****************************************************************************/
49 #if DECIMAL_CALL_BY_REFERENCE
51 bid128_isSigned (int *pres
,
52 UINT128
* px _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
56 bid128_isSigned (UINT128 x _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
60 res
= ((x
.w
[HIGH_128W
] & MASK_SIGN
) == MASK_SIGN
);
64 // return 1 iff x is not zero, nor NaN nor subnormal nor infinity
65 #if DECIMAL_CALL_BY_REFERENCE
67 bid128_isNormal (int *pres
,
68 UINT128
* px _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
72 bid128_isNormal (UINT128 x _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
75 UINT64 x_exp
, C1_hi
, C1_lo
;
77 int exp
, q
, x_nr_bits
;
80 // test for special values - infinity or NaN
81 if ((x
.w
[1] & MASK_SPECIAL
) == MASK_SPECIAL
) {
87 x_exp
= x
.w
[1] & MASK_EXP
; // biased and shifted left 49 bit positions
88 C1_hi
= x
.w
[1] & MASK_COEFF
;
91 if (C1_hi
== 0 && C1_lo
== 0) {
95 // test for non-canonical values of the argument x
96 if ((((C1_hi
> 0x0001ed09bead87c0ull
)
97 || ((C1_hi
== 0x0001ed09bead87c0ull
)
98 && (C1_lo
> 0x378d8e63ffffffffull
)))
99 && ((x
.w
[1] & 0x6000000000000000ull
) != 0x6000000000000000ull
))
100 || ((x
.w
[1] & 0x6000000000000000ull
) == 0x6000000000000000ull
)) {
104 // x is subnormal or normal
105 // determine the number of digits q in the significand
106 // q = nr. of decimal digits in x
107 // determine first the nr. of bits in x
109 if (C1_lo
>= 0x0020000000000000ull
) { // x >= 2^53
110 // split the 64-bit value in two 32-bit halves to avoid rounding errors
111 if (C1_lo
>= 0x0000000100000000ull
) { // x >= 2^32
112 tmp1
.d
= (double) (C1_lo
>> 32); // exact conversion
114 33 + ((((unsigned int) (tmp1
.ui64
>> 52)) & 0x7ff) - 0x3ff);
116 tmp1
.d
= (double) (C1_lo
); // exact conversion
118 1 + ((((unsigned int) (tmp1
.ui64
>> 52)) & 0x7ff) - 0x3ff);
120 } else { // if x < 2^53
121 tmp1
.d
= (double) C1_lo
; // exact conversion
123 1 + ((((unsigned int) (tmp1
.ui64
>> 52)) & 0x7ff) - 0x3ff);
125 } else { // C1_hi != 0 => nr. bits = 64 + nr_bits (C1_hi)
126 tmp1
.d
= (double) C1_hi
; // exact conversion
128 65 + ((((unsigned int) (tmp1
.ui64
>> 52)) & 0x7ff) - 0x3ff);
130 q
= nr_digits
[x_nr_bits
- 1].digits
;
132 q
= nr_digits
[x_nr_bits
- 1].digits1
;
133 if (C1_hi
> nr_digits
[x_nr_bits
- 1].threshold_hi
||
134 (C1_hi
== nr_digits
[x_nr_bits
- 1].threshold_hi
&&
135 C1_lo
>= nr_digits
[x_nr_bits
- 1].threshold_lo
))
138 exp
= (int) (x_exp
>> 49) - 6176;
139 // test for subnormal values of x
140 if (exp
+ q
<= -6143) {
149 // return 1 iff x is not zero, nor NaN nor normal nor infinity
150 #if DECIMAL_CALL_BY_REFERENCE
152 bid128_isSubnormal (int *pres
,
153 UINT128
* px _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
157 bid128_isSubnormal (UINT128 x _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
160 UINT64 x_exp
, C1_hi
, C1_lo
;
162 int exp
, q
, x_nr_bits
;
165 // test for special values - infinity or NaN
166 if ((x
.w
[1] & MASK_SPECIAL
) == MASK_SPECIAL
) {
172 x_exp
= x
.w
[1] & MASK_EXP
; // biased and shifted left 49 bit positions
173 C1_hi
= x
.w
[1] & MASK_COEFF
;
176 if (C1_hi
== 0 && C1_lo
== 0) {
180 // test for non-canonical values of the argument x
181 if ((((C1_hi
> 0x0001ed09bead87c0ull
)
182 || ((C1_hi
== 0x0001ed09bead87c0ull
)
183 && (C1_lo
> 0x378d8e63ffffffffull
)))
184 && ((x
.w
[1] & 0x6000000000000000ull
) != 0x6000000000000000ull
))
185 || ((x
.w
[1] & 0x6000000000000000ull
) == 0x6000000000000000ull
)) {
189 // x is subnormal or normal
190 // determine the number of digits q in the significand
191 // q = nr. of decimal digits in x
192 // determine first the nr. of bits in x
194 if (C1_lo
>= 0x0020000000000000ull
) { // x >= 2^53
195 // split the 64-bit value in two 32-bit halves to avoid rounding errors
196 if (C1_lo
>= 0x0000000100000000ull
) { // x >= 2^32
197 tmp1
.d
= (double) (C1_lo
>> 32); // exact conversion
199 33 + ((((unsigned int) (tmp1
.ui64
>> 52)) & 0x7ff) - 0x3ff);
201 tmp1
.d
= (double) (C1_lo
); // exact conversion
203 1 + ((((unsigned int) (tmp1
.ui64
>> 52)) & 0x7ff) - 0x3ff);
205 } else { // if x < 2^53
206 tmp1
.d
= (double) C1_lo
; // exact conversion
208 1 + ((((unsigned int) (tmp1
.ui64
>> 52)) & 0x7ff) - 0x3ff);
210 } else { // C1_hi != 0 => nr. bits = 64 + nr_bits (C1_hi)
211 tmp1
.d
= (double) C1_hi
; // exact conversion
213 65 + ((((unsigned int) (tmp1
.ui64
>> 52)) & 0x7ff) - 0x3ff);
215 q
= nr_digits
[x_nr_bits
- 1].digits
;
217 q
= nr_digits
[x_nr_bits
- 1].digits1
;
218 if (C1_hi
> nr_digits
[x_nr_bits
- 1].threshold_hi
||
219 (C1_hi
== nr_digits
[x_nr_bits
- 1].threshold_hi
&&
220 C1_lo
>= nr_digits
[x_nr_bits
- 1].threshold_lo
))
223 exp
= (int) (x_exp
>> 49) - 6176;
224 // test for subnormal values of x
225 if (exp
+ q
<= -6143) {
233 #if DECIMAL_CALL_BY_REFERENCE
235 bid128_isFinite (int *pres
,
236 UINT128
* px _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
240 bid128_isFinite (UINT128 x _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
243 res
= ((x
.w
[HIGH_128W
] & MASK_INF
) != MASK_INF
);
247 #if DECIMAL_CALL_BY_REFERENCE
249 bid128_isZero (int *pres
, UINT128
* px _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
253 bid128_isZero (UINT128 x _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
259 if ((x
.w
[1] & MASK_INF
) == MASK_INF
) {
263 sig_x
.w
[1] = x
.w
[1] & 0x0001ffffffffffffull
;
265 if ((sig_x
.w
[1] > 0x0001ed09bead87c0ull
) || // significand is non-canonical
266 ((sig_x
.w
[1] == 0x0001ed09bead87c0ull
) && (sig_x
.w
[0] > 0x378d8e63ffffffffull
)) || // significand is non-canonical
267 ((x
.w
[1] & 0x6000000000000000ull
) == 0x6000000000000000ull
&& (x
.w
[1] & MASK_INF
) != MASK_INF
) || // significand is non-canonical
268 (sig_x
.w
[1] == 0 && sig_x
.w
[0] == 0)) { // significand is 0
276 #if DECIMAL_CALL_BY_REFERENCE
278 bid128_isInf (int *pres
, UINT128
* px _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
282 bid128_isInf (UINT128 x _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
285 res
= ((x
.w
[HIGH_128W
] & MASK_INF
) == MASK_INF
)
286 && ((x
.w
[HIGH_128W
] & MASK_NAN
) != MASK_NAN
);
290 #if DECIMAL_CALL_BY_REFERENCE
292 bid128_isSignaling (int *pres
,
293 UINT128
* px _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
297 bid128_isSignaling (UINT128 x _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
301 res
= ((x
.w
[HIGH_128W
] & MASK_SNAN
) == MASK_SNAN
);
305 // return 1 iff x is a canonical number ,infinity, or NaN.
306 #if DECIMAL_CALL_BY_REFERENCE
308 bid128_isCanonical (int *pres
,
309 UINT128
* px _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
313 bid128_isCanonical (UINT128 x _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
319 if ((x
.w
[1] & MASK_NAN
) == MASK_NAN
) { // NaN
320 if (x
.w
[1] & 0x01ffc00000000000ull
) {
324 sig_x
.w
[1] = x
.w
[1] & 0x00003fffffffffffull
; // 46 bits
325 sig_x
.w
[0] = x
.w
[0]; // 64 bits
326 // payload must be < 10^33 = 0x0000314dc6448d93_38c15b0a00000000
327 if (sig_x
.w
[1] < 0x0000314dc6448d93ull
328 || (sig_x
.w
[1] == 0x0000314dc6448d93ull
329 && sig_x
.w
[0] < 0x38c15b0a00000000ull
)) {
335 } else if ((x
.w
[1] & MASK_INF
) == MASK_INF
) { // infinity
336 if ((x
.w
[1] & 0x03ffffffffffffffull
) || x
.w
[0]) {
343 // not NaN or infinity; extract significand to ensure it is canonical
344 sig_x
.w
[1] = x
.w
[1] & 0x0001ffffffffffffull
;
346 // a canonical number has a coefficient < 10^34
347 // (0x0001ed09_bead87c0_378d8e64_00000000)
348 if ((sig_x
.w
[1] > 0x0001ed09bead87c0ull
) || // significand is non-canonical
349 ((sig_x
.w
[1] == 0x0001ed09bead87c0ull
) && (sig_x
.w
[0] > 0x378d8e63ffffffffull
)) || // significand is non-canonical
350 ((x
.w
[1] & 0x6000000000000000ull
) == 0x6000000000000000ull
)) {
358 #if DECIMAL_CALL_BY_REFERENCE
360 bid128_isNaN (int *pres
, UINT128
* px _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
364 bid128_isNaN (UINT128 x _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
368 res
= ((x
.w
[HIGH_128W
] & MASK_NAN
) == MASK_NAN
);
372 // copies a floating-point operand x to destination y, with no change
373 #if DECIMAL_CALL_BY_REFERENCE
375 bid128_copy (UINT128
* pres
,
376 UINT128
* px _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
380 bid128_copy (UINT128 x _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
388 // copies a floating-point operand x to destination y, reversing the sign
389 #if DECIMAL_CALL_BY_REFERENCE
391 bid128_negate (UINT128
* pres
,
392 UINT128
* px _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
396 bid128_negate (UINT128 x _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
400 x
.w
[HIGH_128W
] ^= MASK_SIGN
;
405 // copies a floating-point operand x to destination y, changing the sign to positive
406 #if DECIMAL_CALL_BY_REFERENCE
408 bid128_abs (UINT128
* pres
,
409 UINT128
* px _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
413 bid128_abs (UINT128 x _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
417 x
.w
[HIGH_128W
] &= ~MASK_SIGN
;
422 // copies operand x to destination in the same format as x, but with the sign of y
423 #if DECIMAL_CALL_BY_REFERENCE
425 bid128_copySign (UINT128
* pres
, UINT128
* px
,
426 UINT128
* py _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
431 bid128_copySign (UINT128 x
, UINT128 y _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
436 (x
.w
[HIGH_128W
] & ~MASK_SIGN
) | (y
.w
[HIGH_128W
] & MASK_SIGN
);
441 #if DECIMAL_CALL_BY_REFERENCE
443 bid128_class (int *pres
, UINT128
* px _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
447 bid128_class (UINT128 x _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
450 UINT256 sig_x_prime256
;
451 UINT192 sig_x_prime192
;
456 if ((x
.w
[1] & MASK_NAN
) == MASK_NAN
) {
457 if ((x
.w
[1] & MASK_SNAN
) == MASK_SNAN
) {
464 if ((x
.w
[1] & MASK_INF
) == MASK_INF
) {
465 if ((x
.w
[1] & MASK_SIGN
) == MASK_SIGN
) {
466 res
= negativeInfinity
;
468 res
= positiveInfinity
;
472 // decode number into exponent and significand
473 sig_x
.w
[1] = x
.w
[1] & 0x0001ffffffffffffull
;
475 // check for zero or non-canonical
476 if ((sig_x
.w
[1] > 0x0001ed09bead87c0ull
)
477 || ((sig_x
.w
[1] == 0x0001ed09bead87c0ull
)
478 && (sig_x
.w
[0] > 0x378d8e63ffffffffull
))
479 || ((x
.w
[1] & 0x6000000000000000ull
) == 0x6000000000000000ull
)
480 || ((sig_x
.w
[1] == 0) && (sig_x
.w
[0] == 0))) {
481 if ((x
.w
[1] & MASK_SIGN
) == MASK_SIGN
) {
488 exp_x
= (x
.w
[1] >> 49) & 0x000000000003fffull
;
489 // if exponent is less than -6176, the number may be subnormal
490 // (less than the smallest normal value)
491 // the smallest normal value is 1 x 10^-6143 = 10^33 x 10^-6176
492 // if (exp_x - 6176 < -6143)
493 if (exp_x
< 33) { // sig_x * 10^exp_x
495 __mul_128x128_to_256 (sig_x_prime256
, sig_x
,
496 ten2k128
[exp_x
- 20]);
497 // 10^33 = 0x0000314dc6448d93_38c15b0a00000000
498 if ((sig_x_prime256
.w
[3] == 0) && (sig_x_prime256
.w
[2] == 0)
499 && ((sig_x_prime256
.w
[1] < 0x0000314dc6448d93ull
)
500 || ((sig_x_prime256
.w
[1] == 0x0000314dc6448d93ull
)
501 && (sig_x_prime256
.w
[0] < 0x38c15b0a00000000ull
)))) {
502 res
= ((x
.w
[1] & MASK_SIGN
) == MASK_SIGN
) ? negativeSubnormal
:
507 __mul_64x128_to_192 (sig_x_prime192
, ten2k64
[exp_x
], sig_x
);
508 // 10^33 = 0x0000314dc6448d93_38c15b0a00000000
509 if ((sig_x_prime192
.w
[2] == 0)
510 && ((sig_x_prime192
.w
[1] < 0x0000314dc6448d93ull
)
511 || ((sig_x_prime192
.w
[1] == 0x0000314dc6448d93ull
)
512 && (sig_x_prime192
.w
[0] < 0x38c15b0a00000000ull
)))) {
513 res
= ((x
.w
[1] & MASK_SIGN
) == MASK_SIGN
) ? negativeSubnormal
:
519 // otherwise, normal number, determine the sign
521 ((x
.w
[1] & MASK_SIGN
) ==
522 MASK_SIGN
) ? negativeNormal
: positiveNormal
;
526 // true if the exponents of x and y are the same, false otherwise.
527 // The special cases of sameQuantum(NaN, NaN) and sameQuantum(Inf, Inf) are true
528 // If exactly one operand is infinite or exactly one operand is NaN, then false
529 #if DECIMAL_CALL_BY_REFERENCE
531 bid128_sameQuantum (int *pres
, UINT128
* px
,
532 UINT128
* py _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
537 bid128_sameQuantum (UINT128 x
,
538 UINT128 y _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
545 // if both operands are NaN, return true
546 if ((x
.w
[1] & MASK_NAN
) == MASK_NAN
547 || ((y
.w
[1] & MASK_NAN
) == MASK_NAN
)) {
548 res
= ((x
.w
[1] & MASK_NAN
) == MASK_NAN
549 && (y
.w
[1] & MASK_NAN
) == MASK_NAN
);
552 // if both operands are INF, return true
553 if ((x
.w
[1] & MASK_INF
) == MASK_INF
554 || (y
.w
[1] & MASK_INF
) == MASK_INF
) {
555 res
= ((x
.w
[1] & MASK_INF
) == MASK_INF
)
556 && ((y
.w
[1] & MASK_INF
) == MASK_INF
);
559 // decode exponents for both numbers, and return true if they match
560 if ((x
.w
[1] & 0x6000000000000000ull
) == 0x6000000000000000ull
) { // G0_G1=11
561 x_exp
= (x
.w
[1] << 2) & MASK_EXP
; // biased and shifted left 49 bits
562 } else { // G0_G1 != 11
563 x_exp
= x
.w
[1] & MASK_EXP
; // biased and shifted left 49 bits
565 if ((y
.w
[1] & 0x6000000000000000ull
) == 0x6000000000000000ull
) { // G0_G1=11
566 y_exp
= (y
.w
[1] << 2) & MASK_EXP
; // biased and shifted left 49 bits
567 } else { // G0_G1 != 11
568 y_exp
= y
.w
[1] & MASK_EXP
; // biased and shifted left 49 bits
570 res
= (x_exp
== y_exp
);
574 #if DECIMAL_CALL_BY_REFERENCE
576 bid128_totalOrder (int *pres
, UINT128
* px
,
577 UINT128
* py _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
582 bid128_totalOrder (UINT128 x
,
583 UINT128 y _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
587 UINT128 sig_x
, sig_y
, pyld_y
, pyld_x
;
588 UINT192 sig_n_prime192
;
589 UINT256 sig_n_prime256
;
590 char x_is_zero
= 0, y_is_zero
= 0;
595 // if x and y are unordered numerically because either operand is NaN
596 // (1) totalOrder(-NaN, number) is true
597 // (2) totalOrder(number, +NaN) is true
598 // (3) if x and y are both NaN:
599 // i) negative sign bit < positive sign bit
600 // ii) signaling < quiet for +NaN, reverse for -NaN
601 // iii) lesser payload < greater payload for +NaN (reverse for -NaN)
602 // iv) else if bitwise identical (in canonical form), return 1
603 if ((x
.w
[1] & MASK_NAN
) == MASK_NAN
) {
605 if ((x
.w
[1] & MASK_SIGN
) == MASK_SIGN
) {
606 // return true, unless y is -NaN also
607 if ((y
.w
[1] & MASK_NAN
) != MASK_NAN
608 || (y
.w
[1] & MASK_SIGN
) != MASK_SIGN
) {
609 res
= 1; // y is a number, return 1
611 } else { // if y and x are both -NaN
612 pyld_x
.w
[1] = x
.w
[1] & 0x00003fffffffffffull
;
613 pyld_x
.w
[0] = x
.w
[0];
614 pyld_y
.w
[1] = y
.w
[1] & 0x00003fffffffffffull
;
615 pyld_y
.w
[0] = y
.w
[0];
616 if ((pyld_x
.w
[1] > 0x0000314dc6448d93ull
)
617 || ((pyld_x
.w
[1] == 0x0000314dc6448d93ull
)
618 && (pyld_x
.w
[0] > 0x38c15b09ffffffffull
))) {
622 if ((pyld_y
.w
[1] > 0x0000314dc6448d93ull
)
623 || ((pyld_y
.w
[1] == 0x0000314dc6448d93ull
)
624 && (pyld_y
.w
[0] > 0x38c15b09ffffffffull
))) {
628 // if x and y are both -SNaN or both -QNaN, we have to compare payloads
629 // this statement evaluates to true if both are SNaN or QNaN
631 (((y
.w
[1] & MASK_SNAN
) == MASK_SNAN
) ^
632 ((x
.w
[1] & MASK_SNAN
) == MASK_SNAN
))) {
633 // it comes down to the payload. we want to return true if x has a
634 // larger payload, or if the payloads are equal (canonical forms
635 // are bitwise identical)
636 if ((pyld_x
.w
[1] > pyld_y
.w
[1]) ||
637 ((pyld_x
.w
[1] == pyld_y
.w
[1])
638 && (pyld_x
.w
[0] >= pyld_y
.w
[0])))
644 // either x = -SNaN and y = -QNaN or x = -QNaN and y = -SNaN
645 res
= ((y
.w
[1] & MASK_SNAN
) == MASK_SNAN
);
646 // totalOrder (-QNaN, -SNaN) == 1
650 } else { // x is +NaN
651 // return false, unless y is +NaN also
652 if ((y
.w
[1] & MASK_NAN
) != MASK_NAN
653 || (y
.w
[1] & MASK_SIGN
) == MASK_SIGN
) {
654 res
= 0; // y is a number, return 1
657 // x and y are both +NaN;
658 pyld_x
.w
[1] = x
.w
[1] & 0x00003fffffffffffull
;
659 pyld_x
.w
[0] = x
.w
[0];
660 pyld_y
.w
[1] = y
.w
[1] & 0x00003fffffffffffull
;
661 pyld_y
.w
[0] = y
.w
[0];
662 if ((pyld_x
.w
[1] > 0x0000314dc6448d93ull
)
663 || ((pyld_x
.w
[1] == 0x0000314dc6448d93ull
)
664 && (pyld_x
.w
[0] > 0x38c15b09ffffffffull
))) {
668 if ((pyld_y
.w
[1] > 0x0000314dc6448d93ull
)
669 || ((pyld_y
.w
[1] == 0x0000314dc6448d93ull
)
670 && (pyld_y
.w
[0] > 0x38c15b09ffffffffull
))) {
674 // if x and y are both +SNaN or both +QNaN, we have to compare payloads
675 // this statement evaluates to true if both are SNaN or QNaN
677 (((y
.w
[1] & MASK_SNAN
) == MASK_SNAN
) ^
678 ((x
.w
[1] & MASK_SNAN
) == MASK_SNAN
))) {
679 // it comes down to the payload. we want to return true if x has a
680 // smaller payload, or if the payloads are equal (canonical forms
681 // are bitwise identical)
682 if ((pyld_x
.w
[1] < pyld_y
.w
[1]) ||
683 ((pyld_x
.w
[1] == pyld_y
.w
[1])
684 && (pyld_x
.w
[0] <= pyld_y
.w
[0])))
690 // either x = SNaN and y = QNaN or x = QNaN and y = SNaN
691 res
= ((x
.w
[1] & MASK_SNAN
) == MASK_SNAN
);
692 // totalOrder (-QNaN, -SNaN) == 1
697 } else if ((y
.w
[1] & MASK_NAN
) == MASK_NAN
) {
698 // x is certainly not NAN in this case.
699 // return true if y is positive
700 res
= ((y
.w
[1] & MASK_SIGN
) != MASK_SIGN
);
704 // if all the bits are the same, the numbers are equal.
705 if ((x
.w
[1] == y
.w
[1]) && (x
.w
[0] == y
.w
[0])) {
709 // OPPOSITE SIGNS (CASE 3)
710 // if signs are opposite, return 1 if x is negative
711 // (if x < y, totalOrder is true)
712 if (((x
.w
[1] & MASK_SIGN
) == MASK_SIGN
) ^ ((y
.w
[1] & MASK_SIGN
) ==
714 res
= ((x
.w
[1] & MASK_SIGN
) == MASK_SIGN
);
718 if ((x
.w
[1] & MASK_INF
) == MASK_INF
) {
719 // if x == neg_inf, return (y == neg_inf);
720 if ((x
.w
[1] & MASK_SIGN
) == MASK_SIGN
) {
724 // x is positive infinity, only return1 if y is positive infinity as well
725 res
= ((y
.w
[1] & MASK_INF
) == MASK_INF
);
727 // && (y & MASK_SIGN) != MASK_SIGN); (we know y has same sign as x)
729 } else if ((y
.w
[1] & MASK_INF
) == MASK_INF
) {
733 res
= ((y
.w
[1] & MASK_SIGN
) != MASK_SIGN
);
737 sig_x
.w
[1] = x
.w
[1] & 0x0001ffffffffffffull
;
739 exp_x
= (x
.w
[1] >> 49) & 0x000000000003fffull
;
741 // CHECK IF x IS CANONICAL
742 // 9999999999999999999999999999999999 (decimal) =
743 // 1ed09_bead87c0_378d8e63_ffffffff(hexadecimal)
744 // [0, 10^34) is the 754r supported canonical range.
745 // If the value exceeds that, it is interpreted as 0.
746 if ((((sig_x
.w
[1] > 0x0001ed09bead87c0ull
) ||
747 ((sig_x
.w
[1] == 0x0001ed09bead87c0ull
) &&
748 (sig_x
.w
[0] > 0x378d8e63ffffffffull
))) &&
749 ((x
.w
[1] & 0x6000000000000000ull
) != 0x6000000000000000ull
)) ||
750 ((x
.w
[1] & 0x6000000000000000ull
) == 0x6000000000000000ull
) ||
751 ((sig_x
.w
[1] == 0) && (sig_x
.w
[0] == 0))) {
753 // check for the case where the exponent is shifted right by 2 bits!
754 if ((x
.w
[1] & 0x6000000000000000ull
) == 0x6000000000000000ull
) {
755 exp_x
= (x
.w
[1] >> 47) & 0x000000000003fffull
;
759 exp_y
= (y
.w
[1] >> 49) & 0x0000000000003fffull
;
760 sig_y
.w
[1] = y
.w
[1] & 0x0001ffffffffffffull
;
763 // CHECK IF y IS CANONICAL
764 // 9999999999999999999999999999999999(decimal) =
765 // 1ed09_bead87c0_378d8e63_ffffffff(hexadecimal)
766 // [0, 10^34) is the 754r supported canonical range.
767 // If the value exceeds that, it is interpreted as 0.
768 if ((((sig_y
.w
[1] > 0x0001ed09bead87c0ull
) ||
769 ((sig_y
.w
[1] == 0x0001ed09bead87c0ull
) &&
770 (sig_y
.w
[0] > 0x378d8e63ffffffffull
))) &&
771 ((y
.w
[1] & 0x6000000000000000ull
) != 0x6000000000000000ull
)) ||
772 ((y
.w
[1] & 0x6000000000000000ull
) == 0x6000000000000000ull
) ||
773 ((sig_y
.w
[1] == 0) && (sig_y
.w
[0] == 0))) {
775 // check for the case where the exponent is shifted right by 2 bits!
776 if ((y
.w
[1] & 0x6000000000000000ull
) == 0x6000000000000000ull
) {
777 exp_y
= (y
.w
[1] >> 47) & 0x000000000003fffull
;
781 // if x and y represent the same entities, and both are negative
782 // return true iff exp_x <= exp_y
783 if (x_is_zero
&& y_is_zero
) {
784 // we know that signs must be the same because we would have caught it
785 // in case3 if signs were different
786 // totalOrder(x,y) iff exp_x >= exp_y for negative numbers
787 // totalOrder(x,y) iff exp_x <= exp_y for positive numbers
788 if (exp_x
== exp_y
) {
792 res
= ((exp_x
<= exp_y
) ^ ((x
.w
[1] & MASK_SIGN
) == MASK_SIGN
));
795 // if x is zero and y isn't, clearly x has the smaller payload
797 res
= ((y
.w
[1] & MASK_SIGN
) != MASK_SIGN
);
800 // if y is zero, and x isn't, clearly y has the smaller payload
802 res
= ((x
.w
[1] & MASK_SIGN
) == MASK_SIGN
);
805 // REDUNDANT REPRESENTATIONS (CASE 6)
806 // if both components are either bigger or smaller
807 if (((sig_x
.w
[1] > sig_y
.w
[1])
808 || (sig_x
.w
[1] == sig_y
.w
[1] && sig_x
.w
[0] > sig_y
.w
[0]))
810 res
= ((x
.w
[1] & MASK_SIGN
) == MASK_SIGN
);
813 if (((sig_x
.w
[1] < sig_y
.w
[1])
814 || (sig_x
.w
[1] == sig_y
.w
[1] && sig_x
.w
[0] < sig_y
.w
[0]))
816 res
= ((x
.w
[1] & MASK_SIGN
) != MASK_SIGN
);
819 // if |exp_x - exp_y| < 33, it comes down to the compensated significand
821 // if exp_x is 33 greater than exp_y, it is definitely larger,
822 // so no need for compensation
823 if (exp_x
- exp_y
> 33) {
824 res
= ((x
.w
[1] & MASK_SIGN
) == MASK_SIGN
);
826 // difference cannot be greater than 10^33
828 // otherwise adjust the x significand upwards
829 if (exp_x
- exp_y
> 19) {
830 __mul_128x128_to_256 (sig_n_prime256
, sig_x
,
831 ten2k128
[exp_x
- exp_y
- 20]);
832 // the compensated significands are equal (ie "x and y represent the same
833 // entities") return 1 if (negative && expx > expy) ||
834 // (positive && expx < expy)
835 if ((sig_n_prime256
.w
[3] == 0) && (sig_n_prime256
.w
[2] == 0)
836 && (sig_n_prime256
.w
[1] == sig_y
.w
[1])
837 && (sig_n_prime256
.w
[0] == sig_y
.w
[0])) {
838 // the case exp_x == exp_y cannot occur, because all bits must be
839 // the same - would have been caught if (x == y)
840 res
= ((exp_x
<= exp_y
) ^ ((x
.w
[1] & MASK_SIGN
) == MASK_SIGN
));
843 // if positive, return 1 if adjusted x is smaller than y
844 res
= (((sig_n_prime256
.w
[3] == 0) && (sig_n_prime256
.w
[2] == 0)
845 && ((sig_n_prime256
.w
[1] < sig_y
.w
[1])
846 || (sig_n_prime256
.w
[1] == sig_y
.w
[1]
847 && sig_n_prime256
.w
[0] <
848 sig_y
.w
[0]))) ^ ((x
.w
[1] & MASK_SIGN
) ==
852 __mul_64x128_to_192 (sig_n_prime192
, ten2k64
[exp_x
- exp_y
], sig_x
);
853 // if positive, return whichever significand is larger
854 // (converse if negative)
855 if ((sig_n_prime192
.w
[2] == 0) && sig_n_prime192
.w
[1] == sig_y
.w
[1]
856 && (sig_n_prime192
.w
[0] == sig_y
.w
[0])) {
857 res
= ((exp_x
<= exp_y
) ^ ((x
.w
[1] & MASK_SIGN
) == MASK_SIGN
));
860 res
= (((sig_n_prime192
.w
[2] == 0)
861 && ((sig_n_prime192
.w
[1] < sig_y
.w
[1])
862 || (sig_n_prime192
.w
[1] == sig_y
.w
[1]
863 && sig_n_prime192
.w
[0] <
864 sig_y
.w
[0]))) ^ ((x
.w
[1] & MASK_SIGN
) ==
868 // if exp_x is 33 less than exp_y, it is definitely smaller,
869 // no need for compensation
870 if (exp_y
- exp_x
> 33) {
871 res
= ((x
.w
[1] & MASK_SIGN
) != MASK_SIGN
);
874 if (exp_y
- exp_x
> 19) {
875 // adjust the y significand upwards
876 __mul_128x128_to_256 (sig_n_prime256
, sig_y
,
877 ten2k128
[exp_y
- exp_x
- 20]);
878 // if x and y represent the same entities and both are negative
879 // return true iff exp_x <= exp_y
880 if ((sig_n_prime256
.w
[3] == 0) && (sig_n_prime256
.w
[2] == 0)
881 && (sig_n_prime256
.w
[1] == sig_x
.w
[1])
882 && (sig_n_prime256
.w
[0] == sig_x
.w
[0])) {
883 res
= (exp_x
<= exp_y
) ^ ((x
.w
[1] & MASK_SIGN
) == MASK_SIGN
);
886 // values are not equal, for positive numbers return 1 if x is less than y
888 res
= (((sig_n_prime256
.w
[3] != 0) ||
889 // if upper128 bits of compensated y are non-zero, y is bigger
890 (sig_n_prime256
.w
[2] != 0) ||
891 // if upper128 bits of compensated y are non-zero, y is bigger
892 (sig_n_prime256
.w
[1] > sig_x
.w
[1]) ||
893 // if compensated y is bigger, y is bigger
894 (sig_n_prime256
.w
[1] == sig_x
.w
[1]
895 && sig_n_prime256
.w
[0] >
896 sig_x
.w
[0])) ^ ((x
.w
[1] & MASK_SIGN
) == MASK_SIGN
));
899 __mul_64x128_to_192 (sig_n_prime192
, ten2k64
[exp_y
- exp_x
], sig_y
);
900 if ((sig_n_prime192
.w
[2] == 0) && (sig_n_prime192
.w
[1] == sig_x
.w
[1])
901 && (sig_n_prime192
.w
[0] == sig_x
.w
[0])) {
902 res
= (exp_x
<= exp_y
) ^ ((x
.w
[1] & MASK_SIGN
) == MASK_SIGN
);
905 res
= (((sig_n_prime192
.w
[2] != 0) ||
906 // if upper128 bits of compensated y are non-zero, y is bigger
907 (sig_n_prime192
.w
[1] > sig_x
.w
[1]) ||
908 // if compensated y is bigger, y is bigger
909 (sig_n_prime192
.w
[1] == sig_x
.w
[1]
910 && sig_n_prime192
.w
[0] >
911 sig_x
.w
[0])) ^ ((x
.w
[1] & MASK_SIGN
) == MASK_SIGN
));
915 #if DECIMAL_CALL_BY_REFERENCE
917 bid128_totalOrderMag (int *pres
, UINT128
* px
,
918 UINT128
* py _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
923 bid128_totalOrderMag (UINT128 x
,
924 UINT128 y _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
928 UINT128 sig_x
, sig_y
, pyld_y
, pyld_x
;
929 UINT192 sig_n_prime192
;
930 UINT256 sig_n_prime256
;
931 char x_is_zero
= 0, y_is_zero
= 0;
935 x
.w
[1] = x
.w
[1] & 0x7fffffffffffffffull
;
936 y
.w
[1] = y
.w
[1] & 0x7fffffffffffffffull
;
939 // if x and y are unordered numerically because either operand is NaN
940 // (1) totalOrder(number, +NaN) is true
941 // (2) if x and y are both NaN:
942 // i) signaling < quiet for +NaN
943 // ii) lesser payload < greater payload for +NaN
944 // iii) else if bitwise identical (in canonical form), return 1
945 if ((x
.w
[1] & MASK_NAN
) == MASK_NAN
) {
947 // return false, unless y is +NaN also
948 if ((y
.w
[1] & MASK_NAN
) != MASK_NAN
) {
949 res
= 0; // y is a number, return 0
952 // x and y are both +NaN;
953 pyld_x
.w
[1] = x
.w
[1] & 0x00003fffffffffffull
;
954 pyld_x
.w
[0] = x
.w
[0];
955 pyld_y
.w
[1] = y
.w
[1] & 0x00003fffffffffffull
;
956 pyld_y
.w
[0] = y
.w
[0];
957 if ((pyld_x
.w
[1] > 0x0000314dc6448d93ull
)
958 || ((pyld_x
.w
[1] == 0x0000314dc6448d93ull
)
959 && (pyld_x
.w
[0] > 0x38c15b09ffffffffull
))) {
963 if ((pyld_y
.w
[1] > 0x0000314dc6448d93ull
)
964 || ((pyld_y
.w
[1] == 0x0000314dc6448d93ull
)
965 && (pyld_y
.w
[0] > 0x38c15b09ffffffffull
))) {
969 // if x and y are both +SNaN or both +QNaN, we have to compare payloads
970 // this statement evaluates to true if both are SNaN or QNaN
972 (((y
.w
[1] & MASK_SNAN
) == MASK_SNAN
) ^
973 ((x
.w
[1] & MASK_SNAN
) == MASK_SNAN
))) {
974 // it comes down to the payload. we want to return true if x has a
975 // smaller payload, or if the payloads are equal (canonical forms
976 // are bitwise identical)
977 if ((pyld_x
.w
[1] < pyld_y
.w
[1]) ||
978 ((pyld_x
.w
[1] == pyld_y
.w
[1])
979 && (pyld_x
.w
[0] <= pyld_y
.w
[0]))) {
986 // either x = SNaN and y = QNaN or x = QNaN and y = SNaN
987 res
= ((x
.w
[1] & MASK_SNAN
) == MASK_SNAN
);
988 // totalOrder (-QNaN, -SNaN) == 1
992 } else if ((y
.w
[1] & MASK_NAN
) == MASK_NAN
) {
993 // x is certainly not NAN in this case.
994 // return true because y is positive
999 // if all the bits are the same, the numbers are equal.
1000 if ((x
.w
[1] == y
.w
[1]) && (x
.w
[0] == y
.w
[0])) {
1004 // INFINITY (CASE 3)
1005 if ((x
.w
[1] & MASK_INF
) == MASK_INF
) {
1006 // x is positive infinity, only return 1 if y is positive infinity as well
1007 res
= ((y
.w
[1] & MASK_INF
) == MASK_INF
);
1009 // (we know y has same sign as x)
1010 } else if ((y
.w
[1] & MASK_INF
) == MASK_INF
) {
1012 // since y is +inf, x<y
1020 sig_x
.w
[1] = x
.w
[1] & 0x0001ffffffffffffull
;
1021 sig_x
.w
[0] = x
.w
[0];
1022 exp_x
= (x
.w
[1] >> 49) & 0x000000000003fffull
;
1024 // CHECK IF x IS CANONICAL
1025 // 9999999999999999999999999999999999 (decimal) =
1026 // 1ed09_bead87c0_378d8e63_ffffffff(hexadecimal)
1027 // [0, 10^34) is the 754r supported canonical range.
1028 // If the value exceeds that, it is interpreted as 0.
1029 if ((((sig_x
.w
[1] > 0x0001ed09bead87c0ull
) ||
1030 ((sig_x
.w
[1] == 0x0001ed09bead87c0ull
) &&
1031 (sig_x
.w
[0] > 0x378d8e63ffffffffull
))) &&
1032 ((x
.w
[1] & 0x6000000000000000ull
) != 0x6000000000000000ull
)) ||
1033 ((x
.w
[1] & 0x6000000000000000ull
) == 0x6000000000000000ull
) ||
1034 ((sig_x
.w
[1] == 0) && (sig_x
.w
[0] == 0))) {
1036 // check for the case where the exponent is shifted right by 2 bits!
1037 if ((x
.w
[1] & 0x6000000000000000ull
) == 0x6000000000000000ull
) {
1038 exp_x
= (x
.w
[1] >> 47) & 0x000000000003fffull
;
1042 exp_y
= (y
.w
[1] >> 49) & 0x0000000000003fffull
;
1043 sig_y
.w
[1] = y
.w
[1] & 0x0001ffffffffffffull
;
1044 sig_y
.w
[0] = y
.w
[0];
1046 // CHECK IF y IS CANONICAL
1047 // 9999999999999999999999999999999999(decimal) =
1048 // 1ed09_bead87c0_378d8e63_ffffffff(hexadecimal)
1049 // [0, 10^34) is the 754r supported canonical range.
1050 // If the value exceeds that, it is interpreted as 0.
1051 if ((((sig_y
.w
[1] > 0x0001ed09bead87c0ull
) ||
1052 ((sig_y
.w
[1] == 0x0001ed09bead87c0ull
) &&
1053 (sig_y
.w
[0] > 0x378d8e63ffffffffull
))) &&
1054 ((y
.w
[1] & 0x6000000000000000ull
) != 0x6000000000000000ull
)) ||
1055 ((y
.w
[1] & 0x6000000000000000ull
) == 0x6000000000000000ull
) ||
1056 ((sig_y
.w
[1] == 0) && (sig_y
.w
[0] == 0))) {
1058 // check for the case where the exponent is shifted right by 2 bits!
1059 if ((y
.w
[1] & 0x6000000000000000ull
) == 0x6000000000000000ull
) {
1060 exp_y
= (y
.w
[1] >> 47) & 0x000000000003fffull
;
1064 if (x_is_zero
&& y_is_zero
) {
1065 // we know that signs must be the same because we would have caught it
1066 // in case3 if signs were different
1067 // totalOrder(x,y) iff exp_x <= exp_y for positive numbers
1068 if (exp_x
== exp_y
) {
1072 res
= (exp_x
<= exp_y
);
1075 // if x is zero and y isn't, clearly x has the smaller payload
1080 // if y is zero, and x isn't, clearly y has the smaller payload
1085 // REDUNDANT REPRESENTATIONS (CASE 5)
1086 // if both components are either bigger or smaller
1087 if (((sig_x
.w
[1] > sig_y
.w
[1])
1088 || (sig_x
.w
[1] == sig_y
.w
[1] && sig_x
.w
[0] > sig_y
.w
[0]))
1089 && exp_x
>= exp_y
) {
1093 if (((sig_x
.w
[1] < sig_y
.w
[1])
1094 || (sig_x
.w
[1] == sig_y
.w
[1] && sig_x
.w
[0] < sig_y
.w
[0]))
1095 && exp_x
<= exp_y
) {
1099 // if |exp_x - exp_y| < 33, it comes down to the compensated significand
1100 if (exp_x
> exp_y
) {
1101 // if exp_x is 33 greater than exp_y, it is definitely larger,
1102 // so no need for compensation
1103 if (exp_x
- exp_y
> 33) {
1104 res
= 0; // difference cannot be greater than 10^33
1107 // otherwise adjust the x significand upwards
1108 if (exp_x
- exp_y
> 19) {
1109 __mul_128x128_to_256 (sig_n_prime256
, sig_x
,
1110 ten2k128
[exp_x
- exp_y
- 20]);
1111 // the compensated significands are equal (ie "x and y represent the same
1112 // entities") return 1 if (negative && expx > expy) ||
1113 // (positive && expx < expy)
1114 if ((sig_n_prime256
.w
[3] == 0) && (sig_n_prime256
.w
[2] == 0)
1115 && (sig_n_prime256
.w
[1] == sig_y
.w
[1])
1116 && (sig_n_prime256
.w
[0] == sig_y
.w
[0])) {
1117 // the case (exp_x == exp_y) cannot occur, because all bits must be
1118 // the same - would have been caught if (x == y)
1119 res
= (exp_x
<= exp_y
);
1122 // since positive, return 1 if adjusted x is smaller than y
1123 res
= ((sig_n_prime256
.w
[3] == 0) && (sig_n_prime256
.w
[2] == 0)
1124 && ((sig_n_prime256
.w
[1] < sig_y
.w
[1])
1125 || (sig_n_prime256
.w
[1] == sig_y
.w
[1]
1126 && sig_n_prime256
.w
[0] < sig_y
.w
[0])));
1129 __mul_64x128_to_192 (sig_n_prime192
, ten2k64
[exp_x
- exp_y
], sig_x
);
1130 // if positive, return whichever significand is larger
1131 // (converse if negative)
1132 if ((sig_n_prime192
.w
[2] == 0) && sig_n_prime192
.w
[1] == sig_y
.w
[1]
1133 && (sig_n_prime192
.w
[0] == sig_y
.w
[0])) {
1134 res
= (exp_x
<= exp_y
);
1137 res
= ((sig_n_prime192
.w
[2] == 0)
1138 && ((sig_n_prime192
.w
[1] < sig_y
.w
[1])
1139 || (sig_n_prime192
.w
[1] == sig_y
.w
[1]
1140 && sig_n_prime192
.w
[0] < sig_y
.w
[0])));
1143 // if exp_x is 33 less than exp_y, it is definitely smaller,
1144 // no need for compensation
1145 if (exp_y
- exp_x
> 33) {
1149 if (exp_y
- exp_x
> 19) {
1150 // adjust the y significand upwards
1151 __mul_128x128_to_256 (sig_n_prime256
, sig_y
,
1152 ten2k128
[exp_y
- exp_x
- 20]);
1153 if ((sig_n_prime256
.w
[3] == 0) && (sig_n_prime256
.w
[2] == 0)
1154 && (sig_n_prime256
.w
[1] == sig_x
.w
[1])
1155 && (sig_n_prime256
.w
[0] == sig_x
.w
[0])) {
1156 res
= (exp_x
<= exp_y
);
1159 // values are not equal, for positive numbers return 1 if x is less than y
1161 res
= ((sig_n_prime256
.w
[3] != 0) ||
1162 // if upper128 bits of compensated y are non-zero, y is bigger
1163 (sig_n_prime256
.w
[2] != 0) ||
1164 // if upper128 bits of compensated y are non-zero, y is bigger
1165 (sig_n_prime256
.w
[1] > sig_x
.w
[1]) ||
1166 // if compensated y is bigger, y is bigger
1167 (sig_n_prime256
.w
[1] == sig_x
.w
[1]
1168 && sig_n_prime256
.w
[0] > sig_x
.w
[0]));
1171 __mul_64x128_to_192 (sig_n_prime192
, ten2k64
[exp_y
- exp_x
], sig_y
);
1172 if ((sig_n_prime192
.w
[2] == 0) && (sig_n_prime192
.w
[1] == sig_x
.w
[1])
1173 && (sig_n_prime192
.w
[0] == sig_x
.w
[0])) {
1174 res
= (exp_x
<= exp_y
);
1177 res
= ((sig_n_prime192
.w
[2] != 0) ||
1178 // if upper128 bits of compensated y are non-zero, y is bigger
1179 (sig_n_prime192
.w
[1] > sig_x
.w
[1]) ||
1180 // if compensated y is bigger, y is bigger
1181 (sig_n_prime192
.w
[1] == sig_x
.w
[1]
1182 && sig_n_prime192
.w
[0] > sig_x
.w
[0]));
1186 #if DECIMAL_CALL_BY_REFERENCE
1188 bid128_radix (int *pres
, UINT128
* px _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
1192 bid128_radix (UINT128 x _EXC_MASKS_PARAM _EXC_INFO_PARAM
) {
1195 if (x
.w
[LOW_128W
]) // dummy test