1 /* real.c - software floating point emulation.
2 Copyright (C) 1993, 1994, 1995, 1996, 1997, 1998, 1999,
3 2000, 2002, 2003, 2004 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 C99 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 the largest
53 supported target floating-point type by at least 2 bits. This gives
54 us proper rounding when we truncate to the target type. In addition,
55 E must be large enough to hold the smallest supported denormal number
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 27.
62 Note that the decimal string conversion routines are sensitive to
63 rounding errors. 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 (REAL_VALUE_TYPE *, int);
83 static void get_canonical_qnan (REAL_VALUE_TYPE *, int);
84 static void get_canonical_snan (REAL_VALUE_TYPE *, int);
85 static void get_inf (REAL_VALUE_TYPE *, int);
86 static bool sticky_rshift_significand (REAL_VALUE_TYPE *,
87 const REAL_VALUE_TYPE *, unsigned int);
88 static void rshift_significand (REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *,
90 static void lshift_significand (REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *,
92 static void lshift_significand_1 (REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *);
93 static bool add_significands (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *,
94 const REAL_VALUE_TYPE *);
95 static bool sub_significands (REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *,
96 const REAL_VALUE_TYPE *, int);
97 static void neg_significand (REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *);
98 static int cmp_significands (const REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *);
99 static int cmp_significand_0 (const REAL_VALUE_TYPE *);
100 static void set_significand_bit (REAL_VALUE_TYPE *, unsigned int);
101 static void clear_significand_bit (REAL_VALUE_TYPE *, unsigned int);
102 static bool test_significand_bit (REAL_VALUE_TYPE *, unsigned int);
103 static void clear_significand_below (REAL_VALUE_TYPE *, unsigned int);
104 static bool div_significands (REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *,
105 const REAL_VALUE_TYPE *);
106 static void normalize (REAL_VALUE_TYPE *);
108 static bool do_add (REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *,
109 const REAL_VALUE_TYPE *, int);
110 static bool do_multiply (REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *,
111 const REAL_VALUE_TYPE *);
112 static bool do_divide (REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *,
113 const REAL_VALUE_TYPE *);
114 static int do_compare (const REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *, int);
115 static void do_fix_trunc (REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *);
117 static unsigned long rtd_divmod (REAL_VALUE_TYPE *, REAL_VALUE_TYPE *);
119 static const REAL_VALUE_TYPE * ten_to_ptwo (int);
120 static const REAL_VALUE_TYPE * ten_to_mptwo (int);
121 static const REAL_VALUE_TYPE * real_digit (int);
122 static void times_pten (REAL_VALUE_TYPE *, int);
124 static void round_for_format (const struct real_format *, REAL_VALUE_TYPE *);
126 /* Initialize R with a positive zero. */
129 get_zero (REAL_VALUE_TYPE *r, int sign)
131 memset (r, 0, sizeof (*r));
135 /* Initialize R with the canonical quiet NaN. */
138 get_canonical_qnan (REAL_VALUE_TYPE *r, int sign)
140 memset (r, 0, sizeof (*r));
147 get_canonical_snan (REAL_VALUE_TYPE *r, int sign)
149 memset (r, 0, sizeof (*r));
157 get_inf (REAL_VALUE_TYPE *r, int sign)
159 memset (r, 0, sizeof (*r));
165 /* Right-shift the significand of A by N bits; put the result in the
166 significand of R. If any one bits are shifted out, return true. */
169 sticky_rshift_significand (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *a,
172 unsigned long sticky = 0;
173 unsigned int i, ofs = 0;
175 if (n >= HOST_BITS_PER_LONG)
177 for (i = 0, ofs = n / HOST_BITS_PER_LONG; i < ofs; ++i)
179 n &= HOST_BITS_PER_LONG - 1;
184 sticky |= a->sig[ofs] & (((unsigned long)1 << n) - 1);
185 for (i = 0; i < SIGSZ; ++i)
188 = (((ofs + i >= SIGSZ ? 0 : a->sig[ofs + i]) >> n)
189 | ((ofs + i + 1 >= SIGSZ ? 0 : a->sig[ofs + i + 1])
190 << (HOST_BITS_PER_LONG - n)));
195 for (i = 0; ofs + i < SIGSZ; ++i)
196 r->sig[i] = a->sig[ofs + i];
197 for (; i < SIGSZ; ++i)
204 /* Right-shift the significand of A by N bits; put the result in the
208 rshift_significand (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *a,
211 unsigned int i, ofs = n / HOST_BITS_PER_LONG;
213 n &= HOST_BITS_PER_LONG - 1;
216 for (i = 0; i < SIGSZ; ++i)
219 = (((ofs + i >= SIGSZ ? 0 : a->sig[ofs + i]) >> n)
220 | ((ofs + i + 1 >= SIGSZ ? 0 : a->sig[ofs + i + 1])
221 << (HOST_BITS_PER_LONG - n)));
226 for (i = 0; ofs + i < SIGSZ; ++i)
227 r->sig[i] = a->sig[ofs + i];
228 for (; i < SIGSZ; ++i)
233 /* Left-shift the significand of A by N bits; put the result in the
237 lshift_significand (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *a,
240 unsigned int i, ofs = n / HOST_BITS_PER_LONG;
242 n &= HOST_BITS_PER_LONG - 1;
245 for (i = 0; ofs + i < SIGSZ; ++i)
246 r->sig[SIGSZ-1-i] = a->sig[SIGSZ-1-i-ofs];
247 for (; i < SIGSZ; ++i)
248 r->sig[SIGSZ-1-i] = 0;
251 for (i = 0; i < SIGSZ; ++i)
254 = (((ofs + i >= SIGSZ ? 0 : a->sig[SIGSZ-1-i-ofs]) << n)
255 | ((ofs + i + 1 >= SIGSZ ? 0 : a->sig[SIGSZ-1-i-ofs-1])
256 >> (HOST_BITS_PER_LONG - n)));
260 /* Likewise, but N is specialized to 1. */
263 lshift_significand_1 (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *a)
267 for (i = SIGSZ - 1; i > 0; --i)
268 r->sig[i] = (a->sig[i] << 1) | (a->sig[i-1] >> (HOST_BITS_PER_LONG - 1));
269 r->sig[0] = a->sig[0] << 1;
272 /* Add the significands of A and B, placing the result in R. Return
273 true if there was carry out of the most significant word. */
276 add_significands (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *a,
277 const REAL_VALUE_TYPE *b)
282 for (i = 0; i < SIGSZ; ++i)
284 unsigned long ai = a->sig[i];
285 unsigned long ri = ai + b->sig[i];
301 /* Subtract the significands of A and B, placing the result in R. CARRY is
302 true if there's a borrow incoming to the least significant word.
303 Return true if there was borrow out of the most significant word. */
306 sub_significands (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *a,
307 const REAL_VALUE_TYPE *b, int carry)
311 for (i = 0; i < SIGSZ; ++i)
313 unsigned long ai = a->sig[i];
314 unsigned long ri = ai - b->sig[i];
330 /* Negate the significand A, placing the result in R. */
333 neg_significand (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *a)
338 for (i = 0; i < SIGSZ; ++i)
340 unsigned long ri, ai = a->sig[i];
359 /* Compare significands. Return tri-state vs zero. */
362 cmp_significands (const REAL_VALUE_TYPE *a, const REAL_VALUE_TYPE *b)
366 for (i = SIGSZ - 1; i >= 0; --i)
368 unsigned long ai = a->sig[i];
369 unsigned long bi = b->sig[i];
380 /* Return true if A is nonzero. */
383 cmp_significand_0 (const REAL_VALUE_TYPE *a)
387 for (i = SIGSZ - 1; i >= 0; --i)
394 /* Set bit N of the significand of R. */
397 set_significand_bit (REAL_VALUE_TYPE *r, unsigned int n)
399 r->sig[n / HOST_BITS_PER_LONG]
400 |= (unsigned long)1 << (n % HOST_BITS_PER_LONG);
403 /* Clear bit N of the significand of R. */
406 clear_significand_bit (REAL_VALUE_TYPE *r, unsigned int n)
408 r->sig[n / HOST_BITS_PER_LONG]
409 &= ~((unsigned long)1 << (n % HOST_BITS_PER_LONG));
412 /* Test bit N of the significand of R. */
415 test_significand_bit (REAL_VALUE_TYPE *r, unsigned int n)
417 /* ??? Compiler bug here if we return this expression directly.
418 The conversion to bool strips the "&1" and we wind up testing
419 e.g. 2 != 0 -> true. Seen in gcc version 3.2 20020520. */
420 int t = (r->sig[n / HOST_BITS_PER_LONG] >> (n % HOST_BITS_PER_LONG)) & 1;
424 /* Clear bits 0..N-1 of the significand of R. */
427 clear_significand_below (REAL_VALUE_TYPE *r, unsigned int n)
429 int i, w = n / HOST_BITS_PER_LONG;
431 for (i = 0; i < w; ++i)
434 r->sig[w] &= ~(((unsigned long)1 << (n % HOST_BITS_PER_LONG)) - 1);
437 /* Divide the significands of A and B, placing the result in R. Return
438 true if the division was inexact. */
441 div_significands (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *a,
442 const REAL_VALUE_TYPE *b)
445 int i, bit = SIGNIFICAND_BITS - 1;
446 unsigned long msb, inexact;
449 memset (r->sig, 0, sizeof (r->sig));
455 msb = u.sig[SIGSZ-1] & SIG_MSB;
456 lshift_significand_1 (&u, &u);
458 if (msb || cmp_significands (&u, b) >= 0)
460 sub_significands (&u, &u, b, 0);
461 set_significand_bit (r, bit);
466 for (i = 0, inexact = 0; i < SIGSZ; i++)
472 /* Adjust the exponent and significand of R such that the most
473 significant bit is set. We underflow to zero and overflow to
474 infinity here, without denormals. (The intermediate representation
475 exponent is large enough to handle target denormals normalized.) */
478 normalize (REAL_VALUE_TYPE *r)
483 /* Find the first word that is nonzero. */
484 for (i = SIGSZ - 1; i >= 0; i--)
486 shift += HOST_BITS_PER_LONG;
490 /* Zero significand flushes to zero. */
498 /* Find the first bit that is nonzero. */
500 if (r->sig[i] & ((unsigned long)1 << (HOST_BITS_PER_LONG - 1 - j)))
506 exp = r->exp - shift;
508 get_inf (r, r->sign);
509 else if (exp < -MAX_EXP)
510 get_zero (r, r->sign);
514 lshift_significand (r, r, shift);
519 /* Calculate R = A + (SUBTRACT_P ? -B : B). Return true if the
520 result may be inexact due to a loss of precision. */
523 do_add (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *a,
524 const REAL_VALUE_TYPE *b, int subtract_p)
528 bool inexact = false;
530 /* Determine if we need to add or subtract. */
532 subtract_p = (sign ^ b->sign) ^ subtract_p;
534 switch (CLASS2 (a->class, b->class))
536 case CLASS2 (rvc_zero, rvc_zero):
537 /* -0 + -0 = -0, -0 - +0 = -0; all other cases yield +0. */
538 get_zero (r, sign & !subtract_p);
541 case CLASS2 (rvc_zero, rvc_normal):
542 case CLASS2 (rvc_zero, rvc_inf):
543 case CLASS2 (rvc_zero, rvc_nan):
545 case CLASS2 (rvc_normal, rvc_nan):
546 case CLASS2 (rvc_inf, rvc_nan):
547 case CLASS2 (rvc_nan, rvc_nan):
548 /* ANY + NaN = NaN. */
549 case CLASS2 (rvc_normal, rvc_inf):
552 r->sign = sign ^ subtract_p;
555 case CLASS2 (rvc_normal, rvc_zero):
556 case CLASS2 (rvc_inf, rvc_zero):
557 case CLASS2 (rvc_nan, rvc_zero):
559 case CLASS2 (rvc_nan, rvc_normal):
560 case CLASS2 (rvc_nan, rvc_inf):
561 /* NaN + ANY = NaN. */
562 case CLASS2 (rvc_inf, rvc_normal):
567 case CLASS2 (rvc_inf, rvc_inf):
569 /* Inf - Inf = NaN. */
570 get_canonical_qnan (r, 0);
572 /* Inf + Inf = Inf. */
576 case CLASS2 (rvc_normal, rvc_normal):
583 /* Swap the arguments such that A has the larger exponent. */
584 dexp = a->exp - b->exp;
587 const REAL_VALUE_TYPE *t;
594 /* If the exponents are not identical, we need to shift the
595 significand of B down. */
598 /* If the exponents are too far apart, the significands
599 do not overlap, which makes the subtraction a noop. */
600 if (dexp >= SIGNIFICAND_BITS)
607 inexact |= sticky_rshift_significand (&t, b, dexp);
613 if (sub_significands (r, a, b, inexact))
615 /* We got a borrow out of the subtraction. That means that
616 A and B had the same exponent, and B had the larger
617 significand. We need to swap the sign and negate the
620 neg_significand (r, r);
625 if (add_significands (r, a, b))
627 /* We got carry out of the addition. This means we need to
628 shift the significand back down one bit and increase the
630 inexact |= sticky_rshift_significand (r, r, 1);
631 r->sig[SIGSZ-1] |= SIG_MSB;
640 r->class = rvc_normal;
643 /* Zero out the remaining fields. */
647 /* Re-normalize the result. */
650 /* Special case: if the subtraction results in zero, the result
652 if (r->class == rvc_zero)
655 r->sig[0] |= inexact;
660 /* Calculate R = A * B. Return true if the result may be inexact. */
663 do_multiply (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *a,
664 const REAL_VALUE_TYPE *b)
666 REAL_VALUE_TYPE u, t, *rr;
667 unsigned int i, j, k;
668 int sign = a->sign ^ b->sign;
669 bool inexact = false;
671 switch (CLASS2 (a->class, b->class))
673 case CLASS2 (rvc_zero, rvc_zero):
674 case CLASS2 (rvc_zero, rvc_normal):
675 case CLASS2 (rvc_normal, rvc_zero):
676 /* +-0 * ANY = 0 with appropriate sign. */
680 case CLASS2 (rvc_zero, rvc_nan):
681 case CLASS2 (rvc_normal, rvc_nan):
682 case CLASS2 (rvc_inf, rvc_nan):
683 case CLASS2 (rvc_nan, rvc_nan):
684 /* ANY * NaN = NaN. */
689 case CLASS2 (rvc_nan, rvc_zero):
690 case CLASS2 (rvc_nan, rvc_normal):
691 case CLASS2 (rvc_nan, rvc_inf):
692 /* NaN * ANY = NaN. */
697 case CLASS2 (rvc_zero, rvc_inf):
698 case CLASS2 (rvc_inf, rvc_zero):
700 get_canonical_qnan (r, sign);
703 case CLASS2 (rvc_inf, rvc_inf):
704 case CLASS2 (rvc_normal, rvc_inf):
705 case CLASS2 (rvc_inf, rvc_normal):
706 /* Inf * Inf = Inf, R * Inf = Inf */
710 case CLASS2 (rvc_normal, rvc_normal):
717 if (r == a || r == b)
723 /* Collect all the partial products. Since we don't have sure access
724 to a widening multiply, we split each long into two half-words.
726 Consider the long-hand form of a four half-word multiplication:
736 We construct partial products of the widened half-word products
737 that are known to not overlap, e.g. DF+DH. Each such partial
738 product is given its proper exponent, which allows us to sum them
739 and obtain the finished product. */
741 for (i = 0; i < SIGSZ * 2; ++i)
743 unsigned long ai = a->sig[i / 2];
745 ai >>= HOST_BITS_PER_LONG / 2;
747 ai &= ((unsigned long)1 << (HOST_BITS_PER_LONG / 2)) - 1;
752 for (j = 0; j < 2; ++j)
754 int exp = (a->exp - (2*SIGSZ-1-i)*(HOST_BITS_PER_LONG/2)
755 + (b->exp - (1-j)*(HOST_BITS_PER_LONG/2)));
764 /* Would underflow to zero, which we shouldn't bother adding. */
769 memset (&u, 0, sizeof (u));
770 u.class = rvc_normal;
773 for (k = j; k < SIGSZ * 2; k += 2)
775 unsigned long bi = b->sig[k / 2];
777 bi >>= HOST_BITS_PER_LONG / 2;
779 bi &= ((unsigned long)1 << (HOST_BITS_PER_LONG / 2)) - 1;
781 u.sig[k / 2] = ai * bi;
785 inexact |= do_add (rr, rr, &u, 0);
796 /* Calculate R = A / B. Return true if the result may be inexact. */
799 do_divide (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *a,
800 const REAL_VALUE_TYPE *b)
802 int exp, sign = a->sign ^ b->sign;
803 REAL_VALUE_TYPE t, *rr;
806 switch (CLASS2 (a->class, b->class))
808 case CLASS2 (rvc_zero, rvc_zero):
810 case CLASS2 (rvc_inf, rvc_inf):
811 /* Inf / Inf = NaN. */
812 get_canonical_qnan (r, sign);
815 case CLASS2 (rvc_zero, rvc_normal):
816 case CLASS2 (rvc_zero, rvc_inf):
818 case CLASS2 (rvc_normal, rvc_inf):
823 case CLASS2 (rvc_normal, rvc_zero):
825 case CLASS2 (rvc_inf, rvc_zero):
830 case CLASS2 (rvc_zero, rvc_nan):
831 case CLASS2 (rvc_normal, rvc_nan):
832 case CLASS2 (rvc_inf, rvc_nan):
833 case CLASS2 (rvc_nan, rvc_nan):
834 /* ANY / NaN = NaN. */
839 case CLASS2 (rvc_nan, rvc_zero):
840 case CLASS2 (rvc_nan, rvc_normal):
841 case CLASS2 (rvc_nan, rvc_inf):
842 /* NaN / ANY = NaN. */
847 case CLASS2 (rvc_inf, rvc_normal):
852 case CLASS2 (rvc_normal, rvc_normal):
859 if (r == a || r == b)
864 /* Make sure all fields in the result are initialized. */
866 rr->class = rvc_normal;
869 exp = a->exp - b->exp + 1;
882 inexact = div_significands (rr, a, b);
884 /* Re-normalize the result. */
886 rr->sig[0] |= inexact;
894 /* Return a tri-state comparison of A vs B. Return NAN_RESULT if
895 one of the two operands is a NaN. */
898 do_compare (const REAL_VALUE_TYPE *a, const REAL_VALUE_TYPE *b,
903 switch (CLASS2 (a->class, b->class))
905 case CLASS2 (rvc_zero, rvc_zero):
906 /* Sign of zero doesn't matter for compares. */
909 case CLASS2 (rvc_inf, rvc_zero):
910 case CLASS2 (rvc_inf, rvc_normal):
911 case CLASS2 (rvc_normal, rvc_zero):
912 return (a->sign ? -1 : 1);
914 case CLASS2 (rvc_inf, rvc_inf):
915 return -a->sign - -b->sign;
917 case CLASS2 (rvc_zero, rvc_normal):
918 case CLASS2 (rvc_zero, rvc_inf):
919 case CLASS2 (rvc_normal, rvc_inf):
920 return (b->sign ? 1 : -1);
922 case CLASS2 (rvc_zero, rvc_nan):
923 case CLASS2 (rvc_normal, rvc_nan):
924 case CLASS2 (rvc_inf, rvc_nan):
925 case CLASS2 (rvc_nan, rvc_nan):
926 case CLASS2 (rvc_nan, rvc_zero):
927 case CLASS2 (rvc_nan, rvc_normal):
928 case CLASS2 (rvc_nan, rvc_inf):
931 case CLASS2 (rvc_normal, rvc_normal):
938 if (a->sign != b->sign)
939 return -a->sign - -b->sign;
943 else if (a->exp < b->exp)
946 ret = cmp_significands (a, b);
948 return (a->sign ? -ret : ret);
951 /* Return A truncated to an integral value toward zero. */
954 do_fix_trunc (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *a)
967 get_zero (r, r->sign);
968 else if (r->exp < SIGNIFICAND_BITS)
969 clear_significand_below (r, SIGNIFICAND_BITS - r->exp);
977 /* Perform the binary or unary operation described by CODE.
978 For a unary operation, leave OP1 NULL. */
981 real_arithmetic (REAL_VALUE_TYPE *r, int icode, const REAL_VALUE_TYPE *op0,
982 const REAL_VALUE_TYPE *op1)
984 enum tree_code code = icode;
989 do_add (r, op0, op1, 0);
993 do_add (r, op0, op1, 1);
997 do_multiply (r, op0, op1);
1001 do_divide (r, op0, op1);
1005 if (op1->class == rvc_nan)
1007 else if (do_compare (op0, op1, -1) < 0)
1014 if (op1->class == rvc_nan)
1016 else if (do_compare (op0, op1, 1) < 0)
1032 case FIX_TRUNC_EXPR:
1033 do_fix_trunc (r, op0);
1041 /* Legacy. Similar, but return the result directly. */
1044 real_arithmetic2 (int icode, const REAL_VALUE_TYPE *op0,
1045 const REAL_VALUE_TYPE *op1)
1048 real_arithmetic (&r, icode, op0, op1);
1053 real_compare (int icode, const REAL_VALUE_TYPE *op0,
1054 const REAL_VALUE_TYPE *op1)
1056 enum tree_code code = icode;
1061 return do_compare (op0, op1, 1) < 0;
1063 return do_compare (op0, op1, 1) <= 0;
1065 return do_compare (op0, op1, -1) > 0;
1067 return do_compare (op0, op1, -1) >= 0;
1069 return do_compare (op0, op1, -1) == 0;
1071 return do_compare (op0, op1, -1) != 0;
1072 case UNORDERED_EXPR:
1073 return op0->class == rvc_nan || op1->class == rvc_nan;
1075 return op0->class != rvc_nan && op1->class != rvc_nan;
1077 return do_compare (op0, op1, -1) < 0;
1079 return do_compare (op0, op1, -1) <= 0;
1081 return do_compare (op0, op1, 1) > 0;
1083 return do_compare (op0, op1, 1) >= 0;
1085 return do_compare (op0, op1, 0) == 0;
1092 /* Return floor log2(R). */
1095 real_exponent (const REAL_VALUE_TYPE *r)
1103 return (unsigned int)-1 >> 1;
1111 /* R = OP0 * 2**EXP. */
1114 real_ldexp (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *op0, int exp)
1127 get_inf (r, r->sign);
1128 else if (exp < -MAX_EXP)
1129 get_zero (r, r->sign);
1139 /* Determine whether a floating-point value X is infinite. */
1142 real_isinf (const REAL_VALUE_TYPE *r)
1144 return (r->class == rvc_inf);
1147 /* Determine whether a floating-point value X is a NaN. */
1150 real_isnan (const REAL_VALUE_TYPE *r)
1152 return (r->class == rvc_nan);
1155 /* Determine whether a floating-point value X is negative. */
1158 real_isneg (const REAL_VALUE_TYPE *r)
1163 /* Determine whether a floating-point value X is minus zero. */
1166 real_isnegzero (const REAL_VALUE_TYPE *r)
1168 return r->sign && r->class == rvc_zero;
1171 /* Compare two floating-point objects for bitwise identity. */
1174 real_identical (const REAL_VALUE_TYPE *a, const REAL_VALUE_TYPE *b)
1178 if (a->class != b->class)
1180 if (a->sign != b->sign)
1190 if (a->exp != b->exp)
1195 if (a->signalling != b->signalling)
1197 /* The significand is ignored for canonical NaNs. */
1198 if (a->canonical || b->canonical)
1199 return a->canonical == b->canonical;
1206 for (i = 0; i < SIGSZ; ++i)
1207 if (a->sig[i] != b->sig[i])
1213 /* Try to change R into its exact multiplicative inverse in machine
1214 mode MODE. Return true if successful. */
1217 exact_real_inverse (enum machine_mode mode, REAL_VALUE_TYPE *r)
1219 const REAL_VALUE_TYPE *one = real_digit (1);
1223 if (r->class != rvc_normal)
1226 /* Check for a power of two: all significand bits zero except the MSB. */
1227 for (i = 0; i < SIGSZ-1; ++i)
1230 if (r->sig[SIGSZ-1] != SIG_MSB)
1233 /* Find the inverse and truncate to the required mode. */
1234 do_divide (&u, one, r);
1235 real_convert (&u, mode, &u);
1237 /* The rounding may have overflowed. */
1238 if (u.class != rvc_normal)
1240 for (i = 0; i < SIGSZ-1; ++i)
1243 if (u.sig[SIGSZ-1] != SIG_MSB)
1250 /* Render R as an integer. */
1253 real_to_integer (const REAL_VALUE_TYPE *r)
1255 unsigned HOST_WIDE_INT i;
1266 i = (unsigned HOST_WIDE_INT) 1 << (HOST_BITS_PER_WIDE_INT - 1);
1274 /* Only force overflow for unsigned overflow. Signed overflow is
1275 undefined, so it doesn't matter what we return, and some callers
1276 expect to be able to use this routine for both signed and
1277 unsigned conversions. */
1278 if (r->exp > HOST_BITS_PER_WIDE_INT)
1281 if (HOST_BITS_PER_WIDE_INT == HOST_BITS_PER_LONG)
1282 i = r->sig[SIGSZ-1];
1283 else if (HOST_BITS_PER_WIDE_INT == 2*HOST_BITS_PER_LONG)
1285 i = r->sig[SIGSZ-1];
1286 i = i << (HOST_BITS_PER_LONG - 1) << 1;
1287 i |= r->sig[SIGSZ-2];
1292 i >>= HOST_BITS_PER_WIDE_INT - r->exp;
1303 /* Likewise, but to an integer pair, HI+LOW. */
1306 real_to_integer2 (HOST_WIDE_INT *plow, HOST_WIDE_INT *phigh,
1307 const REAL_VALUE_TYPE *r)
1310 HOST_WIDE_INT low, high;
1323 high = (unsigned HOST_WIDE_INT) 1 << (HOST_BITS_PER_WIDE_INT - 1);
1337 /* Only force overflow for unsigned overflow. Signed overflow is
1338 undefined, so it doesn't matter what we return, and some callers
1339 expect to be able to use this routine for both signed and
1340 unsigned conversions. */
1341 if (exp > 2*HOST_BITS_PER_WIDE_INT)
1344 rshift_significand (&t, r, 2*HOST_BITS_PER_WIDE_INT - exp);
1345 if (HOST_BITS_PER_WIDE_INT == HOST_BITS_PER_LONG)
1347 high = t.sig[SIGSZ-1];
1348 low = t.sig[SIGSZ-2];
1350 else if (HOST_BITS_PER_WIDE_INT == 2*HOST_BITS_PER_LONG)
1352 high = t.sig[SIGSZ-1];
1353 high = high << (HOST_BITS_PER_LONG - 1) << 1;
1354 high |= t.sig[SIGSZ-2];
1356 low = t.sig[SIGSZ-3];
1357 low = low << (HOST_BITS_PER_LONG - 1) << 1;
1358 low |= t.sig[SIGSZ-4];
1368 low = -low, high = ~high;
1380 /* A subroutine of real_to_decimal. Compute the quotient and remainder
1381 of NUM / DEN. Return the quotient and place the remainder in NUM.
1382 It is expected that NUM / DEN are close enough that the quotient is
1385 static unsigned long
1386 rtd_divmod (REAL_VALUE_TYPE *num, REAL_VALUE_TYPE *den)
1388 unsigned long q, msb;
1389 int expn = num->exp, expd = den->exp;
1398 msb = num->sig[SIGSZ-1] & SIG_MSB;
1400 lshift_significand_1 (num, num);
1402 if (msb || cmp_significands (num, den) >= 0)
1404 sub_significands (num, num, den, 0);
1408 while (--expn >= expd);
1416 /* Render R as a decimal floating point constant. Emit DIGITS significant
1417 digits in the result, bounded by BUF_SIZE. If DIGITS is 0, choose the
1418 maximum for the representation. If CROP_TRAILING_ZEROS, strip trailing
1421 #define M_LOG10_2 0.30102999566398119521
1424 real_to_decimal (char *str, const REAL_VALUE_TYPE *r_orig, size_t buf_size,
1425 size_t digits, int crop_trailing_zeros)
1427 const REAL_VALUE_TYPE *one, *ten;
1428 REAL_VALUE_TYPE r, pten, u, v;
1429 int dec_exp, cmp_one, digit;
1431 char *p, *first, *last;
1438 strcpy (str, (r.sign ? "-0.0" : "0.0"));
1443 strcpy (str, (r.sign ? "-Inf" : "+Inf"));
1446 /* ??? Print the significand as well, if not canonical? */
1447 strcpy (str, (r.sign ? "-NaN" : "+NaN"));
1453 /* Bound the number of digits printed by the size of the representation. */
1454 max_digits = SIGNIFICAND_BITS * M_LOG10_2;
1455 if (digits == 0 || digits > max_digits)
1456 digits = max_digits;
1458 /* Estimate the decimal exponent, and compute the length of the string it
1459 will print as. Be conservative and add one to account for possible
1460 overflow or rounding error. */
1461 dec_exp = r.exp * M_LOG10_2;
1462 for (max_digits = 1; dec_exp ; max_digits++)
1465 /* Bound the number of digits printed by the size of the output buffer. */
1466 max_digits = buf_size - 1 - 1 - 2 - max_digits - 1;
1467 if (max_digits > buf_size)
1469 if (digits > max_digits)
1470 digits = max_digits;
1472 one = real_digit (1);
1473 ten = ten_to_ptwo (0);
1481 cmp_one = do_compare (&r, one, 0);
1486 /* Number is greater than one. Convert significand to an integer
1487 and strip trailing decimal zeros. */
1490 u.exp = SIGNIFICAND_BITS - 1;
1492 /* Largest M, such that 10**2**M fits within SIGNIFICAND_BITS. */
1493 m = floor_log2 (max_digits);
1495 /* Iterate over the bits of the possible powers of 10 that might
1496 be present in U and eliminate them. That is, if we find that
1497 10**2**M divides U evenly, keep the division and increase
1503 do_divide (&t, &u, ten_to_ptwo (m));
1504 do_fix_trunc (&v, &t);
1505 if (cmp_significands (&v, &t) == 0)
1513 /* Revert the scaling to integer that we performed earlier. */
1514 u.exp += r.exp - (SIGNIFICAND_BITS - 1);
1517 /* Find power of 10. Do this by dividing out 10**2**M when
1518 this is larger than the current remainder. Fill PTEN with
1519 the power of 10 that we compute. */
1522 m = floor_log2 ((int)(r.exp * M_LOG10_2)) + 1;
1525 const REAL_VALUE_TYPE *ptentwo = ten_to_ptwo (m);
1526 if (do_compare (&u, ptentwo, 0) >= 0)
1528 do_divide (&u, &u, ptentwo);
1529 do_multiply (&pten, &pten, ptentwo);
1536 /* We managed to divide off enough tens in the above reduction
1537 loop that we've now got a negative exponent. Fall into the
1538 less-than-one code to compute the proper value for PTEN. */
1545 /* Number is less than one. Pad significand with leading
1551 /* Stop if we'd shift bits off the bottom. */
1555 do_multiply (&u, &v, ten);
1557 /* Stop if we're now >= 1. */
1566 /* Find power of 10. Do this by multiplying in P=10**2**M when
1567 the current remainder is smaller than 1/P. Fill PTEN with the
1568 power of 10 that we compute. */
1569 m = floor_log2 ((int)(-r.exp * M_LOG10_2)) + 1;
1572 const REAL_VALUE_TYPE *ptentwo = ten_to_ptwo (m);
1573 const REAL_VALUE_TYPE *ptenmtwo = ten_to_mptwo (m);
1575 if (do_compare (&v, ptenmtwo, 0) <= 0)
1577 do_multiply (&v, &v, ptentwo);
1578 do_multiply (&pten, &pten, ptentwo);
1584 /* Invert the positive power of 10 that we've collected so far. */
1585 do_divide (&pten, one, &pten);
1593 /* At this point, PTEN should contain the nearest power of 10 smaller
1594 than R, such that this division produces the first digit.
1596 Using a divide-step primitive that returns the complete integral
1597 remainder avoids the rounding error that would be produced if
1598 we were to use do_divide here and then simply multiply by 10 for
1599 each subsequent digit. */
1601 digit = rtd_divmod (&r, &pten);
1603 /* Be prepared for error in that division via underflow ... */
1604 if (digit == 0 && cmp_significand_0 (&r))
1606 /* Multiply by 10 and try again. */
1607 do_multiply (&r, &r, ten);
1608 digit = rtd_divmod (&r, &pten);
1614 /* ... or overflow. */
1622 else if (digit > 10)
1627 /* Generate subsequent digits. */
1628 while (--digits > 0)
1630 do_multiply (&r, &r, ten);
1631 digit = rtd_divmod (&r, &pten);
1636 /* Generate one more digit with which to do rounding. */
1637 do_multiply (&r, &r, ten);
1638 digit = rtd_divmod (&r, &pten);
1640 /* Round the result. */
1643 /* Round to nearest. If R is nonzero there are additional
1644 nonzero digits to be extracted. */
1645 if (cmp_significand_0 (&r))
1647 /* Round to even. */
1648 else if ((p[-1] - '0') & 1)
1665 /* Carry out of the first digit. This means we had all 9's and
1666 now have all 0's. "Prepend" a 1 by overwriting the first 0. */
1674 /* Insert the decimal point. */
1675 first[0] = first[1];
1678 /* If requested, drop trailing zeros. Never crop past "1.0". */
1679 if (crop_trailing_zeros)
1680 while (last > first + 3 && last[-1] == '0')
1683 /* Append the exponent. */
1684 sprintf (last, "e%+d", dec_exp);
1687 /* Render R as a hexadecimal floating point constant. Emit DIGITS
1688 significant digits in the result, bounded by BUF_SIZE. If DIGITS is 0,
1689 choose the maximum for the representation. If CROP_TRAILING_ZEROS,
1690 strip trailing zeros. */
1693 real_to_hexadecimal (char *str, const REAL_VALUE_TYPE *r, size_t buf_size,
1694 size_t digits, int crop_trailing_zeros)
1696 int i, j, exp = r->exp;
1709 strcpy (str, (r->sign ? "-Inf" : "+Inf"));
1712 /* ??? Print the significand as well, if not canonical? */
1713 strcpy (str, (r->sign ? "-NaN" : "+NaN"));
1720 digits = SIGNIFICAND_BITS / 4;
1722 /* Bound the number of digits printed by the size of the output buffer. */
1724 sprintf (exp_buf, "p%+d", exp);
1725 max_digits = buf_size - strlen (exp_buf) - r->sign - 4 - 1;
1726 if (max_digits > buf_size)
1728 if (digits > max_digits)
1729 digits = max_digits;
1740 for (i = SIGSZ - 1; i >= 0; --i)
1741 for (j = HOST_BITS_PER_LONG - 4; j >= 0; j -= 4)
1743 *p++ = "0123456789abcdef"[(r->sig[i] >> j) & 15];
1749 if (crop_trailing_zeros)
1750 while (p > first + 1 && p[-1] == '0')
1753 sprintf (p, "p%+d", exp);
1756 /* Initialize R from a decimal or hexadecimal string. The string is
1757 assumed to have been syntax checked already. */
1760 real_from_string (REAL_VALUE_TYPE *r, const char *str)
1772 else if (*str == '+')
1775 if (str[0] == '0' && (str[1] == 'x' || str[1] == 'X'))
1777 /* Hexadecimal floating point. */
1778 int pos = SIGNIFICAND_BITS - 4, d;
1786 d = hex_value (*str);
1791 r->sig[pos / HOST_BITS_PER_LONG]
1792 |= (unsigned long) d << (pos % HOST_BITS_PER_LONG);
1801 if (pos == SIGNIFICAND_BITS - 4)
1808 d = hex_value (*str);
1813 r->sig[pos / HOST_BITS_PER_LONG]
1814 |= (unsigned long) d << (pos % HOST_BITS_PER_LONG);
1820 if (*str == 'p' || *str == 'P')
1822 bool exp_neg = false;
1830 else if (*str == '+')
1834 while (ISDIGIT (*str))
1840 /* Overflowed the exponent. */
1854 r->class = rvc_normal;
1861 /* Decimal floating point. */
1862 const REAL_VALUE_TYPE *ten = ten_to_ptwo (0);
1867 while (ISDIGIT (*str))
1870 do_multiply (r, r, ten);
1872 do_add (r, r, real_digit (d), 0);
1877 if (r->class == rvc_zero)
1882 while (ISDIGIT (*str))
1885 do_multiply (r, r, ten);
1887 do_add (r, r, real_digit (d), 0);
1892 if (*str == 'e' || *str == 'E')
1894 bool exp_neg = false;
1902 else if (*str == '+')
1906 while (ISDIGIT (*str))
1912 /* Overflowed the exponent. */
1926 times_pten (r, exp);
1941 /* Legacy. Similar, but return the result directly. */
1944 real_from_string2 (const char *s, enum machine_mode mode)
1948 real_from_string (&r, s);
1949 if (mode != VOIDmode)
1950 real_convert (&r, mode, &r);
1955 /* Initialize R from the integer pair HIGH+LOW. */
1958 real_from_integer (REAL_VALUE_TYPE *r, enum machine_mode mode,
1959 unsigned HOST_WIDE_INT low, HOST_WIDE_INT high,
1962 if (low == 0 && high == 0)
1966 memset (r, 0, sizeof (*r));
1967 r->class = rvc_normal;
1968 r->sign = high < 0 && !unsigned_p;
1969 r->exp = 2 * HOST_BITS_PER_WIDE_INT;
1980 if (HOST_BITS_PER_LONG == HOST_BITS_PER_WIDE_INT)
1982 r->sig[SIGSZ-1] = high;
1983 r->sig[SIGSZ-2] = low;
1985 else if (HOST_BITS_PER_LONG*2 == HOST_BITS_PER_WIDE_INT)
1987 r->sig[SIGSZ-1] = high >> (HOST_BITS_PER_LONG - 1) >> 1;
1988 r->sig[SIGSZ-2] = high;
1989 r->sig[SIGSZ-3] = low >> (HOST_BITS_PER_LONG - 1) >> 1;
1990 r->sig[SIGSZ-4] = low;
1998 if (mode != VOIDmode)
1999 real_convert (r, mode, r);
2002 /* Returns 10**2**N. */
2004 static const REAL_VALUE_TYPE *
2007 static REAL_VALUE_TYPE tens[EXP_BITS];
2009 if (n < 0 || n >= EXP_BITS)
2012 if (tens[n].class == rvc_zero)
2014 if (n < (HOST_BITS_PER_WIDE_INT == 64 ? 5 : 4))
2016 HOST_WIDE_INT t = 10;
2019 for (i = 0; i < n; ++i)
2022 real_from_integer (&tens[n], VOIDmode, t, 0, 1);
2026 const REAL_VALUE_TYPE *t = ten_to_ptwo (n - 1);
2027 do_multiply (&tens[n], t, t);
2034 /* Returns 10**(-2**N). */
2036 static const REAL_VALUE_TYPE *
2037 ten_to_mptwo (int n)
2039 static REAL_VALUE_TYPE tens[EXP_BITS];
2041 if (n < 0 || n >= EXP_BITS)
2044 if (tens[n].class == rvc_zero)
2045 do_divide (&tens[n], real_digit (1), ten_to_ptwo (n));
2052 static const REAL_VALUE_TYPE *
2055 static REAL_VALUE_TYPE num[10];
2060 if (n > 0 && num[n].class == rvc_zero)
2061 real_from_integer (&num[n], VOIDmode, n, 0, 1);
2066 /* Multiply R by 10**EXP. */
2069 times_pten (REAL_VALUE_TYPE *r, int exp)
2071 REAL_VALUE_TYPE pten, *rr;
2072 bool negative = (exp < 0);
2078 pten = *real_digit (1);
2084 for (i = 0; exp > 0; ++i, exp >>= 1)
2086 do_multiply (rr, rr, ten_to_ptwo (i));
2089 do_divide (r, r, &pten);
2092 /* Fills R with +Inf. */
2095 real_inf (REAL_VALUE_TYPE *r)
2100 /* Fills R with a NaN whose significand is described by STR. If QUIET,
2101 we force a QNaN, else we force an SNaN. The string, if not empty,
2102 is parsed as a number and placed in the significand. Return true
2103 if the string was successfully parsed. */
2106 real_nan (REAL_VALUE_TYPE *r, const char *str, int quiet,
2107 enum machine_mode mode)
2109 const struct real_format *fmt;
2111 fmt = REAL_MODE_FORMAT (mode);
2118 get_canonical_qnan (r, 0);
2120 get_canonical_snan (r, 0);
2127 memset (r, 0, sizeof (*r));
2130 /* Parse akin to strtol into the significand of R. */
2132 while (ISSPACE (*str))
2136 else if (*str == '+')
2146 while ((d = hex_value (*str)) < base)
2153 lshift_significand (r, r, 3);
2156 lshift_significand (r, r, 4);
2159 lshift_significand_1 (&u, r);
2160 lshift_significand (r, r, 3);
2161 add_significands (r, r, &u);
2169 add_significands (r, r, &u);
2174 /* Must have consumed the entire string for success. */
2178 /* Shift the significand into place such that the bits
2179 are in the most significant bits for the format. */
2180 lshift_significand (r, r, SIGNIFICAND_BITS - fmt->pnan);
2182 /* Our MSB is always unset for NaNs. */
2183 r->sig[SIGSZ-1] &= ~SIG_MSB;
2185 /* Force quiet or signalling NaN. */
2186 r->signalling = !quiet;
2192 /* Fills R with the largest finite value representable in mode MODE.
2193 If SIGN is nonzero, R is set to the most negative finite value. */
2196 real_maxval (REAL_VALUE_TYPE *r, int sign, enum machine_mode mode)
2198 const struct real_format *fmt;
2201 fmt = REAL_MODE_FORMAT (mode);
2205 r->class = rvc_normal;
2209 r->exp = fmt->emax * fmt->log2_b;
2211 np2 = SIGNIFICAND_BITS - fmt->p * fmt->log2_b;
2212 memset (r->sig, -1, SIGSZ * sizeof (unsigned long));
2213 clear_significand_below (r, np2);
2216 /* Fills R with 2**N. */
2219 real_2expN (REAL_VALUE_TYPE *r, int n)
2221 memset (r, 0, sizeof (*r));
2226 else if (n < -MAX_EXP)
2230 r->class = rvc_normal;
2232 r->sig[SIGSZ-1] = SIG_MSB;
2238 round_for_format (const struct real_format *fmt, REAL_VALUE_TYPE *r)
2241 unsigned long sticky;
2245 p2 = fmt->p * fmt->log2_b;
2246 emin2m1 = (fmt->emin - 1) * fmt->log2_b;
2247 emax2 = fmt->emax * fmt->log2_b;
2249 np2 = SIGNIFICAND_BITS - p2;
2253 get_zero (r, r->sign);
2255 if (!fmt->has_signed_zero)
2260 get_inf (r, r->sign);
2265 clear_significand_below (r, np2);
2275 /* If we're not base2, normalize the exponent to a multiple of
2277 if (fmt->log2_b != 1)
2279 int shift = r->exp & (fmt->log2_b - 1);
2282 shift = fmt->log2_b - shift;
2283 r->sig[0] |= sticky_rshift_significand (r, r, shift);
2288 /* Check the range of the exponent. If we're out of range,
2289 either underflow or overflow. */
2292 else if (r->exp <= emin2m1)
2296 if (!fmt->has_denorm)
2298 /* Don't underflow completely until we've had a chance to round. */
2299 if (r->exp < emin2m1)
2304 diff = emin2m1 - r->exp + 1;
2308 /* De-normalize the significand. */
2309 r->sig[0] |= sticky_rshift_significand (r, r, diff);
2314 /* There are P2 true significand bits, followed by one guard bit,
2315 followed by one sticky bit, followed by stuff. Fold nonzero
2316 stuff into the sticky bit. */
2319 for (i = 0, w = (np2 - 1) / HOST_BITS_PER_LONG; i < w; ++i)
2320 sticky |= r->sig[i];
2322 r->sig[w] & (((unsigned long)1 << ((np2 - 1) % HOST_BITS_PER_LONG)) - 1);
2324 guard = test_significand_bit (r, np2 - 1);
2325 lsb = test_significand_bit (r, np2);
2327 /* Round to even. */
2328 if (guard && (sticky || lsb))
2332 set_significand_bit (&u, np2);
2334 if (add_significands (r, r, &u))
2336 /* Overflow. Means the significand had been all ones, and
2337 is now all zeros. Need to increase the exponent, and
2338 possibly re-normalize it. */
2339 if (++r->exp > emax2)
2341 r->sig[SIGSZ-1] = SIG_MSB;
2343 if (fmt->log2_b != 1)
2345 int shift = r->exp & (fmt->log2_b - 1);
2348 shift = fmt->log2_b - shift;
2349 rshift_significand (r, r, shift);
2358 /* Catch underflow that we deferred until after rounding. */
2359 if (r->exp <= emin2m1)
2362 /* Clear out trailing garbage. */
2363 clear_significand_below (r, np2);
2366 /* Extend or truncate to a new mode. */
2369 real_convert (REAL_VALUE_TYPE *r, enum machine_mode mode,
2370 const REAL_VALUE_TYPE *a)
2372 const struct real_format *fmt;
2374 fmt = REAL_MODE_FORMAT (mode);
2379 round_for_format (fmt, r);
2381 /* round_for_format de-normalizes denormals. Undo just that part. */
2382 if (r->class == rvc_normal)
2386 /* Legacy. Likewise, except return the struct directly. */
2389 real_value_truncate (enum machine_mode mode, REAL_VALUE_TYPE a)
2392 real_convert (&r, mode, &a);
2396 /* Return true if truncating to MODE is exact. */
2399 exact_real_truncate (enum machine_mode mode, const REAL_VALUE_TYPE *a)
2402 real_convert (&t, mode, a);
2403 return real_identical (&t, a);
2406 /* Write R to the given target format. Place the words of the result
2407 in target word order in BUF. There are always 32 bits in each
2408 long, no matter the size of the host long.
2410 Legacy: return word 0 for implementing REAL_VALUE_TO_TARGET_SINGLE. */
2413 real_to_target_fmt (long *buf, const REAL_VALUE_TYPE *r_orig,
2414 const struct real_format *fmt)
2420 round_for_format (fmt, &r);
2424 (*fmt->encode) (fmt, buf, &r);
2429 /* Similar, but look up the format from MODE. */
2432 real_to_target (long *buf, const REAL_VALUE_TYPE *r, enum machine_mode mode)
2434 const struct real_format *fmt;
2436 fmt = REAL_MODE_FORMAT (mode);
2440 return real_to_target_fmt (buf, r, fmt);
2443 /* Read R from the given target format. Read the words of the result
2444 in target word order in BUF. There are always 32 bits in each
2445 long, no matter the size of the host long. */
2448 real_from_target_fmt (REAL_VALUE_TYPE *r, const long *buf,
2449 const struct real_format *fmt)
2451 (*fmt->decode) (fmt, r, buf);
2454 /* Similar, but look up the format from MODE. */
2457 real_from_target (REAL_VALUE_TYPE *r, const long *buf, enum machine_mode mode)
2459 const struct real_format *fmt;
2461 fmt = REAL_MODE_FORMAT (mode);
2465 (*fmt->decode) (fmt, r, buf);
2468 /* Return the number of bits in the significand for MODE. */
2469 /* ??? Legacy. Should get access to real_format directly. */
2472 significand_size (enum machine_mode mode)
2474 const struct real_format *fmt;
2476 fmt = REAL_MODE_FORMAT (mode);
2480 return fmt->p * fmt->log2_b;
2483 /* Return a hash value for the given real value. */
2484 /* ??? The "unsigned int" return value is intended to be hashval_t,
2485 but I didn't want to pull hashtab.h into real.h. */
2488 real_hash (const REAL_VALUE_TYPE *r)
2493 h = r->class | (r->sign << 2);
2506 h ^= (unsigned int)-1;
2515 if (sizeof(unsigned long) > sizeof(unsigned int))
2516 for (i = 0; i < SIGSZ; ++i)
2518 unsigned long s = r->sig[i];
2519 h ^= s ^ (s >> (HOST_BITS_PER_LONG / 2));
2522 for (i = 0; i < SIGSZ; ++i)
2528 /* IEEE single-precision format. */
2530 static void encode_ieee_single (const struct real_format *fmt,
2531 long *, const REAL_VALUE_TYPE *);
2532 static void decode_ieee_single (const struct real_format *,
2533 REAL_VALUE_TYPE *, const long *);
2536 encode_ieee_single (const struct real_format *fmt, long *buf,
2537 const REAL_VALUE_TYPE *r)
2539 unsigned long image, sig, exp;
2540 unsigned long sign = r->sign;
2541 bool denormal = (r->sig[SIGSZ-1] & SIG_MSB) == 0;
2544 sig = (r->sig[SIGSZ-1] >> (HOST_BITS_PER_LONG - 24)) & 0x7fffff;
2555 image |= 0x7fffffff;
2563 if (r->signalling == fmt->qnan_msb_set)
2567 /* We overload qnan_msb_set here: it's only clear for
2568 mips_ieee_single, which wants all mantissa bits but the
2569 quiet/signalling one set in canonical NaNs (at least
2571 if (r->canonical && !fmt->qnan_msb_set)
2572 sig |= (1 << 22) - 1;
2580 image |= 0x7fffffff;
2584 /* Recall that IEEE numbers are interpreted as 1.F x 2**exp,
2585 whereas the intermediate representation is 0.F x 2**exp.
2586 Which means we're off by one. */
2590 exp = r->exp + 127 - 1;
2603 decode_ieee_single (const struct real_format *fmt, REAL_VALUE_TYPE *r,
2606 unsigned long image = buf[0] & 0xffffffff;
2607 bool sign = (image >> 31) & 1;
2608 int exp = (image >> 23) & 0xff;
2610 memset (r, 0, sizeof (*r));
2611 image <<= HOST_BITS_PER_LONG - 24;
2616 if (image && fmt->has_denorm)
2618 r->class = rvc_normal;
2621 r->sig[SIGSZ-1] = image << 1;
2624 else if (fmt->has_signed_zero)
2627 else if (exp == 255 && (fmt->has_nans || fmt->has_inf))
2633 r->signalling = (((image >> (HOST_BITS_PER_LONG - 2)) & 1)
2634 ^ fmt->qnan_msb_set);
2635 r->sig[SIGSZ-1] = image;
2645 r->class = rvc_normal;
2647 r->exp = exp - 127 + 1;
2648 r->sig[SIGSZ-1] = image | SIG_MSB;
2652 const struct real_format ieee_single_format =
2670 const struct real_format mips_single_format =
2689 /* IEEE double-precision format. */
2691 static void encode_ieee_double (const struct real_format *fmt,
2692 long *, const REAL_VALUE_TYPE *);
2693 static void decode_ieee_double (const struct real_format *,
2694 REAL_VALUE_TYPE *, const long *);
2697 encode_ieee_double (const struct real_format *fmt, long *buf,
2698 const REAL_VALUE_TYPE *r)
2700 unsigned long image_lo, image_hi, sig_lo, sig_hi, exp;
2701 bool denormal = (r->sig[SIGSZ-1] & SIG_MSB) == 0;
2703 image_hi = r->sign << 31;
2706 if (HOST_BITS_PER_LONG == 64)
2708 sig_hi = r->sig[SIGSZ-1];
2709 sig_lo = (sig_hi >> (64 - 53)) & 0xffffffff;
2710 sig_hi = (sig_hi >> (64 - 53 + 1) >> 31) & 0xfffff;
2714 sig_hi = r->sig[SIGSZ-1];
2715 sig_lo = r->sig[SIGSZ-2];
2716 sig_lo = (sig_hi << 21) | (sig_lo >> 11);
2717 sig_hi = (sig_hi >> 11) & 0xfffff;
2727 image_hi |= 2047 << 20;
2730 image_hi |= 0x7fffffff;
2731 image_lo = 0xffffffff;
2739 sig_hi = sig_lo = 0;
2740 if (r->signalling == fmt->qnan_msb_set)
2741 sig_hi &= ~(1 << 19);
2744 /* We overload qnan_msb_set here: it's only clear for
2745 mips_ieee_single, which wants all mantissa bits but the
2746 quiet/signalling one set in canonical NaNs (at least
2748 if (r->canonical && !fmt->qnan_msb_set)
2750 sig_hi |= (1 << 19) - 1;
2751 sig_lo = 0xffffffff;
2753 else if (sig_hi == 0 && sig_lo == 0)
2756 image_hi |= 2047 << 20;
2762 image_hi |= 0x7fffffff;
2763 image_lo = 0xffffffff;
2768 /* Recall that IEEE numbers are interpreted as 1.F x 2**exp,
2769 whereas the intermediate representation is 0.F x 2**exp.
2770 Which means we're off by one. */
2774 exp = r->exp + 1023 - 1;
2775 image_hi |= exp << 20;
2784 if (FLOAT_WORDS_BIG_ENDIAN)
2785 buf[0] = image_hi, buf[1] = image_lo;
2787 buf[0] = image_lo, buf[1] = image_hi;
2791 decode_ieee_double (const struct real_format *fmt, REAL_VALUE_TYPE *r,
2794 unsigned long image_hi, image_lo;
2798 if (FLOAT_WORDS_BIG_ENDIAN)
2799 image_hi = buf[0], image_lo = buf[1];
2801 image_lo = buf[0], image_hi = buf[1];
2802 image_lo &= 0xffffffff;
2803 image_hi &= 0xffffffff;
2805 sign = (image_hi >> 31) & 1;
2806 exp = (image_hi >> 20) & 0x7ff;
2808 memset (r, 0, sizeof (*r));
2810 image_hi <<= 32 - 21;
2811 image_hi |= image_lo >> 21;
2812 image_hi &= 0x7fffffff;
2813 image_lo <<= 32 - 21;
2817 if ((image_hi || image_lo) && fmt->has_denorm)
2819 r->class = rvc_normal;
2822 if (HOST_BITS_PER_LONG == 32)
2824 image_hi = (image_hi << 1) | (image_lo >> 31);
2826 r->sig[SIGSZ-1] = image_hi;
2827 r->sig[SIGSZ-2] = image_lo;
2831 image_hi = (image_hi << 31 << 2) | (image_lo << 1);
2832 r->sig[SIGSZ-1] = image_hi;
2836 else if (fmt->has_signed_zero)
2839 else if (exp == 2047 && (fmt->has_nans || fmt->has_inf))
2841 if (image_hi || image_lo)
2845 r->signalling = ((image_hi >> 30) & 1) ^ fmt->qnan_msb_set;
2846 if (HOST_BITS_PER_LONG == 32)
2848 r->sig[SIGSZ-1] = image_hi;
2849 r->sig[SIGSZ-2] = image_lo;
2852 r->sig[SIGSZ-1] = (image_hi << 31 << 1) | image_lo;
2862 r->class = rvc_normal;
2864 r->exp = exp - 1023 + 1;
2865 if (HOST_BITS_PER_LONG == 32)
2867 r->sig[SIGSZ-1] = image_hi | SIG_MSB;
2868 r->sig[SIGSZ-2] = image_lo;
2871 r->sig[SIGSZ-1] = (image_hi << 31 << 1) | image_lo | SIG_MSB;
2875 const struct real_format ieee_double_format =
2893 const struct real_format mips_double_format =
2912 /* IEEE extended real format. This comes in three flavors: Intel's as
2913 a 12 byte image, Intel's as a 16 byte image, and Motorola's. Intel
2914 12- and 16-byte images may be big- or little endian; Motorola's is
2915 always big endian. */
2917 /* Helper subroutine which converts from the internal format to the
2918 12-byte little-endian Intel format. Functions below adjust this
2919 for the other possible formats. */
2921 encode_ieee_extended (const struct real_format *fmt, long *buf,
2922 const REAL_VALUE_TYPE *r)
2924 unsigned long image_hi, sig_hi, sig_lo;
2925 bool denormal = (r->sig[SIGSZ-1] & SIG_MSB) == 0;
2927 image_hi = r->sign << 15;
2928 sig_hi = sig_lo = 0;
2940 /* Intel requires the explicit integer bit to be set, otherwise
2941 it considers the value a "pseudo-infinity". Motorola docs
2942 say it doesn't care. */
2943 sig_hi = 0x80000000;
2948 sig_lo = sig_hi = 0xffffffff;
2956 if (HOST_BITS_PER_LONG == 32)
2958 sig_hi = r->sig[SIGSZ-1];
2959 sig_lo = r->sig[SIGSZ-2];
2963 sig_lo = r->sig[SIGSZ-1];
2964 sig_hi = sig_lo >> 31 >> 1;
2965 sig_lo &= 0xffffffff;
2967 if (r->signalling == fmt->qnan_msb_set)
2968 sig_hi &= ~(1 << 30);
2971 if ((sig_hi & 0x7fffffff) == 0 && sig_lo == 0)
2974 /* Intel requires the explicit integer bit to be set, otherwise
2975 it considers the value a "pseudo-nan". Motorola docs say it
2977 sig_hi |= 0x80000000;
2982 sig_lo = sig_hi = 0xffffffff;
2990 /* Recall that IEEE numbers are interpreted as 1.F x 2**exp,
2991 whereas the intermediate representation is 0.F x 2**exp.
2992 Which means we're off by one.
2994 Except for Motorola, which consider exp=0 and explicit
2995 integer bit set to continue to be normalized. In theory
2996 this discrepancy has been taken care of by the difference
2997 in fmt->emin in round_for_format. */
3009 if (HOST_BITS_PER_LONG == 32)
3011 sig_hi = r->sig[SIGSZ-1];
3012 sig_lo = r->sig[SIGSZ-2];
3016 sig_lo = r->sig[SIGSZ-1];
3017 sig_hi = sig_lo >> 31 >> 1;
3018 sig_lo &= 0xffffffff;
3027 buf[0] = sig_lo, buf[1] = sig_hi, buf[2] = image_hi;
3030 /* Convert from the internal format to the 12-byte Motorola format
3031 for an IEEE extended real. */
3033 encode_ieee_extended_motorola (const struct real_format *fmt, long *buf,
3034 const REAL_VALUE_TYPE *r)
3037 encode_ieee_extended (fmt, intermed, r);
3039 /* Motorola chips are assumed always to be big-endian. Also, the
3040 padding in a Motorola extended real goes between the exponent and
3041 the mantissa. At this point the mantissa is entirely within
3042 elements 0 and 1 of intermed, and the exponent entirely within
3043 element 2, so all we have to do is swap the order around, and
3044 shift element 2 left 16 bits. */
3045 buf[0] = intermed[2] << 16;
3046 buf[1] = intermed[1];
3047 buf[2] = intermed[0];
3050 /* Convert from the internal format to the 12-byte Intel format for
3051 an IEEE extended real. */
3053 encode_ieee_extended_intel_96 (const struct real_format *fmt, long *buf,
3054 const REAL_VALUE_TYPE *r)
3056 if (FLOAT_WORDS_BIG_ENDIAN)
3058 /* All the padding in an Intel-format extended real goes at the high
3059 end, which in this case is after the mantissa, not the exponent.
3060 Therefore we must shift everything down 16 bits. */
3062 encode_ieee_extended (fmt, intermed, r);
3063 buf[0] = ((intermed[2] << 16) | ((unsigned long)(intermed[1] & 0xFFFF0000) >> 16));
3064 buf[1] = ((intermed[1] << 16) | ((unsigned long)(intermed[0] & 0xFFFF0000) >> 16));
3065 buf[2] = (intermed[0] << 16);
3068 /* encode_ieee_extended produces what we want directly. */
3069 encode_ieee_extended (fmt, buf, r);
3072 /* Convert from the internal format to the 16-byte Intel format for
3073 an IEEE extended real. */
3075 encode_ieee_extended_intel_128 (const struct real_format *fmt, long *buf,
3076 const REAL_VALUE_TYPE *r)
3078 /* All the padding in an Intel-format extended real goes at the high end. */
3079 encode_ieee_extended_intel_96 (fmt, buf, r);
3083 /* As above, we have a helper function which converts from 12-byte
3084 little-endian Intel format to internal format. Functions below
3085 adjust for the other possible formats. */
3087 decode_ieee_extended (const struct real_format *fmt, REAL_VALUE_TYPE *r,
3090 unsigned long image_hi, sig_hi, sig_lo;
3094 sig_lo = buf[0], sig_hi = buf[1], image_hi = buf[2];
3095 sig_lo &= 0xffffffff;
3096 sig_hi &= 0xffffffff;
3097 image_hi &= 0xffffffff;
3099 sign = (image_hi >> 15) & 1;
3100 exp = image_hi & 0x7fff;
3102 memset (r, 0, sizeof (*r));
3106 if ((sig_hi || sig_lo) && fmt->has_denorm)
3108 r->class = rvc_normal;
3111 /* When the IEEE format contains a hidden bit, we know that
3112 it's zero at this point, and so shift up the significand
3113 and decrease the exponent to match. In this case, Motorola
3114 defines the explicit integer bit to be valid, so we don't
3115 know whether the msb is set or not. */
3117 if (HOST_BITS_PER_LONG == 32)
3119 r->sig[SIGSZ-1] = sig_hi;
3120 r->sig[SIGSZ-2] = sig_lo;
3123 r->sig[SIGSZ-1] = (sig_hi << 31 << 1) | sig_lo;
3127 else if (fmt->has_signed_zero)
3130 else if (exp == 32767 && (fmt->has_nans || fmt->has_inf))
3132 /* See above re "pseudo-infinities" and "pseudo-nans".
3133 Short summary is that the MSB will likely always be
3134 set, and that we don't care about it. */
3135 sig_hi &= 0x7fffffff;
3137 if (sig_hi || sig_lo)
3141 r->signalling = ((sig_hi >> 30) & 1) ^ fmt->qnan_msb_set;
3142 if (HOST_BITS_PER_LONG == 32)
3144 r->sig[SIGSZ-1] = sig_hi;
3145 r->sig[SIGSZ-2] = sig_lo;
3148 r->sig[SIGSZ-1] = (sig_hi << 31 << 1) | sig_lo;
3158 r->class = rvc_normal;
3160 r->exp = exp - 16383 + 1;
3161 if (HOST_BITS_PER_LONG == 32)
3163 r->sig[SIGSZ-1] = sig_hi;
3164 r->sig[SIGSZ-2] = sig_lo;
3167 r->sig[SIGSZ-1] = (sig_hi << 31 << 1) | sig_lo;
3171 /* Convert from the internal format to the 12-byte Motorola format
3172 for an IEEE extended real. */
3174 decode_ieee_extended_motorola (const struct real_format *fmt, REAL_VALUE_TYPE *r,
3179 /* Motorola chips are assumed always to be big-endian. Also, the
3180 padding in a Motorola extended real goes between the exponent and
3181 the mantissa; remove it. */
3182 intermed[0] = buf[2];
3183 intermed[1] = buf[1];
3184 intermed[2] = (unsigned long)buf[0] >> 16;
3186 decode_ieee_extended (fmt, r, intermed);
3189 /* Convert from the internal format to the 12-byte Intel format for
3190 an IEEE extended real. */
3192 decode_ieee_extended_intel_96 (const struct real_format *fmt, REAL_VALUE_TYPE *r,
3195 if (FLOAT_WORDS_BIG_ENDIAN)
3197 /* All the padding in an Intel-format extended real goes at the high
3198 end, which in this case is after the mantissa, not the exponent.
3199 Therefore we must shift everything up 16 bits. */
3202 intermed[0] = (((unsigned long)buf[2] >> 16) | (buf[1] << 16));
3203 intermed[1] = (((unsigned long)buf[1] >> 16) | (buf[0] << 16));
3204 intermed[2] = ((unsigned long)buf[0] >> 16);
3206 decode_ieee_extended (fmt, r, intermed);
3209 /* decode_ieee_extended produces what we want directly. */
3210 decode_ieee_extended (fmt, r, buf);
3213 /* Convert from the internal format to the 16-byte Intel format for
3214 an IEEE extended real. */
3216 decode_ieee_extended_intel_128 (const struct real_format *fmt, REAL_VALUE_TYPE *r,
3219 /* All the padding in an Intel-format extended real goes at the high end. */
3220 decode_ieee_extended_intel_96 (fmt, r, buf);
3223 const struct real_format ieee_extended_motorola_format =
3225 encode_ieee_extended_motorola,
3226 decode_ieee_extended_motorola,
3241 const struct real_format ieee_extended_intel_96_format =
3243 encode_ieee_extended_intel_96,
3244 decode_ieee_extended_intel_96,
3259 const struct real_format ieee_extended_intel_128_format =
3261 encode_ieee_extended_intel_128,
3262 decode_ieee_extended_intel_128,
3277 /* The following caters to i386 systems that set the rounding precision
3278 to 53 bits instead of 64, e.g. FreeBSD. */
3279 const struct real_format ieee_extended_intel_96_round_53_format =
3281 encode_ieee_extended_intel_96,
3282 decode_ieee_extended_intel_96,
3297 /* IBM 128-bit extended precision format: a pair of IEEE double precision
3298 numbers whose sum is equal to the extended precision value. The number
3299 with greater magnitude is first. This format has the same magnitude
3300 range as an IEEE double precision value, but effectively 106 bits of
3301 significand precision. Infinity and NaN are represented by their IEEE
3302 double precision value stored in the first number, the second number is
3303 ignored. Zeroes, Infinities, and NaNs are set in both doubles
3304 due to precedent. */
3306 static void encode_ibm_extended (const struct real_format *fmt,
3307 long *, const REAL_VALUE_TYPE *);
3308 static void decode_ibm_extended (const struct real_format *,
3309 REAL_VALUE_TYPE *, const long *);
3312 encode_ibm_extended (const struct real_format *fmt, long *buf,
3313 const REAL_VALUE_TYPE *r)
3315 REAL_VALUE_TYPE u, normr, v;
3316 const struct real_format *base_fmt;
3318 base_fmt = fmt->qnan_msb_set ? &ieee_double_format : &mips_double_format;
3320 /* Renormlize R before doing any arithmetic on it. */
3322 if (normr.class == rvc_normal)
3325 /* u = IEEE double precision portion of significand. */
3327 round_for_format (base_fmt, &u);
3328 encode_ieee_double (base_fmt, &buf[0], &u);
3330 if (u.class == rvc_normal)
3332 do_add (&v, &normr, &u, 1);
3333 /* Call round_for_format since we might need to denormalize. */
3334 round_for_format (base_fmt, &v);
3335 encode_ieee_double (base_fmt, &buf[2], &v);
3339 /* Inf, NaN, 0 are all representable as doubles, so the
3340 least-significant part can be 0.0. */
3347 decode_ibm_extended (const struct real_format *fmt ATTRIBUTE_UNUSED, REAL_VALUE_TYPE *r,
3350 REAL_VALUE_TYPE u, v;
3351 const struct real_format *base_fmt;
3353 base_fmt = fmt->qnan_msb_set ? &ieee_double_format : &mips_double_format;
3354 decode_ieee_double (base_fmt, &u, &buf[0]);
3356 if (u.class != rvc_zero && u.class != rvc_inf && u.class != rvc_nan)
3358 decode_ieee_double (base_fmt, &v, &buf[2]);
3359 do_add (r, &u, &v, 0);
3365 const struct real_format ibm_extended_format =
3367 encode_ibm_extended,
3368 decode_ibm_extended,
3383 const struct real_format mips_extended_format =
3385 encode_ibm_extended,
3386 decode_ibm_extended,
3402 /* IEEE quad precision format. */
3404 static void encode_ieee_quad (const struct real_format *fmt,
3405 long *, const REAL_VALUE_TYPE *);
3406 static void decode_ieee_quad (const struct real_format *,
3407 REAL_VALUE_TYPE *, const long *);
3410 encode_ieee_quad (const struct real_format *fmt, long *buf,
3411 const REAL_VALUE_TYPE *r)
3413 unsigned long image3, image2, image1, image0, exp;
3414 bool denormal = (r->sig[SIGSZ-1] & SIG_MSB) == 0;
3417 image3 = r->sign << 31;
3422 rshift_significand (&u, r, SIGNIFICAND_BITS - 113);
3431 image3 |= 32767 << 16;
3434 image3 |= 0x7fffffff;
3435 image2 = 0xffffffff;
3436 image1 = 0xffffffff;
3437 image0 = 0xffffffff;
3444 image3 |= 32767 << 16;
3448 /* Don't use bits from the significand. The
3449 initialization above is right. */
3451 else if (HOST_BITS_PER_LONG == 32)
3456 image3 |= u.sig[3] & 0xffff;
3461 image1 = image0 >> 31 >> 1;
3463 image3 |= (image2 >> 31 >> 1) & 0xffff;
3464 image0 &= 0xffffffff;
3465 image2 &= 0xffffffff;
3467 if (r->signalling == fmt->qnan_msb_set)
3471 /* We overload qnan_msb_set here: it's only clear for
3472 mips_ieee_single, which wants all mantissa bits but the
3473 quiet/signalling one set in canonical NaNs (at least
3475 if (r->canonical && !fmt->qnan_msb_set)
3478 image2 = image1 = image0 = 0xffffffff;
3480 else if (((image3 & 0xffff) | image2 | image1 | image0) == 0)
3485 image3 |= 0x7fffffff;
3486 image2 = 0xffffffff;
3487 image1 = 0xffffffff;
3488 image0 = 0xffffffff;
3493 /* Recall that IEEE numbers are interpreted as 1.F x 2**exp,
3494 whereas the intermediate representation is 0.F x 2**exp.
3495 Which means we're off by one. */
3499 exp = r->exp + 16383 - 1;
3500 image3 |= exp << 16;
3502 if (HOST_BITS_PER_LONG == 32)
3507 image3 |= u.sig[3] & 0xffff;
3512 image1 = image0 >> 31 >> 1;
3514 image3 |= (image2 >> 31 >> 1) & 0xffff;
3515 image0 &= 0xffffffff;
3516 image2 &= 0xffffffff;
3524 if (FLOAT_WORDS_BIG_ENDIAN)
3541 decode_ieee_quad (const struct real_format *fmt, REAL_VALUE_TYPE *r,
3544 unsigned long image3, image2, image1, image0;
3548 if (FLOAT_WORDS_BIG_ENDIAN)
3562 image0 &= 0xffffffff;
3563 image1 &= 0xffffffff;
3564 image2 &= 0xffffffff;
3566 sign = (image3 >> 31) & 1;
3567 exp = (image3 >> 16) & 0x7fff;
3570 memset (r, 0, sizeof (*r));
3574 if ((image3 | image2 | image1 | image0) && fmt->has_denorm)
3576 r->class = rvc_normal;
3579 r->exp = -16382 + (SIGNIFICAND_BITS - 112);
3580 if (HOST_BITS_PER_LONG == 32)
3589 r->sig[0] = (image1 << 31 << 1) | image0;
3590 r->sig[1] = (image3 << 31 << 1) | image2;
3595 else if (fmt->has_signed_zero)
3598 else if (exp == 32767 && (fmt->has_nans || fmt->has_inf))
3600 if (image3 | image2 | image1 | image0)
3604 r->signalling = ((image3 >> 15) & 1) ^ fmt->qnan_msb_set;
3606 if (HOST_BITS_PER_LONG == 32)
3615 r->sig[0] = (image1 << 31 << 1) | image0;
3616 r->sig[1] = (image3 << 31 << 1) | image2;
3618 lshift_significand (r, r, SIGNIFICAND_BITS - 113);
3628 r->class = rvc_normal;
3630 r->exp = exp - 16383 + 1;
3632 if (HOST_BITS_PER_LONG == 32)
3641 r->sig[0] = (image1 << 31 << 1) | image0;
3642 r->sig[1] = (image3 << 31 << 1) | image2;
3644 lshift_significand (r, r, SIGNIFICAND_BITS - 113);
3645 r->sig[SIGSZ-1] |= SIG_MSB;
3649 const struct real_format ieee_quad_format =
3667 const struct real_format mips_quad_format =
3685 /* Descriptions of VAX floating point formats can be found beginning at
3687 http://h71000.www7.hp.com/doc/73FINAL/4515/4515pro_013.html#f_floating_point_format
3689 The thing to remember is that they're almost IEEE, except for word
3690 order, exponent bias, and the lack of infinities, nans, and denormals.
3692 We don't implement the H_floating format here, simply because neither
3693 the VAX or Alpha ports use it. */
3695 static void encode_vax_f (const struct real_format *fmt,
3696 long *, const REAL_VALUE_TYPE *);
3697 static void decode_vax_f (const struct real_format *,
3698 REAL_VALUE_TYPE *, const long *);
3699 static void encode_vax_d (const struct real_format *fmt,
3700 long *, const REAL_VALUE_TYPE *);
3701 static void decode_vax_d (const struct real_format *,
3702 REAL_VALUE_TYPE *, const long *);
3703 static void encode_vax_g (const struct real_format *fmt,
3704 long *, const REAL_VALUE_TYPE *);
3705 static void decode_vax_g (const struct real_format *,
3706 REAL_VALUE_TYPE *, const long *);
3709 encode_vax_f (const struct real_format *fmt ATTRIBUTE_UNUSED, long *buf,
3710 const REAL_VALUE_TYPE *r)
3712 unsigned long sign, exp, sig, image;
3714 sign = r->sign << 15;
3724 image = 0xffff7fff | sign;
3728 sig = (r->sig[SIGSZ-1] >> (HOST_BITS_PER_LONG - 24)) & 0x7fffff;
3731 image = (sig << 16) & 0xffff0000;
3745 decode_vax_f (const struct real_format *fmt ATTRIBUTE_UNUSED,
3746 REAL_VALUE_TYPE *r, const long *buf)
3748 unsigned long image = buf[0] & 0xffffffff;
3749 int exp = (image >> 7) & 0xff;
3751 memset (r, 0, sizeof (*r));
3755 r->class = rvc_normal;
3756 r->sign = (image >> 15) & 1;
3759 image = ((image & 0x7f) << 16) | ((image >> 16) & 0xffff);
3760 r->sig[SIGSZ-1] = (image << (HOST_BITS_PER_LONG - 24)) | SIG_MSB;
3765 encode_vax_d (const struct real_format *fmt ATTRIBUTE_UNUSED, long *buf,
3766 const REAL_VALUE_TYPE *r)
3768 unsigned long image0, image1, sign = r->sign << 15;
3773 image0 = image1 = 0;
3778 image0 = 0xffff7fff | sign;
3779 image1 = 0xffffffff;
3783 /* Extract the significand into straight hi:lo. */
3784 if (HOST_BITS_PER_LONG == 64)
3786 image0 = r->sig[SIGSZ-1];
3787 image1 = (image0 >> (64 - 56)) & 0xffffffff;
3788 image0 = (image0 >> (64 - 56 + 1) >> 31) & 0x7fffff;
3792 image0 = r->sig[SIGSZ-1];
3793 image1 = r->sig[SIGSZ-2];
3794 image1 = (image0 << 24) | (image1 >> 8);
3795 image0 = (image0 >> 8) & 0xffffff;
3798 /* Rearrange the half-words of the significand to match the
3800 image0 = ((image0 << 16) | (image0 >> 16)) & 0xffff007f;
3801 image1 = ((image1 << 16) | (image1 >> 16)) & 0xffffffff;
3803 /* Add the sign and exponent. */
3805 image0 |= (r->exp + 128) << 7;
3812 if (FLOAT_WORDS_BIG_ENDIAN)
3813 buf[0] = image1, buf[1] = image0;
3815 buf[0] = image0, buf[1] = image1;
3819 decode_vax_d (const struct real_format *fmt ATTRIBUTE_UNUSED,
3820 REAL_VALUE_TYPE *r, const long *buf)
3822 unsigned long image0, image1;
3825 if (FLOAT_WORDS_BIG_ENDIAN)
3826 image1 = buf[0], image0 = buf[1];
3828 image0 = buf[0], image1 = buf[1];
3829 image0 &= 0xffffffff;
3830 image1 &= 0xffffffff;
3832 exp = (image0 >> 7) & 0xff;
3834 memset (r, 0, sizeof (*r));
3838 r->class = rvc_normal;
3839 r->sign = (image0 >> 15) & 1;
3842 /* Rearrange the half-words of the external format into
3843 proper ascending order. */
3844 image0 = ((image0 & 0x7f) << 16) | ((image0 >> 16) & 0xffff);
3845 image1 = ((image1 & 0xffff) << 16) | ((image1 >> 16) & 0xffff);
3847 if (HOST_BITS_PER_LONG == 64)
3849 image0 = (image0 << 31 << 1) | image1;
3852 r->sig[SIGSZ-1] = image0;
3856 r->sig[SIGSZ-1] = image0;
3857 r->sig[SIGSZ-2] = image1;
3858 lshift_significand (r, r, 2*HOST_BITS_PER_LONG - 56);
3859 r->sig[SIGSZ-1] |= SIG_MSB;
3865 encode_vax_g (const struct real_format *fmt ATTRIBUTE_UNUSED, long *buf,
3866 const REAL_VALUE_TYPE *r)
3868 unsigned long image0, image1, sign = r->sign << 15;
3873 image0 = image1 = 0;
3878 image0 = 0xffff7fff | sign;
3879 image1 = 0xffffffff;
3883 /* Extract the significand into straight hi:lo. */
3884 if (HOST_BITS_PER_LONG == 64)
3886 image0 = r->sig[SIGSZ-1];
3887 image1 = (image0 >> (64 - 53)) & 0xffffffff;
3888 image0 = (image0 >> (64 - 53 + 1) >> 31) & 0xfffff;
3892 image0 = r->sig[SIGSZ-1];
3893 image1 = r->sig[SIGSZ-2];
3894 image1 = (image0 << 21) | (image1 >> 11);
3895 image0 = (image0 >> 11) & 0xfffff;
3898 /* Rearrange the half-words of the significand to match the
3900 image0 = ((image0 << 16) | (image0 >> 16)) & 0xffff000f;
3901 image1 = ((image1 << 16) | (image1 >> 16)) & 0xffffffff;
3903 /* Add the sign and exponent. */
3905 image0 |= (r->exp + 1024) << 4;
3912 if (FLOAT_WORDS_BIG_ENDIAN)
3913 buf[0] = image1, buf[1] = image0;
3915 buf[0] = image0, buf[1] = image1;
3919 decode_vax_g (const struct real_format *fmt ATTRIBUTE_UNUSED,
3920 REAL_VALUE_TYPE *r, const long *buf)
3922 unsigned long image0, image1;
3925 if (FLOAT_WORDS_BIG_ENDIAN)
3926 image1 = buf[0], image0 = buf[1];
3928 image0 = buf[0], image1 = buf[1];
3929 image0 &= 0xffffffff;
3930 image1 &= 0xffffffff;
3932 exp = (image0 >> 4) & 0x7ff;
3934 memset (r, 0, sizeof (*r));
3938 r->class = rvc_normal;
3939 r->sign = (image0 >> 15) & 1;
3940 r->exp = exp - 1024;
3942 /* Rearrange the half-words of the external format into
3943 proper ascending order. */
3944 image0 = ((image0 & 0xf) << 16) | ((image0 >> 16) & 0xffff);
3945 image1 = ((image1 & 0xffff) << 16) | ((image1 >> 16) & 0xffff);
3947 if (HOST_BITS_PER_LONG == 64)
3949 image0 = (image0 << 31 << 1) | image1;
3952 r->sig[SIGSZ-1] = image0;
3956 r->sig[SIGSZ-1] = image0;
3957 r->sig[SIGSZ-2] = image1;
3958 lshift_significand (r, r, 64 - 53);
3959 r->sig[SIGSZ-1] |= SIG_MSB;
3964 const struct real_format vax_f_format =
3982 const struct real_format vax_d_format =
4000 const struct real_format vax_g_format =
4018 /* A good reference for these can be found in chapter 9 of
4019 "ESA/390 Principles of Operation", IBM document number SA22-7201-01.
4020 An on-line version can be found here:
4022 http://publibz.boulder.ibm.com/cgi-bin/bookmgr_OS390/BOOKS/DZ9AR001/9.1?DT=19930923083613
4025 static void encode_i370_single (const struct real_format *fmt,
4026 long *, const REAL_VALUE_TYPE *);
4027 static void decode_i370_single (const struct real_format *,
4028 REAL_VALUE_TYPE *, const long *);
4029 static void encode_i370_double (const struct real_format *fmt,
4030 long *, const REAL_VALUE_TYPE *);
4031 static void decode_i370_double (const struct real_format *,
4032 REAL_VALUE_TYPE *, const long *);
4035 encode_i370_single (const struct real_format *fmt ATTRIBUTE_UNUSED,
4036 long *buf, const REAL_VALUE_TYPE *r)
4038 unsigned long sign, exp, sig, image;
4040 sign = r->sign << 31;
4050 image = 0x7fffffff | sign;
4054 sig = (r->sig[SIGSZ-1] >> (HOST_BITS_PER_LONG - 24)) & 0xffffff;
4055 exp = ((r->exp / 4) + 64) << 24;
4056 image = sign | exp | sig;
4067 decode_i370_single (const struct real_format *fmt ATTRIBUTE_UNUSED,
4068 REAL_VALUE_TYPE *r, const long *buf)
4070 unsigned long sign, sig, image = buf[0];
4073 sign = (image >> 31) & 1;
4074 exp = (image >> 24) & 0x7f;
4075 sig = image & 0xffffff;
4077 memset (r, 0, sizeof (*r));
4081 r->class = rvc_normal;
4083 r->exp = (exp - 64) * 4;
4084 r->sig[SIGSZ-1] = sig << (HOST_BITS_PER_LONG - 24);
4090 encode_i370_double (const struct real_format *fmt ATTRIBUTE_UNUSED,
4091 long *buf, const REAL_VALUE_TYPE *r)
4093 unsigned long sign, exp, image_hi, image_lo;
4095 sign = r->sign << 31;
4100 image_hi = image_lo = 0;
4105 image_hi = 0x7fffffff | sign;
4106 image_lo = 0xffffffff;
4110 if (HOST_BITS_PER_LONG == 64)
4112 image_hi = r->sig[SIGSZ-1];
4113 image_lo = (image_hi >> (64 - 56)) & 0xffffffff;
4114 image_hi = (image_hi >> (64 - 56 + 1) >> 31) & 0xffffff;
4118 image_hi = r->sig[SIGSZ-1];
4119 image_lo = r->sig[SIGSZ-2];
4120 image_lo = (image_lo >> 8) | (image_hi << 24);
4124 exp = ((r->exp / 4) + 64) << 24;
4125 image_hi |= sign | exp;
4132 if (FLOAT_WORDS_BIG_ENDIAN)
4133 buf[0] = image_hi, buf[1] = image_lo;
4135 buf[0] = image_lo, buf[1] = image_hi;
4139 decode_i370_double (const struct real_format *fmt ATTRIBUTE_UNUSED,
4140 REAL_VALUE_TYPE *r, const long *buf)
4142 unsigned long sign, image_hi, image_lo;
4145 if (FLOAT_WORDS_BIG_ENDIAN)
4146 image_hi = buf[0], image_lo = buf[1];
4148 image_lo = buf[0], image_hi = buf[1];
4150 sign = (image_hi >> 31) & 1;
4151 exp = (image_hi >> 24) & 0x7f;
4152 image_hi &= 0xffffff;
4153 image_lo &= 0xffffffff;
4155 memset (r, 0, sizeof (*r));
4157 if (exp || image_hi || image_lo)
4159 r->class = rvc_normal;
4161 r->exp = (exp - 64) * 4 + (SIGNIFICAND_BITS - 56);
4163 if (HOST_BITS_PER_LONG == 32)
4165 r->sig[0] = image_lo;
4166 r->sig[1] = image_hi;
4169 r->sig[0] = image_lo | (image_hi << 31 << 1);
4175 const struct real_format i370_single_format =
4188 false, /* ??? The encoding does allow for "unnormals". */
4189 false, /* ??? The encoding does allow for "unnormals". */
4193 const struct real_format i370_double_format =
4206 false, /* ??? The encoding does allow for "unnormals". */
4207 false, /* ??? The encoding does allow for "unnormals". */
4211 /* The "twos-complement" c4x format is officially defined as
4215 This is rather misleading. One must remember that F is signed.
4216 A better description would be
4218 x = -1**s * ((s + 1 + .f) * 2**e
4220 So if we have a (4 bit) fraction of .1000 with a sign bit of 1,
4221 that's -1 * (1+1+(-.5)) == -1.5. I think.
4223 The constructions here are taken from Tables 5-1 and 5-2 of the
4224 TMS320C4x User's Guide wherein step-by-step instructions for
4225 conversion from IEEE are presented. That's close enough to our
4226 internal representation so as to make things easy.
4228 See http://www-s.ti.com/sc/psheets/spru063c/spru063c.pdf */
4230 static void encode_c4x_single (const struct real_format *fmt,
4231 long *, const REAL_VALUE_TYPE *);
4232 static void decode_c4x_single (const struct real_format *,
4233 REAL_VALUE_TYPE *, const long *);
4234 static void encode_c4x_extended (const struct real_format *fmt,
4235 long *, const REAL_VALUE_TYPE *);
4236 static void decode_c4x_extended (const struct real_format *,
4237 REAL_VALUE_TYPE *, const long *);
4240 encode_c4x_single (const struct real_format *fmt ATTRIBUTE_UNUSED,
4241 long *buf, const REAL_VALUE_TYPE *r)
4243 unsigned long image, exp, sig;
4255 sig = 0x800000 - r->sign;
4260 sig = (r->sig[SIGSZ-1] >> (HOST_BITS_PER_LONG - 24)) & 0x7fffff;
4275 image = ((exp & 0xff) << 24) | (sig & 0xffffff);
4280 decode_c4x_single (const struct real_format *fmt ATTRIBUTE_UNUSED,
4281 REAL_VALUE_TYPE *r, const long *buf)
4283 unsigned long image = buf[0];
4287 exp = (((image >> 24) & 0xff) ^ 0x80) - 0x80;
4288 sf = ((image & 0xffffff) ^ 0x800000) - 0x800000;
4290 memset (r, 0, sizeof (*r));
4294 r->class = rvc_normal;
4296 sig = sf & 0x7fffff;
4305 sig = (sig << (HOST_BITS_PER_LONG - 24)) | SIG_MSB;
4308 r->sig[SIGSZ-1] = sig;
4313 encode_c4x_extended (const struct real_format *fmt ATTRIBUTE_UNUSED,
4314 long *buf, const REAL_VALUE_TYPE *r)
4316 unsigned long exp, sig;
4328 sig = 0x80000000 - r->sign;
4334 sig = r->sig[SIGSZ-1];
4335 if (HOST_BITS_PER_LONG == 64)
4336 sig = sig >> 1 >> 31;
4353 exp = (exp & 0xff) << 24;
4356 if (FLOAT_WORDS_BIG_ENDIAN)
4357 buf[0] = exp, buf[1] = sig;
4359 buf[0] = sig, buf[0] = exp;
4363 decode_c4x_extended (const struct real_format *fmt ATTRIBUTE_UNUSED,
4364 REAL_VALUE_TYPE *r, const long *buf)
4369 if (FLOAT_WORDS_BIG_ENDIAN)
4370 exp = buf[0], sf = buf[1];
4372 sf = buf[0], exp = buf[1];
4374 exp = (((exp >> 24) & 0xff) & 0x80) - 0x80;
4375 sf = ((sf & 0xffffffff) ^ 0x80000000) - 0x80000000;
4377 memset (r, 0, sizeof (*r));
4381 r->class = rvc_normal;
4383 sig = sf & 0x7fffffff;
4392 if (HOST_BITS_PER_LONG == 64)
4393 sig = sig << 1 << 31;
4397 r->sig[SIGSZ-1] = sig;
4401 const struct real_format c4x_single_format =
4419 const struct real_format c4x_extended_format =
4421 encode_c4x_extended,
4422 decode_c4x_extended,
4438 /* A synthetic "format" for internal arithmetic. It's the size of the
4439 internal significand minus the two bits needed for proper rounding.
4440 The encode and decode routines exist only to satisfy our paranoia
4443 static void encode_internal (const struct real_format *fmt,
4444 long *, const REAL_VALUE_TYPE *);
4445 static void decode_internal (const struct real_format *,
4446 REAL_VALUE_TYPE *, const long *);
4449 encode_internal (const struct real_format *fmt ATTRIBUTE_UNUSED, long *buf,
4450 const REAL_VALUE_TYPE *r)
4452 memcpy (buf, r, sizeof (*r));
4456 decode_internal (const struct real_format *fmt ATTRIBUTE_UNUSED,
4457 REAL_VALUE_TYPE *r, const long *buf)
4459 memcpy (r, buf, sizeof (*r));
4462 const struct real_format real_internal_format =
4468 SIGNIFICAND_BITS - 2,
4469 SIGNIFICAND_BITS - 2,
4480 /* Calculate the square root of X in mode MODE, and store the result
4481 in R. Return TRUE if the operation does not raise an exception.
4482 For details see "High Precision Division and Square Root",
4483 Alan H. Karp and Peter Markstein, HP Lab Report 93-93-42, June
4484 1993. http://www.hpl.hp.com/techreports/93/HPL-93-42.pdf. */
4487 real_sqrt (REAL_VALUE_TYPE *r, enum machine_mode mode,
4488 const REAL_VALUE_TYPE *x)
4490 static REAL_VALUE_TYPE halfthree;
4491 static bool init = false;
4492 REAL_VALUE_TYPE h, t, i;
4495 /* sqrt(-0.0) is -0.0. */
4496 if (real_isnegzero (x))
4502 /* Negative arguments return NaN. */
4505 get_canonical_qnan (r, 0);
4509 /* Infinity and NaN return themselves. */
4510 if (real_isinf (x) || real_isnan (x))
4518 do_add (&halfthree, &dconst1, &dconsthalf, 0);
4522 /* Initial guess for reciprocal sqrt, i. */
4523 exp = real_exponent (x);
4524 real_ldexp (&i, &dconst1, -exp/2);
4526 /* Newton's iteration for reciprocal sqrt, i. */
4527 for (iter = 0; iter < 16; iter++)
4529 /* i(n+1) = i(n) * (1.5 - 0.5*i(n)*i(n)*x). */
4530 do_multiply (&t, x, &i);
4531 do_multiply (&h, &t, &i);
4532 do_multiply (&t, &h, &dconsthalf);
4533 do_add (&h, &halfthree, &t, 1);
4534 do_multiply (&t, &i, &h);
4536 /* Check for early convergence. */
4537 if (iter >= 6 && real_identical (&i, &t))
4540 /* ??? Unroll loop to avoid copying. */
4544 /* Final iteration: r = i*x + 0.5*i*x*(1.0 - i*(i*x)). */
4545 do_multiply (&t, x, &i);
4546 do_multiply (&h, &t, &i);
4547 do_add (&i, &dconst1, &h, 1);
4548 do_multiply (&h, &t, &i);
4549 do_multiply (&i, &dconsthalf, &h);
4550 do_add (&h, &t, &i, 0);
4552 /* ??? We need a Tuckerman test to get the last bit. */
4554 real_convert (r, mode, &h);
4558 /* Calculate X raised to the integer exponent N in mode MODE and store
4559 the result in R. Return true if the result may be inexact due to
4560 loss of precision. The algorithm is the classic "left-to-right binary
4561 method" described in section 4.6.3 of Donald Knuth's "Seminumerical
4562 Algorithms", "The Art of Computer Programming", Volume 2. */
4565 real_powi (REAL_VALUE_TYPE *r, enum machine_mode mode,
4566 const REAL_VALUE_TYPE *x, HOST_WIDE_INT n)
4568 unsigned HOST_WIDE_INT bit;
4570 bool inexact = false;
4582 /* Don't worry about overflow, from now on n is unsigned. */
4590 bit = (unsigned HOST_WIDE_INT) 1 << (HOST_BITS_PER_WIDE_INT - 1);
4591 for (i = 0; i < HOST_BITS_PER_WIDE_INT; i++)
4595 inexact |= do_multiply (&t, &t, &t);
4597 inexact |= do_multiply (&t, &t, x);
4605 inexact |= do_divide (&t, &dconst1, &t);
4607 real_convert (r, mode, &t);
4611 /* Round X to the nearest integer not larger in absolute value, i.e.
4612 towards zero, placing the result in R in mode MODE. */
4615 real_trunc (REAL_VALUE_TYPE *r, enum machine_mode mode,
4616 const REAL_VALUE_TYPE *x)
4618 do_fix_trunc (r, x);
4619 if (mode != VOIDmode)
4620 real_convert (r, mode, r);
4623 /* Round X to the largest integer not greater in value, i.e. round
4624 down, placing the result in R in mode MODE. */
4627 real_floor (REAL_VALUE_TYPE *r, enum machine_mode mode,
4628 const REAL_VALUE_TYPE *x)
4630 do_fix_trunc (r, x);
4631 if (! real_identical (r, x) && r->sign)
4632 do_add (r, r, &dconstm1, 0);
4633 if (mode != VOIDmode)
4634 real_convert (r, mode, r);
4637 /* Round X to the smallest integer not less then argument, i.e. round
4638 up, placing the result in R in mode MODE. */
4641 real_ceil (REAL_VALUE_TYPE *r, enum machine_mode mode,
4642 const REAL_VALUE_TYPE *x)
4644 do_fix_trunc (r, x);
4645 if (! real_identical (r, x) && ! r->sign)
4646 do_add (r, r, &dconst1, 0);
4647 if (mode != VOIDmode)
4648 real_convert (r, mode, r);