1 /* real.c - software floating point emulation.
2 Copyright (C) 1993, 1994, 1995, 1996, 1997, 1998, 1999,
3 2000, 2002, 2003 Free Software Foundation, Inc.
4 Contributed by Stephen L. Moshier (moshier@world.std.com).
5 Re-written by Richard Henderson <rth@redhat.com>
7 This file is part of GCC.
9 GCC is free software; you can redistribute it and/or modify it under
10 the terms of the GNU General Public License as published by the Free
11 Software Foundation; either version 2, or (at your option) any later
14 GCC is distributed in the hope that it will be useful, but WITHOUT ANY
15 WARRANTY; without even the implied warranty of MERCHANTABILITY or
16 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19 You should have received a copy of the GNU General Public License
20 along with GCC; see the file COPYING. If not, write to the Free
21 Software Foundation, 59 Temple Place - Suite 330, Boston, MA
26 #include "coretypes.h"
33 /* The floating point model used internally is not exactly IEEE 754
34 compliant, and close to the description in the ISO C standard,
35 section 5.2.4.2.2 Characteristics of floating types.
39 x = s * b^e * \sum_{k=1}^p f_k * b^{-k}
43 b = base or radix, here always 2
45 p = precision (the number of base-b digits in the significand)
46 f_k = the digits of the significand.
48 We differ from typical IEEE 754 encodings in that the entire
49 significand is fractional. Normalized significands are in the
52 A requirement of the model is that P be larger than than the
53 largest supported target floating-point type by at least 2 bits.
54 This gives us proper rounding when we truncate to the target type.
55 In addition, E must be large enough to hold the smallest supported
56 denormal number in a normalized form.
58 Both of these requirements are easily satisfied. The largest target
59 significand is 113 bits; we store at least 160. The smallest
60 denormal number fits in 17 exponent bits; we store 29.
62 Note that the decimal string conversion routines are sensitive to
63 rounding error. Since the raw arithmetic routines do not themselves
64 have guard digits or rounding, the computation of 10**exp can
65 accumulate more than a few digits of error. The previous incarnation
66 of real.c successfully used a 144 bit fraction; given the current
67 layout of REAL_VALUE_TYPE we're forced to expand to at least 160 bits.
69 Target floating point models that use base 16 instead of base 2
70 (i.e. IBM 370), are handled during round_for_format, in which we
71 canonicalize the exponent to be a multiple of 4 (log2(16)), and
72 adjust the significand to match. */
75 /* Used to classify two numbers simultaneously. */
76 #define CLASS2(A, B) ((A) << 2 | (B))
78 #if HOST_BITS_PER_LONG != 64 && HOST_BITS_PER_LONG != 32
79 #error "Some constant folding done by hand to avoid shift count warnings"
82 static void get_zero
PARAMS ((REAL_VALUE_TYPE
*, int));
83 static void get_canonical_qnan
PARAMS ((REAL_VALUE_TYPE
*, int));
84 static void get_canonical_snan
PARAMS ((REAL_VALUE_TYPE
*, int));
85 static void get_inf
PARAMS ((REAL_VALUE_TYPE
*, int));
86 static bool sticky_rshift_significand
PARAMS ((REAL_VALUE_TYPE
*,
87 const REAL_VALUE_TYPE
*,
89 static void rshift_significand
PARAMS ((REAL_VALUE_TYPE
*,
90 const REAL_VALUE_TYPE
*,
92 static void lshift_significand
PARAMS ((REAL_VALUE_TYPE
*,
93 const REAL_VALUE_TYPE
*,
95 static void lshift_significand_1
PARAMS ((REAL_VALUE_TYPE
*,
96 const REAL_VALUE_TYPE
*));
97 static bool add_significands
PARAMS ((REAL_VALUE_TYPE
*r
,
98 const REAL_VALUE_TYPE
*,
99 const REAL_VALUE_TYPE
*));
100 static bool sub_significands
PARAMS ((REAL_VALUE_TYPE
*,
101 const REAL_VALUE_TYPE
*,
102 const REAL_VALUE_TYPE
*, int));
103 static void neg_significand
PARAMS ((REAL_VALUE_TYPE
*,
104 const REAL_VALUE_TYPE
*));
105 static int cmp_significands
PARAMS ((const REAL_VALUE_TYPE
*,
106 const REAL_VALUE_TYPE
*));
107 static int cmp_significand_0
PARAMS ((const REAL_VALUE_TYPE
*));
108 static void set_significand_bit
PARAMS ((REAL_VALUE_TYPE
*, unsigned int));
109 static void clear_significand_bit
PARAMS ((REAL_VALUE_TYPE
*, unsigned int));
110 static bool test_significand_bit
PARAMS ((REAL_VALUE_TYPE
*, unsigned int));
111 static void clear_significand_below
PARAMS ((REAL_VALUE_TYPE
*,
113 static bool div_significands
PARAMS ((REAL_VALUE_TYPE
*,
114 const REAL_VALUE_TYPE
*,
115 const REAL_VALUE_TYPE
*));
116 static void normalize
PARAMS ((REAL_VALUE_TYPE
*));
118 static void do_add
PARAMS ((REAL_VALUE_TYPE
*, const REAL_VALUE_TYPE
*,
119 const REAL_VALUE_TYPE
*, int));
120 static void do_multiply
PARAMS ((REAL_VALUE_TYPE
*,
121 const REAL_VALUE_TYPE
*,
122 const REAL_VALUE_TYPE
*));
123 static void do_divide
PARAMS ((REAL_VALUE_TYPE
*, const REAL_VALUE_TYPE
*,
124 const REAL_VALUE_TYPE
*));
125 static int do_compare
PARAMS ((const REAL_VALUE_TYPE
*,
126 const REAL_VALUE_TYPE
*, int));
127 static void do_fix_trunc
PARAMS ((REAL_VALUE_TYPE
*, const REAL_VALUE_TYPE
*));
129 static unsigned long rtd_divmod
PARAMS ((REAL_VALUE_TYPE
*,
132 static const REAL_VALUE_TYPE
* ten_to_ptwo
PARAMS ((int));
133 static const REAL_VALUE_TYPE
* ten_to_mptwo
PARAMS ((int));
134 static const REAL_VALUE_TYPE
* real_digit
PARAMS ((int));
135 static void times_pten
PARAMS ((REAL_VALUE_TYPE
*, int));
137 static void round_for_format
PARAMS ((const struct real_format
*,
140 /* Initialize R with a positive zero. */
147 memset (r
, 0, sizeof (*r
));
151 /* Initialize R with the canonical quiet NaN. */
154 get_canonical_qnan (r
, sign
)
158 memset (r
, 0, sizeof (*r
));
161 r
->sig
[SIGSZ
-1] = SIG_MSB
>> 1;
165 get_canonical_snan (r
, sign
)
169 memset (r
, 0, sizeof (*r
));
172 r
->sig
[SIGSZ
-1] = SIG_MSB
>> 2;
180 memset (r
, 0, sizeof (*r
));
186 /* Right-shift the significand of A by N bits; put the result in the
187 significand of R. If any one bits are shifted out, return true. */
190 sticky_rshift_significand (r
, a
, n
)
192 const REAL_VALUE_TYPE
*a
;
195 unsigned long sticky
= 0;
196 unsigned int i
, ofs
= 0;
198 if (n
>= HOST_BITS_PER_LONG
)
200 for (i
= 0, ofs
= n
/ HOST_BITS_PER_LONG
; i
< ofs
; ++i
)
202 n
&= HOST_BITS_PER_LONG
- 1;
207 sticky
|= a
->sig
[ofs
] & (((unsigned long)1 << n
) - 1);
208 for (i
= 0; i
< SIGSZ
; ++i
)
211 = (((ofs
+ i
>= SIGSZ
? 0 : a
->sig
[ofs
+ i
]) >> n
)
212 | ((ofs
+ i
+ 1 >= SIGSZ
? 0 : a
->sig
[ofs
+ i
+ 1])
213 << (HOST_BITS_PER_LONG
- n
)));
218 for (i
= 0; ofs
+ i
< SIGSZ
; ++i
)
219 r
->sig
[i
] = a
->sig
[ofs
+ i
];
220 for (; i
< SIGSZ
; ++i
)
227 /* Right-shift the significand of A by N bits; put the result in the
231 rshift_significand (r
, a
, n
)
233 const REAL_VALUE_TYPE
*a
;
236 unsigned int i
, ofs
= n
/ HOST_BITS_PER_LONG
;
238 n
&= HOST_BITS_PER_LONG
- 1;
241 for (i
= 0; i
< SIGSZ
; ++i
)
244 = (((ofs
+ i
>= SIGSZ
? 0 : a
->sig
[ofs
+ i
]) >> n
)
245 | ((ofs
+ i
+ 1 >= SIGSZ
? 0 : a
->sig
[ofs
+ i
+ 1])
246 << (HOST_BITS_PER_LONG
- n
)));
251 for (i
= 0; ofs
+ i
< SIGSZ
; ++i
)
252 r
->sig
[i
] = a
->sig
[ofs
+ i
];
253 for (; i
< SIGSZ
; ++i
)
258 /* Left-shift the significand of A by N bits; put the result in the
262 lshift_significand (r
, a
, n
)
264 const REAL_VALUE_TYPE
*a
;
267 unsigned int i
, ofs
= n
/ HOST_BITS_PER_LONG
;
269 n
&= HOST_BITS_PER_LONG
- 1;
272 for (i
= 0; ofs
+ i
< SIGSZ
; ++i
)
273 r
->sig
[SIGSZ
-1-i
] = a
->sig
[SIGSZ
-1-i
-ofs
];
274 for (; i
< SIGSZ
; ++i
)
275 r
->sig
[SIGSZ
-1-i
] = 0;
278 for (i
= 0; i
< SIGSZ
; ++i
)
281 = (((ofs
+ i
>= SIGSZ
? 0 : a
->sig
[SIGSZ
-1-i
-ofs
]) << n
)
282 | ((ofs
+ i
+ 1 >= SIGSZ
? 0 : a
->sig
[SIGSZ
-1-i
-ofs
-1])
283 >> (HOST_BITS_PER_LONG
- n
)));
287 /* Likewise, but N is specialized to 1. */
290 lshift_significand_1 (r
, a
)
292 const REAL_VALUE_TYPE
*a
;
296 for (i
= SIGSZ
- 1; i
> 0; --i
)
297 r
->sig
[i
] = (a
->sig
[i
] << 1) | (a
->sig
[i
-1] >> (HOST_BITS_PER_LONG
- 1));
298 r
->sig
[0] = a
->sig
[0] << 1;
301 /* Add the significands of A and B, placing the result in R. Return
302 true if there was carry out of the most significant word. */
305 add_significands (r
, a
, b
)
307 const REAL_VALUE_TYPE
*a
, *b
;
312 for (i
= 0; i
< SIGSZ
; ++i
)
314 unsigned long ai
= a
->sig
[i
];
315 unsigned long ri
= ai
+ b
->sig
[i
];
331 /* Subtract the significands of A and B, placing the result in R. CARRY is
332 true if there's a borrow incoming to the least significant word.
333 Return true if there was borrow out of the most significant word. */
336 sub_significands (r
, a
, b
, carry
)
338 const REAL_VALUE_TYPE
*a
, *b
;
343 for (i
= 0; i
< SIGSZ
; ++i
)
345 unsigned long ai
= a
->sig
[i
];
346 unsigned long ri
= ai
- b
->sig
[i
];
362 /* Negate the significand A, placing the result in R. */
365 neg_significand (r
, a
)
367 const REAL_VALUE_TYPE
*a
;
372 for (i
= 0; i
< SIGSZ
; ++i
)
374 unsigned long ri
, ai
= a
->sig
[i
];
393 /* Compare significands. Return tri-state vs zero. */
396 cmp_significands (a
, b
)
397 const REAL_VALUE_TYPE
*a
, *b
;
401 for (i
= SIGSZ
- 1; i
>= 0; --i
)
403 unsigned long ai
= a
->sig
[i
];
404 unsigned long bi
= b
->sig
[i
];
415 /* Return true if A is nonzero. */
418 cmp_significand_0 (a
)
419 const REAL_VALUE_TYPE
*a
;
423 for (i
= SIGSZ
- 1; i
>= 0; --i
)
430 /* Set bit N of the significand of R. */
433 set_significand_bit (r
, n
)
437 r
->sig
[n
/ HOST_BITS_PER_LONG
]
438 |= (unsigned long)1 << (n
% HOST_BITS_PER_LONG
);
441 /* Clear bit N of the significand of R. */
444 clear_significand_bit (r
, n
)
448 r
->sig
[n
/ HOST_BITS_PER_LONG
]
449 &= ~((unsigned long)1 << (n
% HOST_BITS_PER_LONG
));
452 /* Test bit N of the significand of R. */
455 test_significand_bit (r
, n
)
459 /* ??? Compiler bug here if we return this expression directly.
460 The conversion to bool strips the "&1" and we wind up testing
461 e.g. 2 != 0 -> true. Seen in gcc version 3.2 20020520. */
462 int t
= (r
->sig
[n
/ HOST_BITS_PER_LONG
] >> (n
% HOST_BITS_PER_LONG
)) & 1;
466 /* Clear bits 0..N-1 of the significand of R. */
469 clear_significand_below (r
, n
)
473 int i
, w
= n
/ HOST_BITS_PER_LONG
;
475 for (i
= 0; i
< w
; ++i
)
478 r
->sig
[w
] &= ~(((unsigned long)1 << (n
% HOST_BITS_PER_LONG
)) - 1);
481 /* Divide the significands of A and B, placing the result in R. Return
482 true if the division was inexact. */
485 div_significands (r
, a
, b
)
487 const REAL_VALUE_TYPE
*a
, *b
;
490 int i
, bit
= SIGNIFICAND_BITS
- 1;
491 unsigned long msb
, inexact
;
494 memset (r
->sig
, 0, sizeof (r
->sig
));
500 msb
= u
.sig
[SIGSZ
-1] & SIG_MSB
;
501 lshift_significand_1 (&u
, &u
);
503 if (msb
|| cmp_significands (&u
, b
) >= 0)
505 sub_significands (&u
, &u
, b
, 0);
506 set_significand_bit (r
, bit
);
511 for (i
= 0, inexact
= 0; i
< SIGSZ
; i
++)
517 /* Adjust the exponent and significand of R such that the most
518 significant bit is set. We underflow to zero and overflow to
519 infinity here, without denormals. (The intermediate representation
520 exponent is large enough to handle target denormals normalized.) */
529 /* Find the first word that is nonzero. */
530 for (i
= SIGSZ
- 1; i
>= 0; i
--)
532 shift
+= HOST_BITS_PER_LONG
;
536 /* Zero significand flushes to zero. */
544 /* Find the first bit that is nonzero. */
546 if (r
->sig
[i
] & ((unsigned long)1 << (HOST_BITS_PER_LONG
- 1 - j
)))
552 exp
= r
->exp
- shift
;
554 get_inf (r
, r
->sign
);
555 else if (exp
< -MAX_EXP
)
556 get_zero (r
, r
->sign
);
560 lshift_significand (r
, r
, shift
);
565 /* Return R = A + (SUBTRACT_P ? -B : B). */
568 do_add (r
, a
, b
, subtract_p
)
570 const REAL_VALUE_TYPE
*a
, *b
;
575 bool inexact
= false;
577 /* Determine if we need to add or subtract. */
579 subtract_p
= (sign
^ b
->sign
) ^ subtract_p
;
581 switch (CLASS2 (a
->class, b
->class))
583 case CLASS2 (rvc_zero
, rvc_zero
):
584 /* -0 + -0 = -0, -0 - +0 = -0; all other cases yield +0. */
585 get_zero (r
, sign
& !subtract_p
);
588 case CLASS2 (rvc_zero
, rvc_normal
):
589 case CLASS2 (rvc_zero
, rvc_inf
):
590 case CLASS2 (rvc_zero
, rvc_nan
):
592 case CLASS2 (rvc_normal
, rvc_nan
):
593 case CLASS2 (rvc_inf
, rvc_nan
):
594 case CLASS2 (rvc_nan
, rvc_nan
):
595 /* ANY + NaN = NaN. */
596 case CLASS2 (rvc_normal
, rvc_inf
):
599 r
->sign
= sign
^ subtract_p
;
602 case CLASS2 (rvc_normal
, rvc_zero
):
603 case CLASS2 (rvc_inf
, rvc_zero
):
604 case CLASS2 (rvc_nan
, rvc_zero
):
606 case CLASS2 (rvc_nan
, rvc_normal
):
607 case CLASS2 (rvc_nan
, rvc_inf
):
608 /* NaN + ANY = NaN. */
609 case CLASS2 (rvc_inf
, rvc_normal
):
614 case CLASS2 (rvc_inf
, rvc_inf
):
616 /* Inf - Inf = NaN. */
617 get_canonical_qnan (r
, 0);
619 /* Inf + Inf = Inf. */
623 case CLASS2 (rvc_normal
, rvc_normal
):
630 /* Swap the arguments such that A has the larger exponent. */
631 dexp
= a
->exp
- b
->exp
;
634 const REAL_VALUE_TYPE
*t
;
641 /* If the exponents are not identical, we need to shift the
642 significand of B down. */
645 /* If the exponents are too far apart, the significands
646 do not overlap, which makes the subtraction a noop. */
647 if (dexp
>= SIGNIFICAND_BITS
)
654 inexact
|= sticky_rshift_significand (&t
, b
, dexp
);
660 if (sub_significands (r
, a
, b
, inexact
))
662 /* We got a borrow out of the subtraction. That means that
663 A and B had the same exponent, and B had the larger
664 significand. We need to swap the sign and negate the
667 neg_significand (r
, r
);
672 if (add_significands (r
, a
, b
))
674 /* We got carry out of the addition. This means we need to
675 shift the significand back down one bit and increase the
677 inexact
|= sticky_rshift_significand (r
, r
, 1);
678 r
->sig
[SIGSZ
-1] |= SIG_MSB
;
687 r
->class = rvc_normal
;
691 /* Re-normalize the result. */
694 /* Special case: if the subtraction results in zero, the result
696 if (r
->class == rvc_zero
)
699 r
->sig
[0] |= inexact
;
702 /* Return R = A * B. */
705 do_multiply (r
, a
, b
)
707 const REAL_VALUE_TYPE
*a
, *b
;
709 REAL_VALUE_TYPE u
, t
, *rr
;
710 unsigned int i
, j
, k
;
711 int sign
= a
->sign
^ b
->sign
;
713 switch (CLASS2 (a
->class, b
->class))
715 case CLASS2 (rvc_zero
, rvc_zero
):
716 case CLASS2 (rvc_zero
, rvc_normal
):
717 case CLASS2 (rvc_normal
, rvc_zero
):
718 /* +-0 * ANY = 0 with appropriate sign. */
722 case CLASS2 (rvc_zero
, rvc_nan
):
723 case CLASS2 (rvc_normal
, rvc_nan
):
724 case CLASS2 (rvc_inf
, rvc_nan
):
725 case CLASS2 (rvc_nan
, rvc_nan
):
726 /* ANY * NaN = NaN. */
731 case CLASS2 (rvc_nan
, rvc_zero
):
732 case CLASS2 (rvc_nan
, rvc_normal
):
733 case CLASS2 (rvc_nan
, rvc_inf
):
734 /* NaN * ANY = NaN. */
739 case CLASS2 (rvc_zero
, rvc_inf
):
740 case CLASS2 (rvc_inf
, rvc_zero
):
742 get_canonical_qnan (r
, sign
);
745 case CLASS2 (rvc_inf
, rvc_inf
):
746 case CLASS2 (rvc_normal
, rvc_inf
):
747 case CLASS2 (rvc_inf
, rvc_normal
):
748 /* Inf * Inf = Inf, R * Inf = Inf */
753 case CLASS2 (rvc_normal
, rvc_normal
):
760 if (r
== a
|| r
== b
)
766 /* Collect all the partial products. Since we don't have sure access
767 to a widening multiply, we split each long into two half-words.
769 Consider the long-hand form of a four half-word multiplication:
779 We construct partial products of the widened half-word products
780 that are known to not overlap, e.g. DF+DH. Each such partial
781 product is given its proper exponent, which allows us to sum them
782 and obtain the finished product. */
784 for (i
= 0; i
< SIGSZ
* 2; ++i
)
786 unsigned long ai
= a
->sig
[i
/ 2];
788 ai
>>= HOST_BITS_PER_LONG
/ 2;
790 ai
&= ((unsigned long)1 << (HOST_BITS_PER_LONG
/ 2)) - 1;
795 for (j
= 0; j
< 2; ++j
)
797 int exp
= (a
->exp
- (2*SIGSZ
-1-i
)*(HOST_BITS_PER_LONG
/2)
798 + (b
->exp
- (1-j
)*(HOST_BITS_PER_LONG
/2)));
803 /* Would underflow to zero, which we shouldn't bother adding. */
806 u
.class = rvc_normal
;
810 for (k
= j
; k
< SIGSZ
* 2; k
+= 2)
812 unsigned long bi
= b
->sig
[k
/ 2];
814 bi
>>= HOST_BITS_PER_LONG
/ 2;
816 bi
&= ((unsigned long)1 << (HOST_BITS_PER_LONG
/ 2)) - 1;
818 u
.sig
[k
/ 2] = ai
* bi
;
822 do_add (rr
, rr
, &u
, 0);
831 /* Return R = A / B. */
836 const REAL_VALUE_TYPE
*a
, *b
;
838 int exp
, sign
= a
->sign
^ b
->sign
;
839 REAL_VALUE_TYPE t
, *rr
;
842 switch (CLASS2 (a
->class, b
->class))
844 case CLASS2 (rvc_zero
, rvc_zero
):
846 case CLASS2 (rvc_inf
, rvc_inf
):
847 /* Inf / Inf = NaN. */
848 get_canonical_qnan (r
, sign
);
851 case CLASS2 (rvc_zero
, rvc_normal
):
852 case CLASS2 (rvc_zero
, rvc_inf
):
854 case CLASS2 (rvc_normal
, rvc_inf
):
860 case CLASS2 (rvc_normal
, rvc_zero
):
862 case CLASS2 (rvc_inf
, rvc_zero
):
867 case CLASS2 (rvc_zero
, rvc_nan
):
868 case CLASS2 (rvc_normal
, rvc_nan
):
869 case CLASS2 (rvc_inf
, rvc_nan
):
870 case CLASS2 (rvc_nan
, rvc_nan
):
871 /* ANY / NaN = NaN. */
876 case CLASS2 (rvc_nan
, rvc_zero
):
877 case CLASS2 (rvc_nan
, rvc_normal
):
878 case CLASS2 (rvc_nan
, rvc_inf
):
879 /* NaN / ANY = NaN. */
884 case CLASS2 (rvc_inf
, rvc_normal
):
890 case CLASS2 (rvc_normal
, rvc_normal
):
897 if (r
== a
|| r
== b
)
902 rr
->class = rvc_normal
;
905 exp
= a
->exp
- b
->exp
+ 1;
912 inexact
= div_significands (rr
, a
, b
);
914 /* Re-normalize the result. */
916 rr
->sig
[0] |= inexact
;
922 /* Return a tri-state comparison of A vs B. Return NAN_RESULT if
923 one of the two operands is a NaN. */
926 do_compare (a
, b
, nan_result
)
927 const REAL_VALUE_TYPE
*a
, *b
;
932 switch (CLASS2 (a
->class, b
->class))
934 case CLASS2 (rvc_zero
, rvc_zero
):
935 /* Sign of zero doesn't matter for compares. */
938 case CLASS2 (rvc_inf
, rvc_zero
):
939 case CLASS2 (rvc_inf
, rvc_normal
):
940 case CLASS2 (rvc_normal
, rvc_zero
):
941 return (a
->sign
? -1 : 1);
943 case CLASS2 (rvc_inf
, rvc_inf
):
944 return -a
->sign
- -b
->sign
;
946 case CLASS2 (rvc_zero
, rvc_normal
):
947 case CLASS2 (rvc_zero
, rvc_inf
):
948 case CLASS2 (rvc_normal
, rvc_inf
):
949 return (b
->sign
? 1 : -1);
951 case CLASS2 (rvc_zero
, rvc_nan
):
952 case CLASS2 (rvc_normal
, rvc_nan
):
953 case CLASS2 (rvc_inf
, rvc_nan
):
954 case CLASS2 (rvc_nan
, rvc_nan
):
955 case CLASS2 (rvc_nan
, rvc_zero
):
956 case CLASS2 (rvc_nan
, rvc_normal
):
957 case CLASS2 (rvc_nan
, rvc_inf
):
960 case CLASS2 (rvc_normal
, rvc_normal
):
967 if (a
->sign
!= b
->sign
)
968 return -a
->sign
- -b
->sign
;
972 else if (a
->exp
< b
->exp
)
975 ret
= cmp_significands (a
, b
);
977 return (a
->sign
? -ret
: ret
);
980 /* Return A truncated to an integral value toward zero. */
985 const REAL_VALUE_TYPE
*a
;
998 get_zero (r
, r
->sign
);
999 else if (r
->exp
< SIGNIFICAND_BITS
)
1000 clear_significand_below (r
, SIGNIFICAND_BITS
- r
->exp
);
1008 /* Perform the binary or unary operation described by CODE.
1009 For a unary operation, leave OP1 NULL. */
1012 real_arithmetic (r
, icode
, op0
, op1
)
1015 const REAL_VALUE_TYPE
*op0
, *op1
;
1017 enum tree_code code
= icode
;
1022 do_add (r
, op0
, op1
, 0);
1026 do_add (r
, op0
, op1
, 1);
1030 do_multiply (r
, op0
, op1
);
1034 do_divide (r
, op0
, op1
);
1038 if (op1
->class == rvc_nan
)
1040 else if (do_compare (op0
, op1
, -1) < 0)
1047 if (op1
->class == rvc_nan
)
1049 else if (do_compare (op0
, op1
, 1) < 0)
1065 case FIX_TRUNC_EXPR
:
1066 do_fix_trunc (r
, op0
);
1074 /* Legacy. Similar, but return the result directly. */
1077 real_arithmetic2 (icode
, op0
, op1
)
1079 const REAL_VALUE_TYPE
*op0
, *op1
;
1082 real_arithmetic (&r
, icode
, op0
, op1
);
1087 real_compare (icode
, op0
, op1
)
1089 const REAL_VALUE_TYPE
*op0
, *op1
;
1091 enum tree_code code
= icode
;
1096 return do_compare (op0
, op1
, 1) < 0;
1098 return do_compare (op0
, op1
, 1) <= 0;
1100 return do_compare (op0
, op1
, -1) > 0;
1102 return do_compare (op0
, op1
, -1) >= 0;
1104 return do_compare (op0
, op1
, -1) == 0;
1106 return do_compare (op0
, op1
, -1) != 0;
1107 case UNORDERED_EXPR
:
1108 return op0
->class == rvc_nan
|| op1
->class == rvc_nan
;
1110 return op0
->class != rvc_nan
&& op1
->class != rvc_nan
;
1112 return do_compare (op0
, op1
, -1) < 0;
1114 return do_compare (op0
, op1
, -1) <= 0;
1116 return do_compare (op0
, op1
, 1) > 0;
1118 return do_compare (op0
, op1
, 1) >= 0;
1120 return do_compare (op0
, op1
, 0) == 0;
1127 /* Return floor log2(R). */
1131 const REAL_VALUE_TYPE
*r
;
1139 return (unsigned int)-1 >> 1;
1147 /* R = OP0 * 2**EXP. */
1150 real_ldexp (r
, op0
, exp
)
1152 const REAL_VALUE_TYPE
*op0
;
1166 get_inf (r
, r
->sign
);
1167 else if (exp
< -MAX_EXP
)
1168 get_zero (r
, r
->sign
);
1178 /* Determine whether a floating-point value X is infinite. */
1182 const REAL_VALUE_TYPE
*r
;
1184 return (r
->class == rvc_inf
);
1187 /* Determine whether a floating-point value X is a NaN. */
1191 const REAL_VALUE_TYPE
*r
;
1193 return (r
->class == rvc_nan
);
1196 /* Determine whether a floating-point value X is negative. */
1200 const REAL_VALUE_TYPE
*r
;
1205 /* Determine whether a floating-point value X is minus zero. */
1209 const REAL_VALUE_TYPE
*r
;
1211 return r
->sign
&& r
->class == rvc_zero
;
1214 /* Compare two floating-point objects for bitwise identity. */
1217 real_identical (a
, b
)
1218 const REAL_VALUE_TYPE
*a
, *b
;
1222 if (a
->class != b
->class)
1224 if (a
->sign
!= b
->sign
)
1234 if (a
->exp
!= b
->exp
)
1238 for (i
= 0; i
< SIGSZ
; ++i
)
1239 if (a
->sig
[i
] != b
->sig
[i
])
1250 /* Try to change R into its exact multiplicative inverse in machine
1251 mode MODE. Return true if successful. */
1254 exact_real_inverse (mode
, r
)
1255 enum machine_mode mode
;
1258 const REAL_VALUE_TYPE
*one
= real_digit (1);
1262 if (r
->class != rvc_normal
)
1265 /* Check for a power of two: all significand bits zero except the MSB. */
1266 for (i
= 0; i
< SIGSZ
-1; ++i
)
1269 if (r
->sig
[SIGSZ
-1] != SIG_MSB
)
1272 /* Find the inverse and truncate to the required mode. */
1273 do_divide (&u
, one
, r
);
1274 real_convert (&u
, mode
, &u
);
1276 /* The rounding may have overflowed. */
1277 if (u
.class != rvc_normal
)
1279 for (i
= 0; i
< SIGSZ
-1; ++i
)
1282 if (u
.sig
[SIGSZ
-1] != SIG_MSB
)
1289 /* Render R as an integer. */
1293 const REAL_VALUE_TYPE
*r
;
1295 unsigned HOST_WIDE_INT i
;
1306 i
= (unsigned HOST_WIDE_INT
) 1 << (HOST_BITS_PER_WIDE_INT
- 1);
1314 /* Only force overflow for unsigned overflow. Signed overflow is
1315 undefined, so it doesn't matter what we return, and some callers
1316 expect to be able to use this routine for both signed and
1317 unsigned conversions. */
1318 if (r
->exp
> HOST_BITS_PER_WIDE_INT
)
1321 if (HOST_BITS_PER_WIDE_INT
== HOST_BITS_PER_LONG
)
1322 i
= r
->sig
[SIGSZ
-1];
1323 else if (HOST_BITS_PER_WIDE_INT
== 2*HOST_BITS_PER_LONG
)
1325 i
= r
->sig
[SIGSZ
-1];
1326 i
= i
<< (HOST_BITS_PER_LONG
- 1) << 1;
1327 i
|= r
->sig
[SIGSZ
-2];
1332 i
>>= HOST_BITS_PER_WIDE_INT
- r
->exp
;
1343 /* Likewise, but to an integer pair, HI+LOW. */
1346 real_to_integer2 (plow
, phigh
, r
)
1347 HOST_WIDE_INT
*plow
, *phigh
;
1348 const REAL_VALUE_TYPE
*r
;
1351 HOST_WIDE_INT low
, high
;
1364 high
= (unsigned HOST_WIDE_INT
) 1 << (HOST_BITS_PER_WIDE_INT
- 1);
1378 /* Only force overflow for unsigned overflow. Signed overflow is
1379 undefined, so it doesn't matter what we return, and some callers
1380 expect to be able to use this routine for both signed and
1381 unsigned conversions. */
1382 if (exp
> 2*HOST_BITS_PER_WIDE_INT
)
1385 rshift_significand (&t
, r
, 2*HOST_BITS_PER_WIDE_INT
- exp
);
1386 if (HOST_BITS_PER_WIDE_INT
== HOST_BITS_PER_LONG
)
1388 high
= t
.sig
[SIGSZ
-1];
1389 low
= t
.sig
[SIGSZ
-2];
1391 else if (HOST_BITS_PER_WIDE_INT
== 2*HOST_BITS_PER_LONG
)
1393 high
= t
.sig
[SIGSZ
-1];
1394 high
= high
<< (HOST_BITS_PER_LONG
- 1) << 1;
1395 high
|= t
.sig
[SIGSZ
-2];
1397 low
= t
.sig
[SIGSZ
-3];
1398 low
= low
<< (HOST_BITS_PER_LONG
- 1) << 1;
1399 low
|= t
.sig
[SIGSZ
-4];
1409 low
= -low
, high
= ~high
;
1421 /* A subroutine of real_to_decimal. Compute the quotient and remainder
1422 of NUM / DEN. Return the quotient and place the remainder in NUM.
1423 It is expected that NUM / DEN are close enough that the quotient is
1426 static unsigned long
1427 rtd_divmod (num
, den
)
1428 REAL_VALUE_TYPE
*num
, *den
;
1430 unsigned long q
, msb
;
1431 int expn
= num
->exp
, expd
= den
->exp
;
1440 msb
= num
->sig
[SIGSZ
-1] & SIG_MSB
;
1442 lshift_significand_1 (num
, num
);
1444 if (msb
|| cmp_significands (num
, den
) >= 0)
1446 sub_significands (num
, num
, den
, 0);
1450 while (--expn
>= expd
);
1458 /* Render R as a decimal floating point constant. Emit DIGITS significant
1459 digits in the result, bounded by BUF_SIZE. If DIGITS is 0, choose the
1460 maximum for the representation. If CROP_TRAILING_ZEROS, strip trailing
1463 #define M_LOG10_2 0.30102999566398119521
1466 real_to_decimal (str
, r_orig
, buf_size
, digits
, crop_trailing_zeros
)
1468 const REAL_VALUE_TYPE
*r_orig
;
1469 size_t buf_size
, digits
;
1470 int crop_trailing_zeros
;
1472 const REAL_VALUE_TYPE
*one
, *ten
;
1473 REAL_VALUE_TYPE r
, pten
, u
, v
;
1474 int dec_exp
, cmp_one
, digit
;
1476 char *p
, *first
, *last
;
1483 strcpy (str
, (r
.sign
? "-0.0" : "0.0"));
1488 strcpy (str
, (r
.sign
? "-Inf" : "+Inf"));
1491 /* ??? Print the significand as well, if not canonical? */
1492 strcpy (str
, (r
.sign
? "-NaN" : "+NaN"));
1498 /* Bound the number of digits printed by the size of the representation. */
1499 max_digits
= SIGNIFICAND_BITS
* M_LOG10_2
;
1500 if (digits
== 0 || digits
> max_digits
)
1501 digits
= max_digits
;
1503 /* Estimate the decimal exponent, and compute the length of the string it
1504 will print as. Be conservative and add one to account for possible
1505 overflow or rounding error. */
1506 dec_exp
= r
.exp
* M_LOG10_2
;
1507 for (max_digits
= 1; dec_exp
; max_digits
++)
1510 /* Bound the number of digits printed by the size of the output buffer. */
1511 max_digits
= buf_size
- 1 - 1 - 2 - max_digits
- 1;
1512 if (max_digits
> buf_size
)
1514 if (digits
> max_digits
)
1515 digits
= max_digits
;
1517 one
= real_digit (1);
1518 ten
= ten_to_ptwo (0);
1526 cmp_one
= do_compare (&r
, one
, 0);
1531 /* Number is greater than one. Convert significand to an integer
1532 and strip trailing decimal zeros. */
1535 u
.exp
= SIGNIFICAND_BITS
- 1;
1537 /* Largest M, such that 10**2**M fits within SIGNIFICAND_BITS. */
1538 m
= floor_log2 (max_digits
);
1540 /* Iterate over the bits of the possible powers of 10 that might
1541 be present in U and eliminate them. That is, if we find that
1542 10**2**M divides U evenly, keep the division and increase
1548 do_divide (&t
, &u
, ten_to_ptwo (m
));
1549 do_fix_trunc (&v
, &t
);
1550 if (cmp_significands (&v
, &t
) == 0)
1558 /* Revert the scaling to integer that we performed earlier. */
1559 u
.exp
+= r
.exp
- (SIGNIFICAND_BITS
- 1);
1562 /* Find power of 10. Do this by dividing out 10**2**M when
1563 this is larger than the current remainder. Fill PTEN with
1564 the power of 10 that we compute. */
1567 m
= floor_log2 ((int)(r
.exp
* M_LOG10_2
)) + 1;
1570 const REAL_VALUE_TYPE
*ptentwo
= ten_to_ptwo (m
);
1571 if (do_compare (&u
, ptentwo
, 0) >= 0)
1573 do_divide (&u
, &u
, ptentwo
);
1574 do_multiply (&pten
, &pten
, ptentwo
);
1581 /* We managed to divide off enough tens in the above reduction
1582 loop that we've now got a negative exponent. Fall into the
1583 less-than-one code to compute the proper value for PTEN. */
1590 /* Number is less than one. Pad significand with leading
1596 /* Stop if we'd shift bits off the bottom. */
1600 do_multiply (&u
, &v
, ten
);
1602 /* Stop if we're now >= 1. */
1611 /* Find power of 10. Do this by multiplying in P=10**2**M when
1612 the current remainder is smaller than 1/P. Fill PTEN with the
1613 power of 10 that we compute. */
1614 m
= floor_log2 ((int)(-r
.exp
* M_LOG10_2
)) + 1;
1617 const REAL_VALUE_TYPE
*ptentwo
= ten_to_ptwo (m
);
1618 const REAL_VALUE_TYPE
*ptenmtwo
= ten_to_mptwo (m
);
1620 if (do_compare (&v
, ptenmtwo
, 0) <= 0)
1622 do_multiply (&v
, &v
, ptentwo
);
1623 do_multiply (&pten
, &pten
, ptentwo
);
1629 /* Invert the positive power of 10 that we've collected so far. */
1630 do_divide (&pten
, one
, &pten
);
1638 /* At this point, PTEN should contain the nearest power of 10 smaller
1639 than R, such that this division produces the first digit.
1641 Using a divide-step primitive that returns the complete integral
1642 remainder avoids the rounding error that would be produced if
1643 we were to use do_divide here and then simply multiply by 10 for
1644 each subsequent digit. */
1646 digit
= rtd_divmod (&r
, &pten
);
1648 /* Be prepared for error in that division via underflow ... */
1649 if (digit
== 0 && cmp_significand_0 (&r
))
1651 /* Multiply by 10 and try again. */
1652 do_multiply (&r
, &r
, ten
);
1653 digit
= rtd_divmod (&r
, &pten
);
1659 /* ... or overflow. */
1667 else if (digit
> 10)
1672 /* Generate subsequent digits. */
1673 while (--digits
> 0)
1675 do_multiply (&r
, &r
, ten
);
1676 digit
= rtd_divmod (&r
, &pten
);
1681 /* Generate one more digit with which to do rounding. */
1682 do_multiply (&r
, &r
, ten
);
1683 digit
= rtd_divmod (&r
, &pten
);
1685 /* Round the result. */
1688 /* Round to nearest. If R is nonzero there are additional
1689 nonzero digits to be extracted. */
1690 if (cmp_significand_0 (&r
))
1692 /* Round to even. */
1693 else if ((p
[-1] - '0') & 1)
1710 /* Carry out of the first digit. This means we had all 9's and
1711 now have all 0's. "Prepend" a 1 by overwriting the first 0. */
1719 /* Insert the decimal point. */
1720 first
[0] = first
[1];
1723 /* If requested, drop trailing zeros. Never crop past "1.0". */
1724 if (crop_trailing_zeros
)
1725 while (last
> first
+ 3 && last
[-1] == '0')
1728 /* Append the exponent. */
1729 sprintf (last
, "e%+d", dec_exp
);
1732 /* Render R as a hexadecimal floating point constant. Emit DIGITS
1733 significant digits in the result, bounded by BUF_SIZE. If DIGITS is 0,
1734 choose the maximum for the representation. If CROP_TRAILING_ZEROS,
1735 strip trailing zeros. */
1738 real_to_hexadecimal (str
, r
, buf_size
, digits
, crop_trailing_zeros
)
1740 const REAL_VALUE_TYPE
*r
;
1741 size_t buf_size
, digits
;
1742 int crop_trailing_zeros
;
1744 int i
, j
, exp
= r
->exp
;
1757 strcpy (str
, (r
->sign
? "-Inf" : "+Inf"));
1760 /* ??? Print the significand as well, if not canonical? */
1761 strcpy (str
, (r
->sign
? "-NaN" : "+NaN"));
1768 digits
= SIGNIFICAND_BITS
/ 4;
1770 /* Bound the number of digits printed by the size of the output buffer. */
1772 sprintf (exp_buf
, "p%+d", exp
);
1773 max_digits
= buf_size
- strlen (exp_buf
) - r
->sign
- 4 - 1;
1774 if (max_digits
> buf_size
)
1776 if (digits
> max_digits
)
1777 digits
= max_digits
;
1788 for (i
= SIGSZ
- 1; i
>= 0; --i
)
1789 for (j
= HOST_BITS_PER_LONG
- 4; j
>= 0; j
-= 4)
1791 *p
++ = "0123456789abcdef"[(r
->sig
[i
] >> j
) & 15];
1797 if (crop_trailing_zeros
)
1798 while (p
> first
+ 1 && p
[-1] == '0')
1801 sprintf (p
, "p%+d", exp
);
1804 /* Initialize R from a decimal or hexadecimal string. The string is
1805 assumed to have been syntax checked already. */
1808 real_from_string (r
, str
)
1822 else if (*str
== '+')
1825 if (str
[0] == '0' && str
[1] == 'x')
1827 /* Hexadecimal floating point. */
1828 int pos
= SIGNIFICAND_BITS
- 4, d
;
1836 d
= hex_value (*str
);
1841 r
->sig
[pos
/ HOST_BITS_PER_LONG
]
1842 |= (unsigned long) d
<< (pos
% HOST_BITS_PER_LONG
);
1851 if (pos
== SIGNIFICAND_BITS
- 4)
1858 d
= hex_value (*str
);
1863 r
->sig
[pos
/ HOST_BITS_PER_LONG
]
1864 |= (unsigned long) d
<< (pos
% HOST_BITS_PER_LONG
);
1870 if (*str
== 'p' || *str
== 'P')
1872 bool exp_neg
= false;
1880 else if (*str
== '+')
1884 while (ISDIGIT (*str
))
1890 /* Overflowed the exponent. */
1904 r
->class = rvc_normal
;
1911 /* Decimal floating point. */
1912 const REAL_VALUE_TYPE
*ten
= ten_to_ptwo (0);
1917 while (ISDIGIT (*str
))
1920 do_multiply (r
, r
, ten
);
1922 do_add (r
, r
, real_digit (d
), 0);
1927 if (r
->class == rvc_zero
)
1932 while (ISDIGIT (*str
))
1935 do_multiply (r
, r
, ten
);
1937 do_add (r
, r
, real_digit (d
), 0);
1942 if (*str
== 'e' || *str
== 'E')
1944 bool exp_neg
= false;
1952 else if (*str
== '+')
1956 while (ISDIGIT (*str
))
1962 /* Overflowed the exponent. */
1976 times_pten (r
, exp
);
1991 /* Legacy. Similar, but return the result directly. */
1994 real_from_string2 (s
, mode
)
1996 enum machine_mode mode
;
2000 real_from_string (&r
, s
);
2001 if (mode
!= VOIDmode
)
2002 real_convert (&r
, mode
, &r
);
2007 /* Initialize R from the integer pair HIGH+LOW. */
2010 real_from_integer (r
, mode
, low
, high
, unsigned_p
)
2012 enum machine_mode mode
;
2013 unsigned HOST_WIDE_INT low
;
2017 if (low
== 0 && high
== 0)
2021 r
->class = rvc_normal
;
2022 r
->sign
= high
< 0 && !unsigned_p
;
2023 r
->exp
= 2 * HOST_BITS_PER_WIDE_INT
;
2034 if (HOST_BITS_PER_LONG
== HOST_BITS_PER_WIDE_INT
)
2036 r
->sig
[SIGSZ
-1] = high
;
2037 r
->sig
[SIGSZ
-2] = low
;
2038 memset (r
->sig
, 0, sizeof(long)*(SIGSZ
-2));
2040 else if (HOST_BITS_PER_LONG
*2 == HOST_BITS_PER_WIDE_INT
)
2042 r
->sig
[SIGSZ
-1] = high
>> (HOST_BITS_PER_LONG
- 1) >> 1;
2043 r
->sig
[SIGSZ
-2] = high
;
2044 r
->sig
[SIGSZ
-3] = low
>> (HOST_BITS_PER_LONG
- 1) >> 1;
2045 r
->sig
[SIGSZ
-4] = low
;
2047 memset (r
->sig
, 0, sizeof(long)*(SIGSZ
-4));
2055 if (mode
!= VOIDmode
)
2056 real_convert (r
, mode
, r
);
2059 /* Returns 10**2**N. */
2061 static const REAL_VALUE_TYPE
*
2065 static REAL_VALUE_TYPE tens
[EXP_BITS
];
2067 if (n
< 0 || n
>= EXP_BITS
)
2070 if (tens
[n
].class == rvc_zero
)
2072 if (n
< (HOST_BITS_PER_WIDE_INT
== 64 ? 5 : 4))
2074 HOST_WIDE_INT t
= 10;
2077 for (i
= 0; i
< n
; ++i
)
2080 real_from_integer (&tens
[n
], VOIDmode
, t
, 0, 1);
2084 const REAL_VALUE_TYPE
*t
= ten_to_ptwo (n
- 1);
2085 do_multiply (&tens
[n
], t
, t
);
2092 /* Returns 10**(-2**N). */
2094 static const REAL_VALUE_TYPE
*
2098 static REAL_VALUE_TYPE tens
[EXP_BITS
];
2100 if (n
< 0 || n
>= EXP_BITS
)
2103 if (tens
[n
].class == rvc_zero
)
2104 do_divide (&tens
[n
], real_digit (1), ten_to_ptwo (n
));
2111 static const REAL_VALUE_TYPE
*
2115 static REAL_VALUE_TYPE num
[10];
2120 if (n
> 0 && num
[n
].class == rvc_zero
)
2121 real_from_integer (&num
[n
], VOIDmode
, n
, 0, 1);
2126 /* Multiply R by 10**EXP. */
2133 REAL_VALUE_TYPE pten
, *rr
;
2134 bool negative
= (exp
< 0);
2140 pten
= *real_digit (1);
2146 for (i
= 0; exp
> 0; ++i
, exp
>>= 1)
2148 do_multiply (rr
, rr
, ten_to_ptwo (i
));
2151 do_divide (r
, r
, &pten
);
2154 /* Fills R with +Inf. */
2163 /* Fills R with a NaN whose significand is described by STR. If QUIET,
2164 we force a QNaN, else we force an SNaN. The string, if not empty,
2165 is parsed as a number and placed in the significand. Return true
2166 if the string was successfully parsed. */
2169 real_nan (r
, str
, quiet
, mode
)
2173 enum machine_mode mode
;
2175 const struct real_format
*fmt
;
2177 fmt
= real_format_for_mode
[mode
- QFmode
];
2184 get_canonical_qnan (r
, 0);
2186 get_canonical_snan (r
, 0);
2193 memset (r
, 0, sizeof (*r
));
2196 /* Parse akin to strtol into the significand of R. */
2198 while (ISSPACE (*str
))
2202 else if (*str
== '+')
2212 while ((d
= hex_value (*str
)) < base
)
2219 lshift_significand (r
, r
, 3);
2222 lshift_significand (r
, r
, 4);
2225 lshift_significand_1 (&u
, r
);
2226 lshift_significand (r
, r
, 3);
2227 add_significands (r
, r
, &u
);
2235 add_significands (r
, r
, &u
);
2240 /* Must have consumed the entire string for success. */
2244 /* Shift the significand into place such that the bits
2245 are in the most significant bits for the format. */
2246 lshift_significand (r
, r
, SIGNIFICAND_BITS
- fmt
->p
);
2248 /* Our MSB is always unset for NaNs. */
2249 r
->sig
[SIGSZ
-1] &= ~SIG_MSB
;
2251 /* Force quiet or signalling NaN. */
2253 r
->sig
[SIGSZ
-1] |= SIG_MSB
>> 1;
2255 r
->sig
[SIGSZ
-1] &= ~(SIG_MSB
>> 1);
2257 /* Force at least one bit of the significand set. */
2258 for (d
= 0; d
< SIGSZ
; ++d
)
2262 r
->sig
[SIGSZ
-1] |= SIG_MSB
>> 2;
2264 /* Our intermediate format forces QNaNs to have MSB-1 set.
2265 If the target format has QNaNs with the top bit unset,
2266 mirror the output routines and invert the top two bits. */
2267 if (!fmt
->qnan_msb_set
)
2268 r
->sig
[SIGSZ
-1] ^= (SIG_MSB
>> 1) | (SIG_MSB
>> 2);
2274 /* Fills R with 2**N. */
2281 memset (r
, 0, sizeof (*r
));
2286 else if (n
< -MAX_EXP
)
2290 r
->class = rvc_normal
;
2292 r
->sig
[SIGSZ
-1] = SIG_MSB
;
2298 round_for_format (fmt
, r
)
2299 const struct real_format
*fmt
;
2303 unsigned long sticky
;
2307 p2
= fmt
->p
* fmt
->log2_b
;
2308 emin2m1
= (fmt
->emin
- 1) * fmt
->log2_b
;
2309 emax2
= fmt
->emax
* fmt
->log2_b
;
2311 np2
= SIGNIFICAND_BITS
- p2
;
2315 get_zero (r
, r
->sign
);
2317 if (!fmt
->has_signed_zero
)
2322 get_inf (r
, r
->sign
);
2327 clear_significand_below (r
, np2
);
2329 /* If we've cleared the entire significand, we need one bit
2330 set for this to continue to be a NaN. */
2331 for (i
= 0; i
< SIGSZ
; ++i
)
2335 r
->sig
[SIGSZ
-1] = SIG_MSB
>> 2;
2345 /* If we're not base2, normalize the exponent to a multiple of
2347 if (fmt
->log2_b
!= 1)
2349 int shift
= r
->exp
& (fmt
->log2_b
- 1);
2352 shift
= fmt
->log2_b
- shift
;
2353 r
->sig
[0] |= sticky_rshift_significand (r
, r
, shift
);
2358 /* Check the range of the exponent. If we're out of range,
2359 either underflow or overflow. */
2362 else if (r
->exp
<= emin2m1
)
2366 if (!fmt
->has_denorm
)
2368 /* Don't underflow completely until we've had a chance to round. */
2369 if (r
->exp
< emin2m1
)
2374 diff
= emin2m1
- r
->exp
+ 1;
2378 /* De-normalize the significand. */
2379 r
->sig
[0] |= sticky_rshift_significand (r
, r
, diff
);
2384 /* There are P2 true significand bits, followed by one guard bit,
2385 followed by one sticky bit, followed by stuff. Fold nonzero
2386 stuff into the sticky bit. */
2389 for (i
= 0, w
= (np2
- 1) / HOST_BITS_PER_LONG
; i
< w
; ++i
)
2390 sticky
|= r
->sig
[i
];
2392 r
->sig
[w
] & (((unsigned long)1 << ((np2
- 1) % HOST_BITS_PER_LONG
)) - 1);
2394 guard
= test_significand_bit (r
, np2
- 1);
2395 lsb
= test_significand_bit (r
, np2
);
2397 /* Round to even. */
2398 if (guard
&& (sticky
|| lsb
))
2402 set_significand_bit (&u
, np2
);
2404 if (add_significands (r
, r
, &u
))
2406 /* Overflow. Means the significand had been all ones, and
2407 is now all zeros. Need to increase the exponent, and
2408 possibly re-normalize it. */
2409 if (++r
->exp
> emax2
)
2411 r
->sig
[SIGSZ
-1] = SIG_MSB
;
2413 if (fmt
->log2_b
!= 1)
2415 int shift
= r
->exp
& (fmt
->log2_b
- 1);
2418 shift
= fmt
->log2_b
- shift
;
2419 rshift_significand (r
, r
, shift
);
2428 /* Catch underflow that we deferred until after rounding. */
2429 if (r
->exp
<= emin2m1
)
2432 /* Clear out trailing garbage. */
2433 clear_significand_below (r
, np2
);
2436 /* Extend or truncate to a new mode. */
2439 real_convert (r
, mode
, a
)
2441 enum machine_mode mode
;
2442 const REAL_VALUE_TYPE
*a
;
2444 const struct real_format
*fmt
;
2446 fmt
= real_format_for_mode
[mode
- QFmode
];
2451 round_for_format (fmt
, r
);
2453 /* round_for_format de-normalizes denormals. Undo just that part. */
2454 if (r
->class == rvc_normal
)
2458 /* Legacy. Likewise, except return the struct directly. */
2461 real_value_truncate (mode
, a
)
2462 enum machine_mode mode
;
2466 real_convert (&r
, mode
, &a
);
2470 /* Return true if truncating to MODE is exact. */
2473 exact_real_truncate (mode
, a
)
2474 enum machine_mode mode
;
2475 const REAL_VALUE_TYPE
*a
;
2478 real_convert (&t
, mode
, a
);
2479 return real_identical (&t
, a
);
2482 /* Write R to the given target format. Place the words of the result
2483 in target word order in BUF. There are always 32 bits in each
2484 long, no matter the size of the host long.
2486 Legacy: return word 0 for implementing REAL_VALUE_TO_TARGET_SINGLE. */
2489 real_to_target_fmt (buf
, r_orig
, fmt
)
2491 const REAL_VALUE_TYPE
*r_orig
;
2492 const struct real_format
*fmt
;
2498 round_for_format (fmt
, &r
);
2502 (*fmt
->encode
) (fmt
, buf
, &r
);
2507 /* Similar, but look up the format from MODE. */
2510 real_to_target (buf
, r
, mode
)
2512 const REAL_VALUE_TYPE
*r
;
2513 enum machine_mode mode
;
2515 const struct real_format
*fmt
;
2517 fmt
= real_format_for_mode
[mode
- QFmode
];
2521 return real_to_target_fmt (buf
, r
, fmt
);
2524 /* Read R from the given target format. Read the words of the result
2525 in target word order in BUF. There are always 32 bits in each
2526 long, no matter the size of the host long. */
2529 real_from_target_fmt (r
, buf
, fmt
)
2532 const struct real_format
*fmt
;
2534 (*fmt
->decode
) (fmt
, r
, buf
);
2537 /* Similar, but look up the format from MODE. */
2540 real_from_target (r
, buf
, mode
)
2543 enum machine_mode mode
;
2545 const struct real_format
*fmt
;
2547 fmt
= real_format_for_mode
[mode
- QFmode
];
2551 (*fmt
->decode
) (fmt
, r
, buf
);
2554 /* Return the number of bits in the significand for MODE. */
2555 /* ??? Legacy. Should get access to real_format directly. */
2558 significand_size (mode
)
2559 enum machine_mode mode
;
2561 const struct real_format
*fmt
;
2563 fmt
= real_format_for_mode
[mode
- QFmode
];
2567 return fmt
->p
* fmt
->log2_b
;
2570 /* Return a hash value for the given real value. */
2571 /* ??? The "unsigned int" return value is intended to be hashval_t,
2572 but I didn't want to pull hashtab.h into real.h. */
2576 const REAL_VALUE_TYPE
*r
;
2581 h
= r
->class | (r
->sign
<< 2);
2593 if (sizeof(unsigned long) > sizeof(unsigned int))
2594 for (i
= 0; i
< SIGSZ
; ++i
)
2596 unsigned long s
= r
->sig
[i
];
2597 h
^= s
^ (s
>> (HOST_BITS_PER_LONG
/ 2));
2600 for (i
= 0; i
< SIGSZ
; ++i
)
2611 /* IEEE single-precision format. */
2613 static void encode_ieee_single
PARAMS ((const struct real_format
*fmt
,
2614 long *, const REAL_VALUE_TYPE
*));
2615 static void decode_ieee_single
PARAMS ((const struct real_format
*,
2616 REAL_VALUE_TYPE
*, const long *));
2619 encode_ieee_single (fmt
, buf
, r
)
2620 const struct real_format
*fmt
;
2622 const REAL_VALUE_TYPE
*r
;
2624 unsigned long image
, sig
, exp
;
2625 bool denormal
= (r
->sig
[SIGSZ
-1] & SIG_MSB
) == 0;
2627 image
= r
->sign
<< 31;
2628 sig
= (r
->sig
[SIGSZ
-1] >> (HOST_BITS_PER_LONG
- 24)) & 0x7fffff;
2639 image
|= 0x7fffffff;
2647 if (!fmt
->qnan_msb_set
)
2648 image
^= 1 << 23 | 1 << 22;
2651 image
|= 0x7fffffff;
2655 /* Recall that IEEE numbers are interpreted as 1.F x 2**exp,
2656 whereas the intermediate representation is 0.F x 2**exp.
2657 Which means we're off by one. */
2661 exp
= r
->exp
+ 127 - 1;
2674 decode_ieee_single (fmt
, r
, buf
)
2675 const struct real_format
*fmt
;
2679 unsigned long image
= buf
[0] & 0xffffffff;
2680 bool sign
= (image
>> 31) & 1;
2681 int exp
= (image
>> 23) & 0xff;
2683 memset (r
, 0, sizeof (*r
));
2684 image
<<= HOST_BITS_PER_LONG
- 24;
2689 if (image
&& fmt
->has_denorm
)
2691 r
->class = rvc_normal
;
2694 r
->sig
[SIGSZ
-1] = image
<< 1;
2697 else if (fmt
->has_signed_zero
)
2700 else if (exp
== 255 && (fmt
->has_nans
|| fmt
->has_inf
))
2706 if (!fmt
->qnan_msb_set
)
2707 image
^= (SIG_MSB
>> 1 | SIG_MSB
>> 2);
2708 r
->sig
[SIGSZ
-1] = image
;
2718 r
->class = rvc_normal
;
2720 r
->exp
= exp
- 127 + 1;
2721 r
->sig
[SIGSZ
-1] = image
| SIG_MSB
;
2725 const struct real_format ieee_single_format
=
2743 /* IEEE double-precision format. */
2745 static void encode_ieee_double
PARAMS ((const struct real_format
*fmt
,
2746 long *, const REAL_VALUE_TYPE
*));
2747 static void decode_ieee_double
PARAMS ((const struct real_format
*,
2748 REAL_VALUE_TYPE
*, const long *));
2751 encode_ieee_double (fmt
, buf
, r
)
2752 const struct real_format
*fmt
;
2754 const REAL_VALUE_TYPE
*r
;
2756 unsigned long image_lo
, image_hi
, sig_lo
, sig_hi
, exp
;
2757 bool denormal
= (r
->sig
[SIGSZ
-1] & SIG_MSB
) == 0;
2759 image_hi
= r
->sign
<< 31;
2762 if (HOST_BITS_PER_LONG
== 64)
2764 sig_hi
= r
->sig
[SIGSZ
-1];
2765 sig_lo
= (sig_hi
>> (64 - 53)) & 0xffffffff;
2766 sig_hi
= (sig_hi
>> (64 - 53 + 1) >> 31) & 0xfffff;
2770 sig_hi
= r
->sig
[SIGSZ
-1];
2771 sig_lo
= r
->sig
[SIGSZ
-2];
2772 sig_lo
= (sig_hi
<< 21) | (sig_lo
>> 11);
2773 sig_hi
= (sig_hi
>> 11) & 0xfffff;
2783 image_hi
|= 2047 << 20;
2786 image_hi
|= 0x7fffffff;
2787 image_lo
= 0xffffffff;
2794 image_hi
|= 2047 << 20;
2796 if (!fmt
->qnan_msb_set
)
2797 image_hi
^= 1 << 19 | 1 << 18;
2802 image_hi
|= 0x7fffffff;
2803 image_lo
= 0xffffffff;
2808 /* Recall that IEEE numbers are interpreted as 1.F x 2**exp,
2809 whereas the intermediate representation is 0.F x 2**exp.
2810 Which means we're off by one. */
2814 exp
= r
->exp
+ 1023 - 1;
2815 image_hi
|= exp
<< 20;
2824 if (FLOAT_WORDS_BIG_ENDIAN
)
2825 buf
[0] = image_hi
, buf
[1] = image_lo
;
2827 buf
[0] = image_lo
, buf
[1] = image_hi
;
2831 decode_ieee_double (fmt
, r
, buf
)
2832 const struct real_format
*fmt
;
2836 unsigned long image_hi
, image_lo
;
2840 if (FLOAT_WORDS_BIG_ENDIAN
)
2841 image_hi
= buf
[0], image_lo
= buf
[1];
2843 image_lo
= buf
[0], image_hi
= buf
[1];
2844 image_lo
&= 0xffffffff;
2845 image_hi
&= 0xffffffff;
2847 sign
= (image_hi
>> 31) & 1;
2848 exp
= (image_hi
>> 20) & 0x7ff;
2850 memset (r
, 0, sizeof (*r
));
2852 image_hi
<<= 32 - 21;
2853 image_hi
|= image_lo
>> 21;
2854 image_hi
&= 0x7fffffff;
2855 image_lo
<<= 32 - 21;
2859 if ((image_hi
|| image_lo
) && fmt
->has_denorm
)
2861 r
->class = rvc_normal
;
2864 if (HOST_BITS_PER_LONG
== 32)
2866 image_hi
= (image_hi
<< 1) | (image_lo
>> 31);
2868 r
->sig
[SIGSZ
-1] = image_hi
;
2869 r
->sig
[SIGSZ
-2] = image_lo
;
2873 image_hi
= (image_hi
<< 31 << 2) | (image_lo
<< 1);
2874 r
->sig
[SIGSZ
-1] = image_hi
;
2878 else if (fmt
->has_signed_zero
)
2881 else if (exp
== 2047 && (fmt
->has_nans
|| fmt
->has_inf
))
2883 if (image_hi
|| image_lo
)
2887 if (HOST_BITS_PER_LONG
== 32)
2889 r
->sig
[SIGSZ
-1] = image_hi
;
2890 r
->sig
[SIGSZ
-2] = image_lo
;
2893 r
->sig
[SIGSZ
-1] = (image_hi
<< 31 << 1) | image_lo
;
2895 if (!fmt
->qnan_msb_set
)
2896 r
->sig
[SIGSZ
-1] ^= (SIG_MSB
>> 1 | SIG_MSB
>> 2);
2906 r
->class = rvc_normal
;
2908 r
->exp
= exp
- 1023 + 1;
2909 if (HOST_BITS_PER_LONG
== 32)
2911 r
->sig
[SIGSZ
-1] = image_hi
| SIG_MSB
;
2912 r
->sig
[SIGSZ
-2] = image_lo
;
2915 r
->sig
[SIGSZ
-1] = (image_hi
<< 31 << 1) | image_lo
| SIG_MSB
;
2919 const struct real_format ieee_double_format
=
2937 /* IEEE extended double precision format. This comes in three
2938 flavours: Intel's as a 12 byte image, Intel's as a 16 byte image,
2941 static void encode_ieee_extended
PARAMS ((const struct real_format
*fmt
,
2942 long *, const REAL_VALUE_TYPE
*));
2943 static void decode_ieee_extended
PARAMS ((const struct real_format
*,
2944 REAL_VALUE_TYPE
*, const long *));
2946 static void encode_ieee_extended_128
PARAMS ((const struct real_format
*fmt
,
2948 const REAL_VALUE_TYPE
*));
2949 static void decode_ieee_extended_128
PARAMS ((const struct real_format
*,
2954 encode_ieee_extended (fmt
, buf
, r
)
2955 const struct real_format
*fmt
;
2957 const REAL_VALUE_TYPE
*r
;
2959 unsigned long image_hi
, sig_hi
, sig_lo
;
2960 bool denormal
= (r
->sig
[SIGSZ
-1] & SIG_MSB
) == 0;
2962 image_hi
= r
->sign
<< 15;
2963 sig_hi
= sig_lo
= 0;
2975 /* Intel requires the explicit integer bit to be set, otherwise
2976 it considers the value a "pseudo-infinity". Motorola docs
2977 say it doesn't care. */
2978 sig_hi
= 0x80000000;
2983 sig_lo
= sig_hi
= 0xffffffff;
2991 if (HOST_BITS_PER_LONG
== 32)
2993 sig_hi
= r
->sig
[SIGSZ
-1];
2994 sig_lo
= r
->sig
[SIGSZ
-2];
2998 sig_lo
= r
->sig
[SIGSZ
-1];
2999 sig_hi
= sig_lo
>> 31 >> 1;
3000 sig_lo
&= 0xffffffff;
3002 if (!fmt
->qnan_msb_set
)
3003 sig_hi
^= 1 << 30 | 1 << 29;
3005 /* Intel requires the explicit integer bit to be set, otherwise
3006 it considers the value a "pseudo-nan". Motorola docs say it
3008 sig_hi
|= 0x80000000;
3013 sig_lo
= sig_hi
= 0xffffffff;
3021 /* Recall that IEEE numbers are interpreted as 1.F x 2**exp,
3022 whereas the intermediate representation is 0.F x 2**exp.
3023 Which means we're off by one.
3025 Except for Motorola, which consider exp=0 and explicit
3026 integer bit set to continue to be normalized. In theory
3027 this discrepancy has been taken care of by the difference
3028 in fmt->emin in round_for_format. */
3040 if (HOST_BITS_PER_LONG
== 32)
3042 sig_hi
= r
->sig
[SIGSZ
-1];
3043 sig_lo
= r
->sig
[SIGSZ
-2];
3047 sig_lo
= r
->sig
[SIGSZ
-1];
3048 sig_hi
= sig_lo
>> 31 >> 1;
3049 sig_lo
&= 0xffffffff;
3058 if (FLOAT_WORDS_BIG_ENDIAN
)
3059 buf
[0] = image_hi
<< 16, buf
[1] = sig_hi
, buf
[2] = sig_lo
;
3061 buf
[0] = sig_lo
, buf
[1] = sig_hi
, buf
[2] = image_hi
;
3065 encode_ieee_extended_128 (fmt
, buf
, r
)
3066 const struct real_format
*fmt
;
3068 const REAL_VALUE_TYPE
*r
;
3070 buf
[3 * !FLOAT_WORDS_BIG_ENDIAN
] = 0;
3071 encode_ieee_extended (fmt
, buf
+!!FLOAT_WORDS_BIG_ENDIAN
, r
);
3075 decode_ieee_extended (fmt
, r
, buf
)
3076 const struct real_format
*fmt
;
3080 unsigned long image_hi
, sig_hi
, sig_lo
;
3084 if (FLOAT_WORDS_BIG_ENDIAN
)
3085 image_hi
= buf
[0] >> 16, sig_hi
= buf
[1], sig_lo
= buf
[2];
3087 sig_lo
= buf
[0], sig_hi
= buf
[1], image_hi
= buf
[2];
3088 sig_lo
&= 0xffffffff;
3089 sig_hi
&= 0xffffffff;
3090 image_hi
&= 0xffffffff;
3092 sign
= (image_hi
>> 15) & 1;
3093 exp
= image_hi
& 0x7fff;
3095 memset (r
, 0, sizeof (*r
));
3099 if ((sig_hi
|| sig_lo
) && fmt
->has_denorm
)
3101 r
->class = rvc_normal
;
3104 /* When the IEEE format contains a hidden bit, we know that
3105 it's zero at this point, and so shift up the significand
3106 and decrease the exponent to match. In this case, Motorola
3107 defines the explicit integer bit to be valid, so we don't
3108 know whether the msb is set or not. */
3110 if (HOST_BITS_PER_LONG
== 32)
3112 r
->sig
[SIGSZ
-1] = sig_hi
;
3113 r
->sig
[SIGSZ
-2] = sig_lo
;
3116 r
->sig
[SIGSZ
-1] = (sig_hi
<< 31 << 1) | sig_lo
;
3120 else if (fmt
->has_signed_zero
)
3123 else if (exp
== 32767 && (fmt
->has_nans
|| fmt
->has_inf
))
3125 /* See above re "pseudo-infinities" and "pseudo-nans".
3126 Short summary is that the MSB will likely always be
3127 set, and that we don't care about it. */
3128 sig_hi
&= 0x7fffffff;
3130 if (sig_hi
|| sig_lo
)
3134 if (HOST_BITS_PER_LONG
== 32)
3136 r
->sig
[SIGSZ
-1] = sig_hi
;
3137 r
->sig
[SIGSZ
-2] = sig_lo
;
3140 r
->sig
[SIGSZ
-1] = (sig_hi
<< 31 << 1) | sig_lo
;
3142 if (!fmt
->qnan_msb_set
)
3143 r
->sig
[SIGSZ
-1] ^= (SIG_MSB
>> 1 | SIG_MSB
>> 2);
3153 r
->class = rvc_normal
;
3155 r
->exp
= exp
- 16383 + 1;
3156 if (HOST_BITS_PER_LONG
== 32)
3158 r
->sig
[SIGSZ
-1] = sig_hi
;
3159 r
->sig
[SIGSZ
-2] = sig_lo
;
3162 r
->sig
[SIGSZ
-1] = (sig_hi
<< 31 << 1) | sig_lo
;
3167 decode_ieee_extended_128 (fmt
, r
, buf
)
3168 const struct real_format
*fmt
;
3172 decode_ieee_extended (fmt
, r
, buf
+!!FLOAT_WORDS_BIG_ENDIAN
);
3175 const struct real_format ieee_extended_motorola_format
=
3177 encode_ieee_extended
,
3178 decode_ieee_extended
,
3192 const struct real_format ieee_extended_intel_96_format
=
3194 encode_ieee_extended
,
3195 decode_ieee_extended
,
3209 const struct real_format ieee_extended_intel_128_format
=
3211 encode_ieee_extended_128
,
3212 decode_ieee_extended_128
,
3227 /* IBM 128-bit extended precision format: a pair of IEEE double precision
3228 numbers whose sum is equal to the extended precision value. The number
3229 with greater magnitude is first. This format has the same magnitude
3230 range as an IEEE double precision value, but effectively 106 bits of
3231 significand precision. Infinity and NaN are represented by their IEEE
3232 double precision value stored in the first number, the second number is
3233 ignored. Zeroes, Infinities, and NaNs are set in both doubles
3234 due to precedent. */
3236 static void encode_ibm_extended
PARAMS ((const struct real_format
*fmt
,
3237 long *, const REAL_VALUE_TYPE
*));
3238 static void decode_ibm_extended
PARAMS ((const struct real_format
*,
3239 REAL_VALUE_TYPE
*, const long *));
3242 encode_ibm_extended (fmt
, buf
, r
)
3243 const struct real_format
*fmt ATTRIBUTE_UNUSED
;
3245 const REAL_VALUE_TYPE
*r
;
3247 REAL_VALUE_TYPE u
, v
;
3252 /* Both doubles have sign bit set. */
3253 buf
[0] = FLOAT_WORDS_BIG_ENDIAN
? r
->sign
<< 31 : 0;
3254 buf
[1] = FLOAT_WORDS_BIG_ENDIAN
? 0 : r
->sign
<< 31;
3261 /* Both doubles set to Inf / NaN. */
3262 encode_ieee_double (&ieee_double_format
, &buf
[0], r
);
3268 /* u = IEEE double precision portion of significand. */
3270 clear_significand_below (&u
, SIGNIFICAND_BITS
- 53);
3273 /* If the upper double is zero, we have a denormal double, so
3274 move it to the first double and leave the second as zero. */
3275 if (u
.class == rvc_zero
)
3283 /* v = remainder containing additional 53 bits of significand. */
3284 do_add (&v
, r
, &u
, 1);
3285 round_for_format (&ieee_double_format
, &v
);
3288 round_for_format (&ieee_double_format
, &u
);
3290 encode_ieee_double (&ieee_double_format
, &buf
[0], &u
);
3291 encode_ieee_double (&ieee_double_format
, &buf
[2], &v
);
3300 decode_ibm_extended (fmt
, r
, buf
)
3301 const struct real_format
*fmt ATTRIBUTE_UNUSED
;
3305 REAL_VALUE_TYPE u
, v
;
3307 decode_ieee_double (&ieee_double_format
, &u
, &buf
[0]);
3309 if (u
.class != rvc_zero
&& u
.class != rvc_inf
&& u
.class != rvc_nan
)
3311 decode_ieee_double (&ieee_double_format
, &v
, &buf
[2]);
3312 do_add (r
, &u
, &v
, 0);
3318 const struct real_format ibm_extended_format
=
3320 encode_ibm_extended
,
3321 decode_ibm_extended
,
3336 /* IEEE quad precision format. */
3338 static void encode_ieee_quad
PARAMS ((const struct real_format
*fmt
,
3339 long *, const REAL_VALUE_TYPE
*));
3340 static void decode_ieee_quad
PARAMS ((const struct real_format
*,
3341 REAL_VALUE_TYPE
*, const long *));
3344 encode_ieee_quad (fmt
, buf
, r
)
3345 const struct real_format
*fmt
;
3347 const REAL_VALUE_TYPE
*r
;
3349 unsigned long image3
, image2
, image1
, image0
, exp
;
3350 bool denormal
= (r
->sig
[SIGSZ
-1] & SIG_MSB
) == 0;
3353 image3
= r
->sign
<< 31;
3358 rshift_significand (&u
, r
, SIGNIFICAND_BITS
- 113);
3367 image3
|= 32767 << 16;
3370 image3
|= 0x7fffffff;
3371 image2
= 0xffffffff;
3372 image1
= 0xffffffff;
3373 image0
= 0xffffffff;
3380 image3
|= 32767 << 16;
3382 if (HOST_BITS_PER_LONG
== 32)
3387 image3
|= u
.sig
[3] & 0xffff;
3392 image1
= image0
>> 31 >> 1;
3394 image3
|= (image2
>> 31 >> 1) & 0xffff;
3395 image0
&= 0xffffffff;
3396 image2
&= 0xffffffff;
3399 if (!fmt
->qnan_msb_set
)
3400 image3
^= 1 << 15 | 1 << 14;
3404 image3
|= 0x7fffffff;
3405 image2
= 0xffffffff;
3406 image1
= 0xffffffff;
3407 image0
= 0xffffffff;
3412 /* Recall that IEEE numbers are interpreted as 1.F x 2**exp,
3413 whereas the intermediate representation is 0.F x 2**exp.
3414 Which means we're off by one. */
3418 exp
= r
->exp
+ 16383 - 1;
3419 image3
|= exp
<< 16;
3421 if (HOST_BITS_PER_LONG
== 32)
3426 image3
|= u
.sig
[3] & 0xffff;
3431 image1
= image0
>> 31 >> 1;
3433 image3
|= (image2
>> 31 >> 1) & 0xffff;
3434 image0
&= 0xffffffff;
3435 image2
&= 0xffffffff;
3443 if (FLOAT_WORDS_BIG_ENDIAN
)
3460 decode_ieee_quad (fmt
, r
, buf
)
3461 const struct real_format
*fmt
;
3465 unsigned long image3
, image2
, image1
, image0
;
3469 if (FLOAT_WORDS_BIG_ENDIAN
)
3483 image0
&= 0xffffffff;
3484 image1
&= 0xffffffff;
3485 image2
&= 0xffffffff;
3487 sign
= (image3
>> 31) & 1;
3488 exp
= (image3
>> 16) & 0x7fff;
3491 memset (r
, 0, sizeof (*r
));
3495 if ((image3
| image2
| image1
| image0
) && fmt
->has_denorm
)
3497 r
->class = rvc_normal
;
3500 r
->exp
= -16382 + (SIGNIFICAND_BITS
- 112);
3501 if (HOST_BITS_PER_LONG
== 32)
3510 r
->sig
[0] = (image1
<< 31 << 1) | image0
;
3511 r
->sig
[1] = (image3
<< 31 << 1) | image2
;
3516 else if (fmt
->has_signed_zero
)
3519 else if (exp
== 32767 && (fmt
->has_nans
|| fmt
->has_inf
))
3521 if (image3
| image2
| image1
| image0
)
3526 if (HOST_BITS_PER_LONG
== 32)
3535 r
->sig
[0] = (image1
<< 31 << 1) | image0
;
3536 r
->sig
[1] = (image3
<< 31 << 1) | image2
;
3538 lshift_significand (r
, r
, SIGNIFICAND_BITS
- 113);
3540 if (!fmt
->qnan_msb_set
)
3541 r
->sig
[SIGSZ
-1] ^= (SIG_MSB
>> 1 | SIG_MSB
>> 2);
3551 r
->class = rvc_normal
;
3553 r
->exp
= exp
- 16383 + 1;
3555 if (HOST_BITS_PER_LONG
== 32)
3564 r
->sig
[0] = (image1
<< 31 << 1) | image0
;
3565 r
->sig
[1] = (image3
<< 31 << 1) | image2
;
3567 lshift_significand (r
, r
, SIGNIFICAND_BITS
- 113);
3568 r
->sig
[SIGSZ
-1] |= SIG_MSB
;
3572 const struct real_format ieee_quad_format
=
3589 /* Descriptions of VAX floating point formats can be found beginning at
3591 http://www.openvms.compaq.com:8000/73final/4515/4515pro_013.html#f_floating_point_format
3593 The thing to remember is that they're almost IEEE, except for word
3594 order, exponent bias, and the lack of infinities, nans, and denormals.
3596 We don't implement the H_floating format here, simply because neither
3597 the VAX or Alpha ports use it. */
3599 static void encode_vax_f
PARAMS ((const struct real_format
*fmt
,
3600 long *, const REAL_VALUE_TYPE
*));
3601 static void decode_vax_f
PARAMS ((const struct real_format
*,
3602 REAL_VALUE_TYPE
*, const long *));
3603 static void encode_vax_d
PARAMS ((const struct real_format
*fmt
,
3604 long *, const REAL_VALUE_TYPE
*));
3605 static void decode_vax_d
PARAMS ((const struct real_format
*,
3606 REAL_VALUE_TYPE
*, const long *));
3607 static void encode_vax_g
PARAMS ((const struct real_format
*fmt
,
3608 long *, const REAL_VALUE_TYPE
*));
3609 static void decode_vax_g
PARAMS ((const struct real_format
*,
3610 REAL_VALUE_TYPE
*, const long *));
3613 encode_vax_f (fmt
, buf
, r
)
3614 const struct real_format
*fmt ATTRIBUTE_UNUSED
;
3616 const REAL_VALUE_TYPE
*r
;
3618 unsigned long sign
, exp
, sig
, image
;
3620 sign
= r
->sign
<< 15;
3630 image
= 0xffff7fff | sign
;
3634 sig
= (r
->sig
[SIGSZ
-1] >> (HOST_BITS_PER_LONG
- 24)) & 0x7fffff;
3637 image
= (sig
<< 16) & 0xffff0000;
3651 decode_vax_f (fmt
, r
, buf
)
3652 const struct real_format
*fmt ATTRIBUTE_UNUSED
;
3656 unsigned long image
= buf
[0] & 0xffffffff;
3657 int exp
= (image
>> 7) & 0xff;
3659 memset (r
, 0, sizeof (*r
));
3663 r
->class = rvc_normal
;
3664 r
->sign
= (image
>> 15) & 1;
3667 image
= ((image
& 0x7f) << 16) | ((image
>> 16) & 0xffff);
3668 r
->sig
[SIGSZ
-1] = (image
<< (HOST_BITS_PER_LONG
- 24)) | SIG_MSB
;
3673 encode_vax_d (fmt
, buf
, r
)
3674 const struct real_format
*fmt ATTRIBUTE_UNUSED
;
3676 const REAL_VALUE_TYPE
*r
;
3678 unsigned long image0
, image1
, sign
= r
->sign
<< 15;
3683 image0
= image1
= 0;
3688 image0
= 0xffff7fff | sign
;
3689 image1
= 0xffffffff;
3693 /* Extract the significand into straight hi:lo. */
3694 if (HOST_BITS_PER_LONG
== 64)
3696 image0
= r
->sig
[SIGSZ
-1];
3697 image1
= (image0
>> (64 - 56)) & 0xffffffff;
3698 image0
= (image0
>> (64 - 56 + 1) >> 31) & 0x7fffff;
3702 image0
= r
->sig
[SIGSZ
-1];
3703 image1
= r
->sig
[SIGSZ
-2];
3704 image1
= (image0
<< 24) | (image1
>> 8);
3705 image0
= (image0
>> 8) & 0xffffff;
3708 /* Rearrange the half-words of the significand to match the
3710 image0
= ((image0
<< 16) | (image0
>> 16)) & 0xffff007f;
3711 image1
= ((image1
<< 16) | (image1
>> 16)) & 0xffffffff;
3713 /* Add the sign and exponent. */
3715 image0
|= (r
->exp
+ 128) << 7;
3722 if (FLOAT_WORDS_BIG_ENDIAN
)
3723 buf
[0] = image1
, buf
[1] = image0
;
3725 buf
[0] = image0
, buf
[1] = image1
;
3729 decode_vax_d (fmt
, r
, buf
)
3730 const struct real_format
*fmt ATTRIBUTE_UNUSED
;
3734 unsigned long image0
, image1
;
3737 if (FLOAT_WORDS_BIG_ENDIAN
)
3738 image1
= buf
[0], image0
= buf
[1];
3740 image0
= buf
[0], image1
= buf
[1];
3741 image0
&= 0xffffffff;
3742 image1
&= 0xffffffff;
3744 exp
= (image0
>> 7) & 0x7f;
3746 memset (r
, 0, sizeof (*r
));
3750 r
->class = rvc_normal
;
3751 r
->sign
= (image0
>> 15) & 1;
3754 /* Rearrange the half-words of the external format into
3755 proper ascending order. */
3756 image0
= ((image0
& 0x7f) << 16) | ((image0
>> 16) & 0xffff);
3757 image1
= ((image1
& 0xffff) << 16) | ((image1
>> 16) & 0xffff);
3759 if (HOST_BITS_PER_LONG
== 64)
3761 image0
= (image0
<< 31 << 1) | image1
;
3764 r
->sig
[SIGSZ
-1] = image0
;
3768 r
->sig
[SIGSZ
-1] = image0
;
3769 r
->sig
[SIGSZ
-2] = image1
;
3770 lshift_significand (r
, r
, 2*HOST_BITS_PER_LONG
- 56);
3771 r
->sig
[SIGSZ
-1] |= SIG_MSB
;
3777 encode_vax_g (fmt
, buf
, r
)
3778 const struct real_format
*fmt ATTRIBUTE_UNUSED
;
3780 const REAL_VALUE_TYPE
*r
;
3782 unsigned long image0
, image1
, sign
= r
->sign
<< 15;
3787 image0
= image1
= 0;
3792 image0
= 0xffff7fff | sign
;
3793 image1
= 0xffffffff;
3797 /* Extract the significand into straight hi:lo. */
3798 if (HOST_BITS_PER_LONG
== 64)
3800 image0
= r
->sig
[SIGSZ
-1];
3801 image1
= (image0
>> (64 - 53)) & 0xffffffff;
3802 image0
= (image0
>> (64 - 53 + 1) >> 31) & 0xfffff;
3806 image0
= r
->sig
[SIGSZ
-1];
3807 image1
= r
->sig
[SIGSZ
-2];
3808 image1
= (image0
<< 21) | (image1
>> 11);
3809 image0
= (image0
>> 11) & 0xfffff;
3812 /* Rearrange the half-words of the significand to match the
3814 image0
= ((image0
<< 16) | (image0
>> 16)) & 0xffff000f;
3815 image1
= ((image1
<< 16) | (image1
>> 16)) & 0xffffffff;
3817 /* Add the sign and exponent. */
3819 image0
|= (r
->exp
+ 1024) << 4;
3826 if (FLOAT_WORDS_BIG_ENDIAN
)
3827 buf
[0] = image1
, buf
[1] = image0
;
3829 buf
[0] = image0
, buf
[1] = image1
;
3833 decode_vax_g (fmt
, r
, buf
)
3834 const struct real_format
*fmt ATTRIBUTE_UNUSED
;
3838 unsigned long image0
, image1
;
3841 if (FLOAT_WORDS_BIG_ENDIAN
)
3842 image1
= buf
[0], image0
= buf
[1];
3844 image0
= buf
[0], image1
= buf
[1];
3845 image0
&= 0xffffffff;
3846 image1
&= 0xffffffff;
3848 exp
= (image0
>> 4) & 0x7ff;
3850 memset (r
, 0, sizeof (*r
));
3854 r
->class = rvc_normal
;
3855 r
->sign
= (image0
>> 15) & 1;
3856 r
->exp
= exp
- 1024;
3858 /* Rearrange the half-words of the external format into
3859 proper ascending order. */
3860 image0
= ((image0
& 0xf) << 16) | ((image0
>> 16) & 0xffff);
3861 image1
= ((image1
& 0xffff) << 16) | ((image1
>> 16) & 0xffff);
3863 if (HOST_BITS_PER_LONG
== 64)
3865 image0
= (image0
<< 31 << 1) | image1
;
3868 r
->sig
[SIGSZ
-1] = image0
;
3872 r
->sig
[SIGSZ
-1] = image0
;
3873 r
->sig
[SIGSZ
-2] = image1
;
3874 lshift_significand (r
, r
, 64 - 53);
3875 r
->sig
[SIGSZ
-1] |= SIG_MSB
;
3880 const struct real_format vax_f_format
=
3897 const struct real_format vax_d_format
=
3914 const struct real_format vax_g_format
=
3931 /* A good reference for these can be found in chapter 9 of
3932 "ESA/390 Principles of Operation", IBM document number SA22-7201-01.
3933 An on-line version can be found here:
3935 http://publibz.boulder.ibm.com/cgi-bin/bookmgr_OS390/BOOKS/DZ9AR001/9.1?DT=19930923083613
3938 static void encode_i370_single
PARAMS ((const struct real_format
*fmt
,
3939 long *, const REAL_VALUE_TYPE
*));
3940 static void decode_i370_single
PARAMS ((const struct real_format
*,
3941 REAL_VALUE_TYPE
*, const long *));
3942 static void encode_i370_double
PARAMS ((const struct real_format
*fmt
,
3943 long *, const REAL_VALUE_TYPE
*));
3944 static void decode_i370_double
PARAMS ((const struct real_format
*,
3945 REAL_VALUE_TYPE
*, const long *));
3948 encode_i370_single (fmt
, buf
, r
)
3949 const struct real_format
*fmt ATTRIBUTE_UNUSED
;
3951 const REAL_VALUE_TYPE
*r
;
3953 unsigned long sign
, exp
, sig
, image
;
3955 sign
= r
->sign
<< 31;
3965 image
= 0x7fffffff | sign
;
3969 sig
= (r
->sig
[SIGSZ
-1] >> (HOST_BITS_PER_LONG
- 24)) & 0xffffff;
3970 exp
= ((r
->exp
/ 4) + 64) << 24;
3971 image
= sign
| exp
| sig
;
3982 decode_i370_single (fmt
, r
, buf
)
3983 const struct real_format
*fmt ATTRIBUTE_UNUSED
;
3987 unsigned long sign
, sig
, image
= buf
[0];
3990 sign
= (image
>> 31) & 1;
3991 exp
= (image
>> 24) & 0x7f;
3992 sig
= image
& 0xffffff;
3994 memset (r
, 0, sizeof (*r
));
3998 r
->class = rvc_normal
;
4000 r
->exp
= (exp
- 64) * 4;
4001 r
->sig
[SIGSZ
-1] = sig
<< (HOST_BITS_PER_LONG
- 24);
4007 encode_i370_double (fmt
, buf
, r
)
4008 const struct real_format
*fmt ATTRIBUTE_UNUSED
;
4010 const REAL_VALUE_TYPE
*r
;
4012 unsigned long sign
, exp
, image_hi
, image_lo
;
4014 sign
= r
->sign
<< 31;
4019 image_hi
= image_lo
= 0;
4024 image_hi
= 0x7fffffff | sign
;
4025 image_lo
= 0xffffffff;
4029 if (HOST_BITS_PER_LONG
== 64)
4031 image_hi
= r
->sig
[SIGSZ
-1];
4032 image_lo
= (image_hi
>> (64 - 56)) & 0xffffffff;
4033 image_hi
= (image_hi
>> (64 - 56 + 1) >> 31) & 0xffffff;
4037 image_hi
= r
->sig
[SIGSZ
-1];
4038 image_lo
= r
->sig
[SIGSZ
-2];
4039 image_lo
= (image_lo
>> 8) | (image_hi
<< 24);
4043 exp
= ((r
->exp
/ 4) + 64) << 24;
4044 image_hi
|= sign
| exp
;
4051 if (FLOAT_WORDS_BIG_ENDIAN
)
4052 buf
[0] = image_hi
, buf
[1] = image_lo
;
4054 buf
[0] = image_lo
, buf
[1] = image_hi
;
4058 decode_i370_double (fmt
, r
, buf
)
4059 const struct real_format
*fmt ATTRIBUTE_UNUSED
;
4063 unsigned long sign
, image_hi
, image_lo
;
4066 if (FLOAT_WORDS_BIG_ENDIAN
)
4067 image_hi
= buf
[0], image_lo
= buf
[1];
4069 image_lo
= buf
[0], image_hi
= buf
[1];
4071 sign
= (image_hi
>> 31) & 1;
4072 exp
= (image_hi
>> 24) & 0x7f;
4073 image_hi
&= 0xffffff;
4074 image_lo
&= 0xffffffff;
4076 memset (r
, 0, sizeof (*r
));
4078 if (exp
|| image_hi
|| image_lo
)
4080 r
->class = rvc_normal
;
4082 r
->exp
= (exp
- 64) * 4 + (SIGNIFICAND_BITS
- 56);
4084 if (HOST_BITS_PER_LONG
== 32)
4086 r
->sig
[0] = image_lo
;
4087 r
->sig
[1] = image_hi
;
4090 r
->sig
[0] = image_lo
| (image_hi
<< 31 << 1);
4096 const struct real_format i370_single_format
=
4108 false, /* ??? The encoding does allow for "unnormals". */
4109 false, /* ??? The encoding does allow for "unnormals". */
4113 const struct real_format i370_double_format
=
4125 false, /* ??? The encoding does allow for "unnormals". */
4126 false, /* ??? The encoding does allow for "unnormals". */
4130 /* The "twos-complement" c4x format is officially defined as
4134 This is rather misleading. One must remember that F is signed.
4135 A better description would be
4137 x = -1**s * ((s + 1 + .f) * 2**e
4139 So if we have a (4 bit) fraction of .1000 with a sign bit of 1,
4140 that's -1 * (1+1+(-.5)) == -1.5. I think.
4142 The constructions here are taken from Tables 5-1 and 5-2 of the
4143 TMS320C4x User's Guide wherein step-by-step instructions for
4144 conversion from IEEE are presented. That's close enough to our
4145 internal representation so as to make things easy.
4147 See http://www-s.ti.com/sc/psheets/spru063c/spru063c.pdf */
4149 static void encode_c4x_single
PARAMS ((const struct real_format
*fmt
,
4150 long *, const REAL_VALUE_TYPE
*));
4151 static void decode_c4x_single
PARAMS ((const struct real_format
*,
4152 REAL_VALUE_TYPE
*, const long *));
4153 static void encode_c4x_extended
PARAMS ((const struct real_format
*fmt
,
4154 long *, const REAL_VALUE_TYPE
*));
4155 static void decode_c4x_extended
PARAMS ((const struct real_format
*,
4156 REAL_VALUE_TYPE
*, const long *));
4159 encode_c4x_single (fmt
, buf
, r
)
4160 const struct real_format
*fmt ATTRIBUTE_UNUSED
;
4162 const REAL_VALUE_TYPE
*r
;
4164 unsigned long image
, exp
, sig
;
4176 sig
= 0x800000 - r
->sign
;
4181 sig
= (r
->sig
[SIGSZ
-1] >> (HOST_BITS_PER_LONG
- 24)) & 0x7fffff;
4196 image
= ((exp
& 0xff) << 24) | (sig
& 0xffffff);
4201 decode_c4x_single (fmt
, r
, buf
)
4202 const struct real_format
*fmt ATTRIBUTE_UNUSED
;
4206 unsigned long image
= buf
[0];
4210 exp
= (((image
>> 24) & 0xff) ^ 0x80) - 0x80;
4211 sf
= ((image
& 0xffffff) ^ 0x800000) - 0x800000;
4213 memset (r
, 0, sizeof (*r
));
4217 r
->class = rvc_normal
;
4219 sig
= sf
& 0x7fffff;
4228 sig
= (sig
<< (HOST_BITS_PER_LONG
- 24)) | SIG_MSB
;
4231 r
->sig
[SIGSZ
-1] = sig
;
4236 encode_c4x_extended (fmt
, buf
, r
)
4237 const struct real_format
*fmt ATTRIBUTE_UNUSED
;
4239 const REAL_VALUE_TYPE
*r
;
4241 unsigned long exp
, sig
;
4253 sig
= 0x80000000 - r
->sign
;
4259 sig
= r
->sig
[SIGSZ
-1];
4260 if (HOST_BITS_PER_LONG
== 64)
4261 sig
= sig
>> 1 >> 31;
4278 exp
= (exp
& 0xff) << 24;
4281 if (FLOAT_WORDS_BIG_ENDIAN
)
4282 buf
[0] = exp
, buf
[1] = sig
;
4284 buf
[0] = sig
, buf
[0] = exp
;
4288 decode_c4x_extended (fmt
, r
, buf
)
4289 const struct real_format
*fmt ATTRIBUTE_UNUSED
;
4296 if (FLOAT_WORDS_BIG_ENDIAN
)
4297 exp
= buf
[0], sf
= buf
[1];
4299 sf
= buf
[0], exp
= buf
[1];
4301 exp
= (((exp
>> 24) & 0xff) & 0x80) - 0x80;
4302 sf
= ((sf
& 0xffffffff) ^ 0x80000000) - 0x80000000;
4304 memset (r
, 0, sizeof (*r
));
4308 r
->class = rvc_normal
;
4310 sig
= sf
& 0x7fffffff;
4319 if (HOST_BITS_PER_LONG
== 64)
4320 sig
= sig
<< 1 << 31;
4324 r
->sig
[SIGSZ
-1] = sig
;
4328 const struct real_format c4x_single_format
=
4345 const struct real_format c4x_extended_format
=
4347 encode_c4x_extended
,
4348 decode_c4x_extended
,
4363 /* A synthetic "format" for internal arithmetic. It's the size of the
4364 internal significand minus the two bits needed for proper rounding.
4365 The encode and decode routines exist only to satisfy our paranoia
4368 static void encode_internal
PARAMS ((const struct real_format
*fmt
,
4369 long *, const REAL_VALUE_TYPE
*));
4370 static void decode_internal
PARAMS ((const struct real_format
*,
4371 REAL_VALUE_TYPE
*, const long *));
4374 encode_internal (fmt
, buf
, r
)
4375 const struct real_format
*fmt ATTRIBUTE_UNUSED
;
4377 const REAL_VALUE_TYPE
*r
;
4379 memcpy (buf
, r
, sizeof (*r
));
4383 decode_internal (fmt
, r
, buf
)
4384 const struct real_format
*fmt ATTRIBUTE_UNUSED
;
4388 memcpy (r
, buf
, sizeof (*r
));
4391 const struct real_format real_internal_format
=
4397 SIGNIFICAND_BITS
- 2,
4408 /* Set up default mode to format mapping for IEEE. Everyone else has
4409 to set these values in OVERRIDE_OPTIONS. */
4411 const struct real_format
*real_format_for_mode
[TFmode
- QFmode
+ 1] =
4416 &ieee_single_format
, /* SFmode */
4417 &ieee_double_format
, /* DFmode */
4419 /* We explicitly don't handle XFmode. There are two formats,
4420 pretty much equally common. Choose one in OVERRIDE_OPTIONS. */
4422 &ieee_quad_format
/* TFmode */
4426 /* Calculate the square root of X in mode MODE, and store the result
4427 in R. Return TRUE if the operation does not raise an exception.
4428 For details see "High Precision Division and Square Root",
4429 Alan H. Karp and Peter Markstein, HP Lab Report 93-93-42, June
4430 1993. http://www.hpl.hp.com/techreports/93/HPL-93-42.pdf. */
4433 real_sqrt (r
, mode
, x
)
4435 enum machine_mode mode
;
4436 const REAL_VALUE_TYPE
*x
;
4438 static REAL_VALUE_TYPE halfthree
;
4439 static REAL_VALUE_TYPE half
;
4440 static bool init
= false;
4441 REAL_VALUE_TYPE h
, t
, i
;
4444 /* sqrt(-0.0) is -0.0. */
4445 if (real_isnegzero (x
))
4451 /* Negative arguments return NaN. */
4454 /* Mode is ignored for canonical NaN. */
4455 real_nan (r
, "", 1, SFmode
);
4459 /* Infinity and NaN return themselves. */
4460 if (real_isinf (x
) || real_isnan (x
))
4468 real_arithmetic (&half
, RDIV_EXPR
, &dconst1
, &dconst2
);
4469 real_arithmetic (&halfthree
, PLUS_EXPR
, &dconst1
, &half
);
4473 /* Initial guess for reciprocal sqrt, i. */
4474 exp
= real_exponent (x
);
4475 real_ldexp (&i
, &dconst1
, -exp
/2);
4477 /* Newton's iteration for reciprocal sqrt, i. */
4478 for (iter
= 0; iter
< 16; iter
++)
4480 /* i(n+1) = i(n) * (1.5 - 0.5*i(n)*i(n)*x). */
4481 real_arithmetic (&t
, MULT_EXPR
, x
, &i
);
4482 real_arithmetic (&h
, MULT_EXPR
, &t
, &i
);
4483 real_arithmetic (&t
, MULT_EXPR
, &h
, &half
);
4484 real_arithmetic (&h
, MINUS_EXPR
, &halfthree
, &t
);
4485 real_arithmetic (&t
, MULT_EXPR
, &i
, &h
);
4487 /* Check for early convergence. */
4488 if (iter
>= 6 && real_identical (&i
, &t
))
4491 /* ??? Unroll loop to avoid copying. */
4495 /* Final iteration: r = i*x + 0.5*i*x*(1.0 - i*(i*x)). */
4496 real_arithmetic (&t
, MULT_EXPR
, x
, &i
);
4497 real_arithmetic (&h
, MULT_EXPR
, &t
, &i
);
4498 real_arithmetic (&i
, MINUS_EXPR
, &dconst1
, &h
);
4499 real_arithmetic (&h
, MULT_EXPR
, &t
, &i
);
4500 real_arithmetic (&i
, MULT_EXPR
, &half
, &h
);
4501 real_arithmetic (&h
, PLUS_EXPR
, &t
, &i
);
4503 /* ??? We need a Tuckerman test to get the last bit. */
4505 real_convert (r
, mode
, &h
);