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 bool do_add PARAMS ((REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *,
119 const REAL_VALUE_TYPE *, int));
120 static bool do_multiply PARAMS ((REAL_VALUE_TYPE *,
121 const REAL_VALUE_TYPE *,
122 const REAL_VALUE_TYPE *));
123 static bool 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));
165 get_canonical_snan (r, sign)
169 memset (r, 0, sizeof (*r));
181 memset (r, 0, sizeof (*r));
187 /* Right-shift the significand of A by N bits; put the result in the
188 significand of R. If any one bits are shifted out, return true. */
191 sticky_rshift_significand (r, a, n)
193 const REAL_VALUE_TYPE *a;
196 unsigned long sticky = 0;
197 unsigned int i, ofs = 0;
199 if (n >= HOST_BITS_PER_LONG)
201 for (i = 0, ofs = n / HOST_BITS_PER_LONG; i < ofs; ++i)
203 n &= HOST_BITS_PER_LONG - 1;
208 sticky |= a->sig[ofs] & (((unsigned long)1 << n) - 1);
209 for (i = 0; i < SIGSZ; ++i)
212 = (((ofs + i >= SIGSZ ? 0 : a->sig[ofs + i]) >> n)
213 | ((ofs + i + 1 >= SIGSZ ? 0 : a->sig[ofs + i + 1])
214 << (HOST_BITS_PER_LONG - n)));
219 for (i = 0; ofs + i < SIGSZ; ++i)
220 r->sig[i] = a->sig[ofs + i];
221 for (; i < SIGSZ; ++i)
228 /* Right-shift the significand of A by N bits; put the result in the
232 rshift_significand (r, a, n)
234 const REAL_VALUE_TYPE *a;
237 unsigned int i, ofs = n / HOST_BITS_PER_LONG;
239 n &= HOST_BITS_PER_LONG - 1;
242 for (i = 0; i < SIGSZ; ++i)
245 = (((ofs + i >= SIGSZ ? 0 : a->sig[ofs + i]) >> n)
246 | ((ofs + i + 1 >= SIGSZ ? 0 : a->sig[ofs + i + 1])
247 << (HOST_BITS_PER_LONG - n)));
252 for (i = 0; ofs + i < SIGSZ; ++i)
253 r->sig[i] = a->sig[ofs + i];
254 for (; i < SIGSZ; ++i)
259 /* Left-shift the significand of A by N bits; put the result in the
263 lshift_significand (r, a, n)
265 const REAL_VALUE_TYPE *a;
268 unsigned int i, ofs = n / HOST_BITS_PER_LONG;
270 n &= HOST_BITS_PER_LONG - 1;
273 for (i = 0; ofs + i < SIGSZ; ++i)
274 r->sig[SIGSZ-1-i] = a->sig[SIGSZ-1-i-ofs];
275 for (; i < SIGSZ; ++i)
276 r->sig[SIGSZ-1-i] = 0;
279 for (i = 0; i < SIGSZ; ++i)
282 = (((ofs + i >= SIGSZ ? 0 : a->sig[SIGSZ-1-i-ofs]) << n)
283 | ((ofs + i + 1 >= SIGSZ ? 0 : a->sig[SIGSZ-1-i-ofs-1])
284 >> (HOST_BITS_PER_LONG - n)));
288 /* Likewise, but N is specialized to 1. */
291 lshift_significand_1 (r, a)
293 const REAL_VALUE_TYPE *a;
297 for (i = SIGSZ - 1; i > 0; --i)
298 r->sig[i] = (a->sig[i] << 1) | (a->sig[i-1] >> (HOST_BITS_PER_LONG - 1));
299 r->sig[0] = a->sig[0] << 1;
302 /* Add the significands of A and B, placing the result in R. Return
303 true if there was carry out of the most significant word. */
306 add_significands (r, a, b)
308 const REAL_VALUE_TYPE *a, *b;
313 for (i = 0; i < SIGSZ; ++i)
315 unsigned long ai = a->sig[i];
316 unsigned long ri = ai + b->sig[i];
332 /* Subtract the significands of A and B, placing the result in R. CARRY is
333 true if there's a borrow incoming to the least significant word.
334 Return true if there was borrow out of the most significant word. */
337 sub_significands (r, a, b, carry)
339 const REAL_VALUE_TYPE *a, *b;
344 for (i = 0; i < SIGSZ; ++i)
346 unsigned long ai = a->sig[i];
347 unsigned long ri = ai - b->sig[i];
363 /* Negate the significand A, placing the result in R. */
366 neg_significand (r, a)
368 const REAL_VALUE_TYPE *a;
373 for (i = 0; i < SIGSZ; ++i)
375 unsigned long ri, ai = a->sig[i];
394 /* Compare significands. Return tri-state vs zero. */
397 cmp_significands (a, b)
398 const REAL_VALUE_TYPE *a, *b;
402 for (i = SIGSZ - 1; i >= 0; --i)
404 unsigned long ai = a->sig[i];
405 unsigned long bi = b->sig[i];
416 /* Return true if A is nonzero. */
419 cmp_significand_0 (a)
420 const REAL_VALUE_TYPE *a;
424 for (i = SIGSZ - 1; i >= 0; --i)
431 /* Set bit N of the significand of R. */
434 set_significand_bit (r, n)
438 r->sig[n / HOST_BITS_PER_LONG]
439 |= (unsigned long)1 << (n % HOST_BITS_PER_LONG);
442 /* Clear bit N of the significand of R. */
445 clear_significand_bit (r, n)
449 r->sig[n / HOST_BITS_PER_LONG]
450 &= ~((unsigned long)1 << (n % HOST_BITS_PER_LONG));
453 /* Test bit N of the significand of R. */
456 test_significand_bit (r, n)
460 /* ??? Compiler bug here if we return this expression directly.
461 The conversion to bool strips the "&1" and we wind up testing
462 e.g. 2 != 0 -> true. Seen in gcc version 3.2 20020520. */
463 int t = (r->sig[n / HOST_BITS_PER_LONG] >> (n % HOST_BITS_PER_LONG)) & 1;
467 /* Clear bits 0..N-1 of the significand of R. */
470 clear_significand_below (r, n)
474 int i, w = n / HOST_BITS_PER_LONG;
476 for (i = 0; i < w; ++i)
479 r->sig[w] &= ~(((unsigned long)1 << (n % HOST_BITS_PER_LONG)) - 1);
482 /* Divide the significands of A and B, placing the result in R. Return
483 true if the division was inexact. */
486 div_significands (r, a, b)
488 const REAL_VALUE_TYPE *a, *b;
491 int i, bit = SIGNIFICAND_BITS - 1;
492 unsigned long msb, inexact;
495 memset (r->sig, 0, sizeof (r->sig));
501 msb = u.sig[SIGSZ-1] & SIG_MSB;
502 lshift_significand_1 (&u, &u);
504 if (msb || cmp_significands (&u, b) >= 0)
506 sub_significands (&u, &u, b, 0);
507 set_significand_bit (r, bit);
512 for (i = 0, inexact = 0; i < SIGSZ; i++)
518 /* Adjust the exponent and significand of R such that the most
519 significant bit is set. We underflow to zero and overflow to
520 infinity here, without denormals. (The intermediate representation
521 exponent is large enough to handle target denormals normalized.) */
530 /* Find the first word that is nonzero. */
531 for (i = SIGSZ - 1; i >= 0; i--)
533 shift += HOST_BITS_PER_LONG;
537 /* Zero significand flushes to zero. */
545 /* Find the first bit that is nonzero. */
547 if (r->sig[i] & ((unsigned long)1 << (HOST_BITS_PER_LONG - 1 - j)))
553 exp = r->exp - shift;
555 get_inf (r, r->sign);
556 else if (exp < -MAX_EXP)
557 get_zero (r, r->sign);
561 lshift_significand (r, r, shift);
566 /* Calculate R = A + (SUBTRACT_P ? -B : B). Return true if the
567 result may be inexact due to a loss of precision. */
570 do_add (r, a, b, subtract_p)
572 const REAL_VALUE_TYPE *a, *b;
577 bool inexact = false;
579 /* Determine if we need to add or subtract. */
581 subtract_p = (sign ^ b->sign) ^ subtract_p;
583 switch (CLASS2 (a->class, b->class))
585 case CLASS2 (rvc_zero, rvc_zero):
586 /* -0 + -0 = -0, -0 - +0 = -0; all other cases yield +0. */
587 get_zero (r, sign & !subtract_p);
590 case CLASS2 (rvc_zero, rvc_normal):
591 case CLASS2 (rvc_zero, rvc_inf):
592 case CLASS2 (rvc_zero, rvc_nan):
594 case CLASS2 (rvc_normal, rvc_nan):
595 case CLASS2 (rvc_inf, rvc_nan):
596 case CLASS2 (rvc_nan, rvc_nan):
597 /* ANY + NaN = NaN. */
598 case CLASS2 (rvc_normal, rvc_inf):
601 r->sign = sign ^ subtract_p;
604 case CLASS2 (rvc_normal, rvc_zero):
605 case CLASS2 (rvc_inf, rvc_zero):
606 case CLASS2 (rvc_nan, rvc_zero):
608 case CLASS2 (rvc_nan, rvc_normal):
609 case CLASS2 (rvc_nan, rvc_inf):
610 /* NaN + ANY = NaN. */
611 case CLASS2 (rvc_inf, rvc_normal):
616 case CLASS2 (rvc_inf, rvc_inf):
618 /* Inf - Inf = NaN. */
619 get_canonical_qnan (r, 0);
621 /* Inf + Inf = Inf. */
625 case CLASS2 (rvc_normal, rvc_normal):
632 /* Swap the arguments such that A has the larger exponent. */
633 dexp = a->exp - b->exp;
636 const REAL_VALUE_TYPE *t;
643 /* If the exponents are not identical, we need to shift the
644 significand of B down. */
647 /* If the exponents are too far apart, the significands
648 do not overlap, which makes the subtraction a noop. */
649 if (dexp >= SIGNIFICAND_BITS)
656 inexact |= sticky_rshift_significand (&t, b, dexp);
662 if (sub_significands (r, a, b, inexact))
664 /* We got a borrow out of the subtraction. That means that
665 A and B had the same exponent, and B had the larger
666 significand. We need to swap the sign and negate the
669 neg_significand (r, r);
674 if (add_significands (r, a, b))
676 /* We got carry out of the addition. This means we need to
677 shift the significand back down one bit and increase the
679 inexact |= sticky_rshift_significand (r, r, 1);
680 r->sig[SIGSZ-1] |= SIG_MSB;
689 r->class = rvc_normal;
693 /* Re-normalize the result. */
696 /* Special case: if the subtraction results in zero, the result
698 if (r->class == rvc_zero)
701 r->sig[0] |= inexact;
706 /* Calculate R = A * B. Return true if the result may be inexact. */
709 do_multiply (r, a, b)
711 const REAL_VALUE_TYPE *a, *b;
713 REAL_VALUE_TYPE u, t, *rr;
714 unsigned int i, j, k;
715 int sign = a->sign ^ b->sign;
716 bool inexact = false;
718 switch (CLASS2 (a->class, b->class))
720 case CLASS2 (rvc_zero, rvc_zero):
721 case CLASS2 (rvc_zero, rvc_normal):
722 case CLASS2 (rvc_normal, rvc_zero):
723 /* +-0 * ANY = 0 with appropriate sign. */
727 case CLASS2 (rvc_zero, rvc_nan):
728 case CLASS2 (rvc_normal, rvc_nan):
729 case CLASS2 (rvc_inf, rvc_nan):
730 case CLASS2 (rvc_nan, rvc_nan):
731 /* ANY * NaN = NaN. */
736 case CLASS2 (rvc_nan, rvc_zero):
737 case CLASS2 (rvc_nan, rvc_normal):
738 case CLASS2 (rvc_nan, rvc_inf):
739 /* NaN * ANY = NaN. */
744 case CLASS2 (rvc_zero, rvc_inf):
745 case CLASS2 (rvc_inf, rvc_zero):
747 get_canonical_qnan (r, sign);
750 case CLASS2 (rvc_inf, rvc_inf):
751 case CLASS2 (rvc_normal, rvc_inf):
752 case CLASS2 (rvc_inf, rvc_normal):
753 /* Inf * Inf = Inf, R * Inf = Inf */
757 case CLASS2 (rvc_normal, rvc_normal):
764 if (r == a || r == b)
770 /* Collect all the partial products. Since we don't have sure access
771 to a widening multiply, we split each long into two half-words.
773 Consider the long-hand form of a four half-word multiplication:
783 We construct partial products of the widened half-word products
784 that are known to not overlap, e.g. DF+DH. Each such partial
785 product is given its proper exponent, which allows us to sum them
786 and obtain the finished product. */
788 for (i = 0; i < SIGSZ * 2; ++i)
790 unsigned long ai = a->sig[i / 2];
792 ai >>= HOST_BITS_PER_LONG / 2;
794 ai &= ((unsigned long)1 << (HOST_BITS_PER_LONG / 2)) - 1;
799 for (j = 0; j < 2; ++j)
801 int exp = (a->exp - (2*SIGSZ-1-i)*(HOST_BITS_PER_LONG/2)
802 + (b->exp - (1-j)*(HOST_BITS_PER_LONG/2)));
811 /* Would underflow to zero, which we shouldn't bother adding. */
816 u.class = rvc_normal;
820 for (k = j; k < SIGSZ * 2; k += 2)
822 unsigned long bi = b->sig[k / 2];
824 bi >>= HOST_BITS_PER_LONG / 2;
826 bi &= ((unsigned long)1 << (HOST_BITS_PER_LONG / 2)) - 1;
828 u.sig[k / 2] = ai * bi;
832 inexact |= do_add (rr, rr, &u, 0);
843 /* Calculate R = A / B. Return true if the result may be inexact. */
848 const REAL_VALUE_TYPE *a, *b;
850 int exp, sign = a->sign ^ b->sign;
851 REAL_VALUE_TYPE t, *rr;
854 switch (CLASS2 (a->class, b->class))
856 case CLASS2 (rvc_zero, rvc_zero):
858 case CLASS2 (rvc_inf, rvc_inf):
859 /* Inf / Inf = NaN. */
860 get_canonical_qnan (r, sign);
863 case CLASS2 (rvc_zero, rvc_normal):
864 case CLASS2 (rvc_zero, rvc_inf):
866 case CLASS2 (rvc_normal, rvc_inf):
871 case CLASS2 (rvc_normal, rvc_zero):
873 case CLASS2 (rvc_inf, rvc_zero):
878 case CLASS2 (rvc_zero, rvc_nan):
879 case CLASS2 (rvc_normal, rvc_nan):
880 case CLASS2 (rvc_inf, rvc_nan):
881 case CLASS2 (rvc_nan, rvc_nan):
882 /* ANY / NaN = NaN. */
887 case CLASS2 (rvc_nan, rvc_zero):
888 case CLASS2 (rvc_nan, rvc_normal):
889 case CLASS2 (rvc_nan, rvc_inf):
890 /* NaN / ANY = NaN. */
895 case CLASS2 (rvc_inf, rvc_normal):
900 case CLASS2 (rvc_normal, rvc_normal):
907 if (r == a || r == b)
912 rr->class = rvc_normal;
915 exp = a->exp - b->exp + 1;
928 inexact = div_significands (rr, a, b);
930 /* Re-normalize the result. */
932 rr->sig[0] |= inexact;
940 /* Return a tri-state comparison of A vs B. Return NAN_RESULT if
941 one of the two operands is a NaN. */
944 do_compare (a, b, nan_result)
945 const REAL_VALUE_TYPE *a, *b;
950 switch (CLASS2 (a->class, b->class))
952 case CLASS2 (rvc_zero, rvc_zero):
953 /* Sign of zero doesn't matter for compares. */
956 case CLASS2 (rvc_inf, rvc_zero):
957 case CLASS2 (rvc_inf, rvc_normal):
958 case CLASS2 (rvc_normal, rvc_zero):
959 return (a->sign ? -1 : 1);
961 case CLASS2 (rvc_inf, rvc_inf):
962 return -a->sign - -b->sign;
964 case CLASS2 (rvc_zero, rvc_normal):
965 case CLASS2 (rvc_zero, rvc_inf):
966 case CLASS2 (rvc_normal, rvc_inf):
967 return (b->sign ? 1 : -1);
969 case CLASS2 (rvc_zero, rvc_nan):
970 case CLASS2 (rvc_normal, rvc_nan):
971 case CLASS2 (rvc_inf, rvc_nan):
972 case CLASS2 (rvc_nan, rvc_nan):
973 case CLASS2 (rvc_nan, rvc_zero):
974 case CLASS2 (rvc_nan, rvc_normal):
975 case CLASS2 (rvc_nan, rvc_inf):
978 case CLASS2 (rvc_normal, rvc_normal):
985 if (a->sign != b->sign)
986 return -a->sign - -b->sign;
990 else if (a->exp < b->exp)
993 ret = cmp_significands (a, b);
995 return (a->sign ? -ret : ret);
998 /* Return A truncated to an integral value toward zero. */
1003 const REAL_VALUE_TYPE *a;
1016 get_zero (r, r->sign);
1017 else if (r->exp < SIGNIFICAND_BITS)
1018 clear_significand_below (r, SIGNIFICAND_BITS - r->exp);
1026 /* Perform the binary or unary operation described by CODE.
1027 For a unary operation, leave OP1 NULL. */
1030 real_arithmetic (r, icode, op0, op1)
1033 const REAL_VALUE_TYPE *op0, *op1;
1035 enum tree_code code = icode;
1040 do_add (r, op0, op1, 0);
1044 do_add (r, op0, op1, 1);
1048 do_multiply (r, op0, op1);
1052 do_divide (r, op0, op1);
1056 if (op1->class == rvc_nan)
1058 else if (do_compare (op0, op1, -1) < 0)
1065 if (op1->class == rvc_nan)
1067 else if (do_compare (op0, op1, 1) < 0)
1083 case FIX_TRUNC_EXPR:
1084 do_fix_trunc (r, op0);
1092 /* Legacy. Similar, but return the result directly. */
1095 real_arithmetic2 (icode, op0, op1)
1097 const REAL_VALUE_TYPE *op0, *op1;
1100 real_arithmetic (&r, icode, op0, op1);
1105 real_compare (icode, op0, op1)
1107 const REAL_VALUE_TYPE *op0, *op1;
1109 enum tree_code code = icode;
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, -1) >= 0;
1122 return do_compare (op0, op1, -1) == 0;
1124 return do_compare (op0, op1, -1) != 0;
1125 case UNORDERED_EXPR:
1126 return op0->class == rvc_nan || op1->class == rvc_nan;
1128 return op0->class != rvc_nan && op1->class != rvc_nan;
1130 return do_compare (op0, op1, -1) < 0;
1132 return do_compare (op0, op1, -1) <= 0;
1134 return do_compare (op0, op1, 1) > 0;
1136 return do_compare (op0, op1, 1) >= 0;
1138 return do_compare (op0, op1, 0) == 0;
1145 /* Return floor log2(R). */
1149 const REAL_VALUE_TYPE *r;
1157 return (unsigned int)-1 >> 1;
1165 /* R = OP0 * 2**EXP. */
1168 real_ldexp (r, op0, exp)
1170 const REAL_VALUE_TYPE *op0;
1184 get_inf (r, r->sign);
1185 else if (exp < -MAX_EXP)
1186 get_zero (r, r->sign);
1196 /* Determine whether a floating-point value X is infinite. */
1200 const REAL_VALUE_TYPE *r;
1202 return (r->class == rvc_inf);
1205 /* Determine whether a floating-point value X is a NaN. */
1209 const REAL_VALUE_TYPE *r;
1211 return (r->class == rvc_nan);
1214 /* Determine whether a floating-point value X is negative. */
1218 const REAL_VALUE_TYPE *r;
1223 /* Determine whether a floating-point value X is minus zero. */
1227 const REAL_VALUE_TYPE *r;
1229 return r->sign && r->class == rvc_zero;
1232 /* Compare two floating-point objects for bitwise identity. */
1235 real_identical (a, b)
1236 const REAL_VALUE_TYPE *a, *b;
1240 if (a->class != b->class)
1242 if (a->sign != b->sign)
1252 if (a->exp != b->exp)
1257 if (a->signalling != b->signalling)
1259 /* The significand is ignored for canonical NaNs. */
1260 if (a->canonical || b->canonical)
1261 return a->canonical == b->canonical;
1268 for (i = 0; i < SIGSZ; ++i)
1269 if (a->sig[i] != b->sig[i])
1275 /* Try to change R into its exact multiplicative inverse in machine
1276 mode MODE. Return true if successful. */
1279 exact_real_inverse (mode, r)
1280 enum machine_mode mode;
1283 const REAL_VALUE_TYPE *one = real_digit (1);
1287 if (r->class != rvc_normal)
1290 /* Check for a power of two: all significand bits zero except the MSB. */
1291 for (i = 0; i < SIGSZ-1; ++i)
1294 if (r->sig[SIGSZ-1] != SIG_MSB)
1297 /* Find the inverse and truncate to the required mode. */
1298 do_divide (&u, one, r);
1299 real_convert (&u, mode, &u);
1301 /* The rounding may have overflowed. */
1302 if (u.class != rvc_normal)
1304 for (i = 0; i < SIGSZ-1; ++i)
1307 if (u.sig[SIGSZ-1] != SIG_MSB)
1314 /* Render R as an integer. */
1318 const REAL_VALUE_TYPE *r;
1320 unsigned HOST_WIDE_INT i;
1331 i = (unsigned HOST_WIDE_INT) 1 << (HOST_BITS_PER_WIDE_INT - 1);
1339 /* Only force overflow for unsigned overflow. Signed overflow is
1340 undefined, so it doesn't matter what we return, and some callers
1341 expect to be able to use this routine for both signed and
1342 unsigned conversions. */
1343 if (r->exp > HOST_BITS_PER_WIDE_INT)
1346 if (HOST_BITS_PER_WIDE_INT == HOST_BITS_PER_LONG)
1347 i = r->sig[SIGSZ-1];
1348 else if (HOST_BITS_PER_WIDE_INT == 2*HOST_BITS_PER_LONG)
1350 i = r->sig[SIGSZ-1];
1351 i = i << (HOST_BITS_PER_LONG - 1) << 1;
1352 i |= r->sig[SIGSZ-2];
1357 i >>= HOST_BITS_PER_WIDE_INT - r->exp;
1368 /* Likewise, but to an integer pair, HI+LOW. */
1371 real_to_integer2 (plow, phigh, r)
1372 HOST_WIDE_INT *plow, *phigh;
1373 const REAL_VALUE_TYPE *r;
1376 HOST_WIDE_INT low, high;
1389 high = (unsigned HOST_WIDE_INT) 1 << (HOST_BITS_PER_WIDE_INT - 1);
1403 /* Only force overflow for unsigned overflow. Signed overflow is
1404 undefined, so it doesn't matter what we return, and some callers
1405 expect to be able to use this routine for both signed and
1406 unsigned conversions. */
1407 if (exp > 2*HOST_BITS_PER_WIDE_INT)
1410 rshift_significand (&t, r, 2*HOST_BITS_PER_WIDE_INT - exp);
1411 if (HOST_BITS_PER_WIDE_INT == HOST_BITS_PER_LONG)
1413 high = t.sig[SIGSZ-1];
1414 low = t.sig[SIGSZ-2];
1416 else if (HOST_BITS_PER_WIDE_INT == 2*HOST_BITS_PER_LONG)
1418 high = t.sig[SIGSZ-1];
1419 high = high << (HOST_BITS_PER_LONG - 1) << 1;
1420 high |= t.sig[SIGSZ-2];
1422 low = t.sig[SIGSZ-3];
1423 low = low << (HOST_BITS_PER_LONG - 1) << 1;
1424 low |= t.sig[SIGSZ-4];
1434 low = -low, high = ~high;
1446 /* A subroutine of real_to_decimal. Compute the quotient and remainder
1447 of NUM / DEN. Return the quotient and place the remainder in NUM.
1448 It is expected that NUM / DEN are close enough that the quotient is
1451 static unsigned long
1452 rtd_divmod (num, den)
1453 REAL_VALUE_TYPE *num, *den;
1455 unsigned long q, msb;
1456 int expn = num->exp, expd = den->exp;
1465 msb = num->sig[SIGSZ-1] & SIG_MSB;
1467 lshift_significand_1 (num, num);
1469 if (msb || cmp_significands (num, den) >= 0)
1471 sub_significands (num, num, den, 0);
1475 while (--expn >= expd);
1483 /* Render R as a decimal floating point constant. Emit DIGITS significant
1484 digits in the result, bounded by BUF_SIZE. If DIGITS is 0, choose the
1485 maximum for the representation. If CROP_TRAILING_ZEROS, strip trailing
1488 #define M_LOG10_2 0.30102999566398119521
1491 real_to_decimal (str, r_orig, buf_size, digits, crop_trailing_zeros)
1493 const REAL_VALUE_TYPE *r_orig;
1494 size_t buf_size, digits;
1495 int crop_trailing_zeros;
1497 const REAL_VALUE_TYPE *one, *ten;
1498 REAL_VALUE_TYPE r, pten, u, v;
1499 int dec_exp, cmp_one, digit;
1501 char *p, *first, *last;
1508 strcpy (str, (r.sign ? "-0.0" : "0.0"));
1513 strcpy (str, (r.sign ? "-Inf" : "+Inf"));
1516 /* ??? Print the significand as well, if not canonical? */
1517 strcpy (str, (r.sign ? "-NaN" : "+NaN"));
1523 /* Bound the number of digits printed by the size of the representation. */
1524 max_digits = SIGNIFICAND_BITS * M_LOG10_2;
1525 if (digits == 0 || digits > max_digits)
1526 digits = max_digits;
1528 /* Estimate the decimal exponent, and compute the length of the string it
1529 will print as. Be conservative and add one to account for possible
1530 overflow or rounding error. */
1531 dec_exp = r.exp * M_LOG10_2;
1532 for (max_digits = 1; dec_exp ; max_digits++)
1535 /* Bound the number of digits printed by the size of the output buffer. */
1536 max_digits = buf_size - 1 - 1 - 2 - max_digits - 1;
1537 if (max_digits > buf_size)
1539 if (digits > max_digits)
1540 digits = max_digits;
1542 one = real_digit (1);
1543 ten = ten_to_ptwo (0);
1551 cmp_one = do_compare (&r, one, 0);
1556 /* Number is greater than one. Convert significand to an integer
1557 and strip trailing decimal zeros. */
1560 u.exp = SIGNIFICAND_BITS - 1;
1562 /* Largest M, such that 10**2**M fits within SIGNIFICAND_BITS. */
1563 m = floor_log2 (max_digits);
1565 /* Iterate over the bits of the possible powers of 10 that might
1566 be present in U and eliminate them. That is, if we find that
1567 10**2**M divides U evenly, keep the division and increase
1573 do_divide (&t, &u, ten_to_ptwo (m));
1574 do_fix_trunc (&v, &t);
1575 if (cmp_significands (&v, &t) == 0)
1583 /* Revert the scaling to integer that we performed earlier. */
1584 u.exp += r.exp - (SIGNIFICAND_BITS - 1);
1587 /* Find power of 10. Do this by dividing out 10**2**M when
1588 this is larger than the current remainder. Fill PTEN with
1589 the power of 10 that we compute. */
1592 m = floor_log2 ((int)(r.exp * M_LOG10_2)) + 1;
1595 const REAL_VALUE_TYPE *ptentwo = ten_to_ptwo (m);
1596 if (do_compare (&u, ptentwo, 0) >= 0)
1598 do_divide (&u, &u, ptentwo);
1599 do_multiply (&pten, &pten, ptentwo);
1606 /* We managed to divide off enough tens in the above reduction
1607 loop that we've now got a negative exponent. Fall into the
1608 less-than-one code to compute the proper value for PTEN. */
1615 /* Number is less than one. Pad significand with leading
1621 /* Stop if we'd shift bits off the bottom. */
1625 do_multiply (&u, &v, ten);
1627 /* Stop if we're now >= 1. */
1636 /* Find power of 10. Do this by multiplying in P=10**2**M when
1637 the current remainder is smaller than 1/P. Fill PTEN with the
1638 power of 10 that we compute. */
1639 m = floor_log2 ((int)(-r.exp * M_LOG10_2)) + 1;
1642 const REAL_VALUE_TYPE *ptentwo = ten_to_ptwo (m);
1643 const REAL_VALUE_TYPE *ptenmtwo = ten_to_mptwo (m);
1645 if (do_compare (&v, ptenmtwo, 0) <= 0)
1647 do_multiply (&v, &v, ptentwo);
1648 do_multiply (&pten, &pten, ptentwo);
1654 /* Invert the positive power of 10 that we've collected so far. */
1655 do_divide (&pten, one, &pten);
1663 /* At this point, PTEN should contain the nearest power of 10 smaller
1664 than R, such that this division produces the first digit.
1666 Using a divide-step primitive that returns the complete integral
1667 remainder avoids the rounding error that would be produced if
1668 we were to use do_divide here and then simply multiply by 10 for
1669 each subsequent digit. */
1671 digit = rtd_divmod (&r, &pten);
1673 /* Be prepared for error in that division via underflow ... */
1674 if (digit == 0 && cmp_significand_0 (&r))
1676 /* Multiply by 10 and try again. */
1677 do_multiply (&r, &r, ten);
1678 digit = rtd_divmod (&r, &pten);
1684 /* ... or overflow. */
1692 else if (digit > 10)
1697 /* Generate subsequent digits. */
1698 while (--digits > 0)
1700 do_multiply (&r, &r, ten);
1701 digit = rtd_divmod (&r, &pten);
1706 /* Generate one more digit with which to do rounding. */
1707 do_multiply (&r, &r, ten);
1708 digit = rtd_divmod (&r, &pten);
1710 /* Round the result. */
1713 /* Round to nearest. If R is nonzero there are additional
1714 nonzero digits to be extracted. */
1715 if (cmp_significand_0 (&r))
1717 /* Round to even. */
1718 else if ((p[-1] - '0') & 1)
1735 /* Carry out of the first digit. This means we had all 9's and
1736 now have all 0's. "Prepend" a 1 by overwriting the first 0. */
1744 /* Insert the decimal point. */
1745 first[0] = first[1];
1748 /* If requested, drop trailing zeros. Never crop past "1.0". */
1749 if (crop_trailing_zeros)
1750 while (last > first + 3 && last[-1] == '0')
1753 /* Append the exponent. */
1754 sprintf (last, "e%+d", dec_exp);
1757 /* Render R as a hexadecimal floating point constant. Emit DIGITS
1758 significant digits in the result, bounded by BUF_SIZE. If DIGITS is 0,
1759 choose the maximum for the representation. If CROP_TRAILING_ZEROS,
1760 strip trailing zeros. */
1763 real_to_hexadecimal (str, r, buf_size, digits, crop_trailing_zeros)
1765 const REAL_VALUE_TYPE *r;
1766 size_t buf_size, digits;
1767 int crop_trailing_zeros;
1769 int i, j, exp = r->exp;
1782 strcpy (str, (r->sign ? "-Inf" : "+Inf"));
1785 /* ??? Print the significand as well, if not canonical? */
1786 strcpy (str, (r->sign ? "-NaN" : "+NaN"));
1793 digits = SIGNIFICAND_BITS / 4;
1795 /* Bound the number of digits printed by the size of the output buffer. */
1797 sprintf (exp_buf, "p%+d", exp);
1798 max_digits = buf_size - strlen (exp_buf) - r->sign - 4 - 1;
1799 if (max_digits > buf_size)
1801 if (digits > max_digits)
1802 digits = max_digits;
1813 for (i = SIGSZ - 1; i >= 0; --i)
1814 for (j = HOST_BITS_PER_LONG - 4; j >= 0; j -= 4)
1816 *p++ = "0123456789abcdef"[(r->sig[i] >> j) & 15];
1822 if (crop_trailing_zeros)
1823 while (p > first + 1 && p[-1] == '0')
1826 sprintf (p, "p%+d", exp);
1829 /* Initialize R from a decimal or hexadecimal string. The string is
1830 assumed to have been syntax checked already. */
1833 real_from_string (r, str)
1847 else if (*str == '+')
1850 if (str[0] == '0' && str[1] == 'x')
1852 /* Hexadecimal floating point. */
1853 int pos = SIGNIFICAND_BITS - 4, d;
1861 d = hex_value (*str);
1866 r->sig[pos / HOST_BITS_PER_LONG]
1867 |= (unsigned long) d << (pos % HOST_BITS_PER_LONG);
1876 if (pos == SIGNIFICAND_BITS - 4)
1883 d = hex_value (*str);
1888 r->sig[pos / HOST_BITS_PER_LONG]
1889 |= (unsigned long) d << (pos % HOST_BITS_PER_LONG);
1895 if (*str == 'p' || *str == 'P')
1897 bool exp_neg = false;
1905 else if (*str == '+')
1909 while (ISDIGIT (*str))
1915 /* Overflowed the exponent. */
1929 r->class = rvc_normal;
1936 /* Decimal floating point. */
1937 const REAL_VALUE_TYPE *ten = ten_to_ptwo (0);
1942 while (ISDIGIT (*str))
1945 do_multiply (r, r, ten);
1947 do_add (r, r, real_digit (d), 0);
1952 if (r->class == rvc_zero)
1957 while (ISDIGIT (*str))
1960 do_multiply (r, r, ten);
1962 do_add (r, r, real_digit (d), 0);
1967 if (*str == 'e' || *str == 'E')
1969 bool exp_neg = false;
1977 else if (*str == '+')
1981 while (ISDIGIT (*str))
1987 /* Overflowed the exponent. */
2001 times_pten (r, exp);
2016 /* Legacy. Similar, but return the result directly. */
2019 real_from_string2 (s, mode)
2021 enum machine_mode mode;
2025 real_from_string (&r, s);
2026 if (mode != VOIDmode)
2027 real_convert (&r, mode, &r);
2032 /* Initialize R from the integer pair HIGH+LOW. */
2035 real_from_integer (r, mode, low, high, unsigned_p)
2037 enum machine_mode mode;
2038 unsigned HOST_WIDE_INT low;
2042 if (low == 0 && high == 0)
2046 r->class = rvc_normal;
2047 r->sign = high < 0 && !unsigned_p;
2048 r->exp = 2 * HOST_BITS_PER_WIDE_INT;
2059 if (HOST_BITS_PER_LONG == HOST_BITS_PER_WIDE_INT)
2061 r->sig[SIGSZ-1] = high;
2062 r->sig[SIGSZ-2] = low;
2063 memset (r->sig, 0, sizeof(long)*(SIGSZ-2));
2065 else if (HOST_BITS_PER_LONG*2 == HOST_BITS_PER_WIDE_INT)
2067 r->sig[SIGSZ-1] = high >> (HOST_BITS_PER_LONG - 1) >> 1;
2068 r->sig[SIGSZ-2] = high;
2069 r->sig[SIGSZ-3] = low >> (HOST_BITS_PER_LONG - 1) >> 1;
2070 r->sig[SIGSZ-4] = low;
2072 memset (r->sig, 0, sizeof(long)*(SIGSZ-4));
2080 if (mode != VOIDmode)
2081 real_convert (r, mode, r);
2084 /* Returns 10**2**N. */
2086 static const REAL_VALUE_TYPE *
2090 static REAL_VALUE_TYPE tens[EXP_BITS];
2092 if (n < 0 || n >= EXP_BITS)
2095 if (tens[n].class == rvc_zero)
2097 if (n < (HOST_BITS_PER_WIDE_INT == 64 ? 5 : 4))
2099 HOST_WIDE_INT t = 10;
2102 for (i = 0; i < n; ++i)
2105 real_from_integer (&tens[n], VOIDmode, t, 0, 1);
2109 const REAL_VALUE_TYPE *t = ten_to_ptwo (n - 1);
2110 do_multiply (&tens[n], t, t);
2117 /* Returns 10**(-2**N). */
2119 static const REAL_VALUE_TYPE *
2123 static REAL_VALUE_TYPE tens[EXP_BITS];
2125 if (n < 0 || n >= EXP_BITS)
2128 if (tens[n].class == rvc_zero)
2129 do_divide (&tens[n], real_digit (1), ten_to_ptwo (n));
2136 static const REAL_VALUE_TYPE *
2140 static REAL_VALUE_TYPE num[10];
2145 if (n > 0 && num[n].class == rvc_zero)
2146 real_from_integer (&num[n], VOIDmode, n, 0, 1);
2151 /* Multiply R by 10**EXP. */
2158 REAL_VALUE_TYPE pten, *rr;
2159 bool negative = (exp < 0);
2165 pten = *real_digit (1);
2171 for (i = 0; exp > 0; ++i, exp >>= 1)
2173 do_multiply (rr, rr, ten_to_ptwo (i));
2176 do_divide (r, r, &pten);
2179 /* Fills R with +Inf. */
2188 /* Fills R with a NaN whose significand is described by STR. If QUIET,
2189 we force a QNaN, else we force an SNaN. The string, if not empty,
2190 is parsed as a number and placed in the significand. Return true
2191 if the string was successfully parsed. */
2194 real_nan (r, str, quiet, mode)
2198 enum machine_mode mode;
2200 const struct real_format *fmt;
2202 fmt = real_format_for_mode[mode - QFmode];
2209 get_canonical_qnan (r, 0);
2211 get_canonical_snan (r, 0);
2218 memset (r, 0, sizeof (*r));
2221 /* Parse akin to strtol into the significand of R. */
2223 while (ISSPACE (*str))
2227 else if (*str == '+')
2237 while ((d = hex_value (*str)) < base)
2244 lshift_significand (r, r, 3);
2247 lshift_significand (r, r, 4);
2250 lshift_significand_1 (&u, r);
2251 lshift_significand (r, r, 3);
2252 add_significands (r, r, &u);
2260 add_significands (r, r, &u);
2265 /* Must have consumed the entire string for success. */
2269 /* Shift the significand into place such that the bits
2270 are in the most significant bits for the format. */
2271 lshift_significand (r, r, SIGNIFICAND_BITS - fmt->pnan);
2273 /* Our MSB is always unset for NaNs. */
2274 r->sig[SIGSZ-1] &= ~SIG_MSB;
2276 /* Force quiet or signalling NaN. */
2277 r->signalling = !quiet;
2283 /* Fills R with the largest finite value representable in mode MODE.
2284 If SIGN is nonzero, R is set to the most negative finite value. */
2287 real_maxval (r, sign, mode)
2290 enum machine_mode mode;
2292 const struct real_format *fmt;
2295 fmt = real_format_for_mode[mode - QFmode];
2299 r->class = rvc_normal;
2303 r->exp = fmt->emax * fmt->log2_b;
2305 np2 = SIGNIFICAND_BITS - fmt->p * fmt->log2_b;
2306 memset (r->sig, -1, SIGSZ * sizeof (unsigned long));
2307 clear_significand_below (r, np2);
2310 /* Fills R with 2**N. */
2317 memset (r, 0, sizeof (*r));
2322 else if (n < -MAX_EXP)
2326 r->class = rvc_normal;
2328 r->sig[SIGSZ-1] = SIG_MSB;
2334 round_for_format (fmt, r)
2335 const struct real_format *fmt;
2339 unsigned long sticky;
2343 p2 = fmt->p * fmt->log2_b;
2344 emin2m1 = (fmt->emin - 1) * fmt->log2_b;
2345 emax2 = fmt->emax * fmt->log2_b;
2347 np2 = SIGNIFICAND_BITS - p2;
2351 get_zero (r, r->sign);
2353 if (!fmt->has_signed_zero)
2358 get_inf (r, r->sign);
2363 clear_significand_below (r, np2);
2373 /* If we're not base2, normalize the exponent to a multiple of
2375 if (fmt->log2_b != 1)
2377 int shift = r->exp & (fmt->log2_b - 1);
2380 shift = fmt->log2_b - shift;
2381 r->sig[0] |= sticky_rshift_significand (r, r, shift);
2386 /* Check the range of the exponent. If we're out of range,
2387 either underflow or overflow. */
2390 else if (r->exp <= emin2m1)
2394 if (!fmt->has_denorm)
2396 /* Don't underflow completely until we've had a chance to round. */
2397 if (r->exp < emin2m1)
2402 diff = emin2m1 - r->exp + 1;
2406 /* De-normalize the significand. */
2407 r->sig[0] |= sticky_rshift_significand (r, r, diff);
2412 /* There are P2 true significand bits, followed by one guard bit,
2413 followed by one sticky bit, followed by stuff. Fold nonzero
2414 stuff into the sticky bit. */
2417 for (i = 0, w = (np2 - 1) / HOST_BITS_PER_LONG; i < w; ++i)
2418 sticky |= r->sig[i];
2420 r->sig[w] & (((unsigned long)1 << ((np2 - 1) % HOST_BITS_PER_LONG)) - 1);
2422 guard = test_significand_bit (r, np2 - 1);
2423 lsb = test_significand_bit (r, np2);
2425 /* Round to even. */
2426 if (guard && (sticky || lsb))
2430 set_significand_bit (&u, np2);
2432 if (add_significands (r, r, &u))
2434 /* Overflow. Means the significand had been all ones, and
2435 is now all zeros. Need to increase the exponent, and
2436 possibly re-normalize it. */
2437 if (++r->exp > emax2)
2439 r->sig[SIGSZ-1] = SIG_MSB;
2441 if (fmt->log2_b != 1)
2443 int shift = r->exp & (fmt->log2_b - 1);
2446 shift = fmt->log2_b - shift;
2447 rshift_significand (r, r, shift);
2456 /* Catch underflow that we deferred until after rounding. */
2457 if (r->exp <= emin2m1)
2460 /* Clear out trailing garbage. */
2461 clear_significand_below (r, np2);
2464 /* Extend or truncate to a new mode. */
2467 real_convert (r, mode, a)
2469 enum machine_mode mode;
2470 const REAL_VALUE_TYPE *a;
2472 const struct real_format *fmt;
2474 fmt = real_format_for_mode[mode - QFmode];
2479 round_for_format (fmt, r);
2481 /* round_for_format de-normalizes denormals. Undo just that part. */
2482 if (r->class == rvc_normal)
2486 /* Legacy. Likewise, except return the struct directly. */
2489 real_value_truncate (mode, a)
2490 enum machine_mode mode;
2494 real_convert (&r, mode, &a);
2498 /* Return true if truncating to MODE is exact. */
2501 exact_real_truncate (mode, a)
2502 enum machine_mode mode;
2503 const REAL_VALUE_TYPE *a;
2506 real_convert (&t, mode, a);
2507 return real_identical (&t, a);
2510 /* Write R to the given target format. Place the words of the result
2511 in target word order in BUF. There are always 32 bits in each
2512 long, no matter the size of the host long.
2514 Legacy: return word 0 for implementing REAL_VALUE_TO_TARGET_SINGLE. */
2517 real_to_target_fmt (buf, r_orig, fmt)
2519 const REAL_VALUE_TYPE *r_orig;
2520 const struct real_format *fmt;
2526 round_for_format (fmt, &r);
2530 (*fmt->encode) (fmt, buf, &r);
2535 /* Similar, but look up the format from MODE. */
2538 real_to_target (buf, r, mode)
2540 const REAL_VALUE_TYPE *r;
2541 enum machine_mode mode;
2543 const struct real_format *fmt;
2545 fmt = real_format_for_mode[mode - QFmode];
2549 return real_to_target_fmt (buf, r, fmt);
2552 /* Read R from the given target format. Read the words of the result
2553 in target word order in BUF. There are always 32 bits in each
2554 long, no matter the size of the host long. */
2557 real_from_target_fmt (r, buf, fmt)
2560 const struct real_format *fmt;
2562 (*fmt->decode) (fmt, r, buf);
2565 /* Similar, but look up the format from MODE. */
2568 real_from_target (r, buf, mode)
2571 enum machine_mode mode;
2573 const struct real_format *fmt;
2575 fmt = real_format_for_mode[mode - QFmode];
2579 (*fmt->decode) (fmt, r, buf);
2582 /* Return the number of bits in the significand for MODE. */
2583 /* ??? Legacy. Should get access to real_format directly. */
2586 significand_size (mode)
2587 enum machine_mode mode;
2589 const struct real_format *fmt;
2591 fmt = real_format_for_mode[mode - QFmode];
2595 return fmt->p * fmt->log2_b;
2598 /* Return a hash value for the given real value. */
2599 /* ??? The "unsigned int" return value is intended to be hashval_t,
2600 but I didn't want to pull hashtab.h into real.h. */
2604 const REAL_VALUE_TYPE *r;
2609 h = r->class | (r->sign << 2);
2622 h ^= (unsigned int)-1;
2631 if (sizeof(unsigned long) > sizeof(unsigned int))
2632 for (i = 0; i < SIGSZ; ++i)
2634 unsigned long s = r->sig[i];
2635 h ^= s ^ (s >> (HOST_BITS_PER_LONG / 2));
2638 for (i = 0; i < SIGSZ; ++i)
2644 /* IEEE single-precision format. */
2646 static void encode_ieee_single PARAMS ((const struct real_format *fmt,
2647 long *, const REAL_VALUE_TYPE *));
2648 static void decode_ieee_single PARAMS ((const struct real_format *,
2649 REAL_VALUE_TYPE *, const long *));
2652 encode_ieee_single (fmt, buf, r)
2653 const struct real_format *fmt;
2655 const REAL_VALUE_TYPE *r;
2657 unsigned long image, sig, exp;
2658 bool denormal = (r->sig[SIGSZ-1] & SIG_MSB) == 0;
2660 image = r->sign << 31;
2661 sig = (r->sig[SIGSZ-1] >> (HOST_BITS_PER_LONG - 24)) & 0x7fffff;
2672 image |= 0x7fffffff;
2680 if (r->signalling == fmt->qnan_msb_set)
2684 /* We overload qnan_msb_set here: it's only clear for
2685 mips_ieee_single, which wants all mantissa bits but the
2686 quiet/signalling one set in canonical NaNs (at least
2688 if (r->canonical && !fmt->qnan_msb_set)
2689 sig |= (1 << 22) - 1;
2697 image |= 0x7fffffff;
2701 /* Recall that IEEE numbers are interpreted as 1.F x 2**exp,
2702 whereas the intermediate representation is 0.F x 2**exp.
2703 Which means we're off by one. */
2707 exp = r->exp + 127 - 1;
2720 decode_ieee_single (fmt, r, buf)
2721 const struct real_format *fmt;
2725 unsigned long image = buf[0] & 0xffffffff;
2726 bool sign = (image >> 31) & 1;
2727 int exp = (image >> 23) & 0xff;
2729 memset (r, 0, sizeof (*r));
2730 image <<= HOST_BITS_PER_LONG - 24;
2735 if (image && fmt->has_denorm)
2737 r->class = rvc_normal;
2740 r->sig[SIGSZ-1] = image << 1;
2743 else if (fmt->has_signed_zero)
2746 else if (exp == 255 && (fmt->has_nans || fmt->has_inf))
2752 r->signalling = (((image >> (HOST_BITS_PER_LONG - 2)) & 1)
2753 ^ fmt->qnan_msb_set);
2754 r->sig[SIGSZ-1] = image;
2764 r->class = rvc_normal;
2766 r->exp = exp - 127 + 1;
2767 r->sig[SIGSZ-1] = image | SIG_MSB;
2771 const struct real_format ieee_single_format =
2789 const struct real_format mips_single_format =
2808 /* IEEE double-precision format. */
2810 static void encode_ieee_double PARAMS ((const struct real_format *fmt,
2811 long *, const REAL_VALUE_TYPE *));
2812 static void decode_ieee_double PARAMS ((const struct real_format *,
2813 REAL_VALUE_TYPE *, const long *));
2816 encode_ieee_double (fmt, buf, r)
2817 const struct real_format *fmt;
2819 const REAL_VALUE_TYPE *r;
2821 unsigned long image_lo, image_hi, sig_lo, sig_hi, exp;
2822 bool denormal = (r->sig[SIGSZ-1] & SIG_MSB) == 0;
2824 image_hi = r->sign << 31;
2827 if (HOST_BITS_PER_LONG == 64)
2829 sig_hi = r->sig[SIGSZ-1];
2830 sig_lo = (sig_hi >> (64 - 53)) & 0xffffffff;
2831 sig_hi = (sig_hi >> (64 - 53 + 1) >> 31) & 0xfffff;
2835 sig_hi = r->sig[SIGSZ-1];
2836 sig_lo = r->sig[SIGSZ-2];
2837 sig_lo = (sig_hi << 21) | (sig_lo >> 11);
2838 sig_hi = (sig_hi >> 11) & 0xfffff;
2848 image_hi |= 2047 << 20;
2851 image_hi |= 0x7fffffff;
2852 image_lo = 0xffffffff;
2860 sig_hi = sig_lo = 0;
2861 if (r->signalling == fmt->qnan_msb_set)
2862 sig_hi &= ~(1 << 19);
2865 /* We overload qnan_msb_set here: it's only clear for
2866 mips_ieee_single, which wants all mantissa bits but the
2867 quiet/signalling one set in canonical NaNs (at least
2869 if (r->canonical && !fmt->qnan_msb_set)
2871 sig_hi |= (1 << 19) - 1;
2872 sig_lo = 0xffffffff;
2874 else if (sig_hi == 0 && sig_lo == 0)
2877 image_hi |= 2047 << 20;
2883 image_hi |= 0x7fffffff;
2884 image_lo = 0xffffffff;
2889 /* Recall that IEEE numbers are interpreted as 1.F x 2**exp,
2890 whereas the intermediate representation is 0.F x 2**exp.
2891 Which means we're off by one. */
2895 exp = r->exp + 1023 - 1;
2896 image_hi |= exp << 20;
2905 if (FLOAT_WORDS_BIG_ENDIAN)
2906 buf[0] = image_hi, buf[1] = image_lo;
2908 buf[0] = image_lo, buf[1] = image_hi;
2912 decode_ieee_double (fmt, r, buf)
2913 const struct real_format *fmt;
2917 unsigned long image_hi, image_lo;
2921 if (FLOAT_WORDS_BIG_ENDIAN)
2922 image_hi = buf[0], image_lo = buf[1];
2924 image_lo = buf[0], image_hi = buf[1];
2925 image_lo &= 0xffffffff;
2926 image_hi &= 0xffffffff;
2928 sign = (image_hi >> 31) & 1;
2929 exp = (image_hi >> 20) & 0x7ff;
2931 memset (r, 0, sizeof (*r));
2933 image_hi <<= 32 - 21;
2934 image_hi |= image_lo >> 21;
2935 image_hi &= 0x7fffffff;
2936 image_lo <<= 32 - 21;
2940 if ((image_hi || image_lo) && fmt->has_denorm)
2942 r->class = rvc_normal;
2945 if (HOST_BITS_PER_LONG == 32)
2947 image_hi = (image_hi << 1) | (image_lo >> 31);
2949 r->sig[SIGSZ-1] = image_hi;
2950 r->sig[SIGSZ-2] = image_lo;
2954 image_hi = (image_hi << 31 << 2) | (image_lo << 1);
2955 r->sig[SIGSZ-1] = image_hi;
2959 else if (fmt->has_signed_zero)
2962 else if (exp == 2047 && (fmt->has_nans || fmt->has_inf))
2964 if (image_hi || image_lo)
2968 r->signalling = ((image_hi >> 30) & 1) ^ fmt->qnan_msb_set;
2969 if (HOST_BITS_PER_LONG == 32)
2971 r->sig[SIGSZ-1] = image_hi;
2972 r->sig[SIGSZ-2] = image_lo;
2975 r->sig[SIGSZ-1] = (image_hi << 31 << 1) | image_lo;
2985 r->class = rvc_normal;
2987 r->exp = exp - 1023 + 1;
2988 if (HOST_BITS_PER_LONG == 32)
2990 r->sig[SIGSZ-1] = image_hi | SIG_MSB;
2991 r->sig[SIGSZ-2] = image_lo;
2994 r->sig[SIGSZ-1] = (image_hi << 31 << 1) | image_lo | SIG_MSB;
2998 const struct real_format ieee_double_format =
3016 const struct real_format mips_double_format =
3035 /* IEEE extended double precision format. This comes in three
3036 flavors: Intel's as a 12 byte image, Intel's as a 16 byte image,
3039 static void encode_ieee_extended PARAMS ((const struct real_format *fmt,
3040 long *, const REAL_VALUE_TYPE *));
3041 static void decode_ieee_extended PARAMS ((const struct real_format *,
3042 REAL_VALUE_TYPE *, const long *));
3044 static void encode_ieee_extended_128 PARAMS ((const struct real_format *fmt,
3046 const REAL_VALUE_TYPE *));
3047 static void decode_ieee_extended_128 PARAMS ((const struct real_format *,
3052 encode_ieee_extended (fmt, buf, r)
3053 const struct real_format *fmt;
3055 const REAL_VALUE_TYPE *r;
3057 unsigned long image_hi, sig_hi, sig_lo;
3058 bool denormal = (r->sig[SIGSZ-1] & SIG_MSB) == 0;
3060 image_hi = r->sign << 15;
3061 sig_hi = sig_lo = 0;
3073 /* Intel requires the explicit integer bit to be set, otherwise
3074 it considers the value a "pseudo-infinity". Motorola docs
3075 say it doesn't care. */
3076 sig_hi = 0x80000000;
3081 sig_lo = sig_hi = 0xffffffff;
3089 if (HOST_BITS_PER_LONG == 32)
3091 sig_hi = r->sig[SIGSZ-1];
3092 sig_lo = r->sig[SIGSZ-2];
3096 sig_lo = r->sig[SIGSZ-1];
3097 sig_hi = sig_lo >> 31 >> 1;
3098 sig_lo &= 0xffffffff;
3100 if (r->signalling == fmt->qnan_msb_set)
3101 sig_hi &= ~(1 << 30);
3104 if ((sig_hi & 0x7fffffff) == 0 && sig_lo == 0)
3107 /* Intel requires the explicit integer bit to be set, otherwise
3108 it considers the value a "pseudo-nan". Motorola docs say it
3110 sig_hi |= 0x80000000;
3115 sig_lo = sig_hi = 0xffffffff;
3123 /* Recall that IEEE numbers are interpreted as 1.F x 2**exp,
3124 whereas the intermediate representation is 0.F x 2**exp.
3125 Which means we're off by one.
3127 Except for Motorola, which consider exp=0 and explicit
3128 integer bit set to continue to be normalized. In theory
3129 this discrepancy has been taken care of by the difference
3130 in fmt->emin in round_for_format. */
3142 if (HOST_BITS_PER_LONG == 32)
3144 sig_hi = r->sig[SIGSZ-1];
3145 sig_lo = r->sig[SIGSZ-2];
3149 sig_lo = r->sig[SIGSZ-1];
3150 sig_hi = sig_lo >> 31 >> 1;
3151 sig_lo &= 0xffffffff;
3160 if (FLOAT_WORDS_BIG_ENDIAN)
3161 buf[0] = image_hi << 16, buf[1] = sig_hi, buf[2] = sig_lo;
3163 buf[0] = sig_lo, buf[1] = sig_hi, buf[2] = image_hi;
3167 encode_ieee_extended_128 (fmt, buf, r)
3168 const struct real_format *fmt;
3170 const REAL_VALUE_TYPE *r;
3172 buf[3 * !FLOAT_WORDS_BIG_ENDIAN] = 0;
3173 encode_ieee_extended (fmt, buf+!!FLOAT_WORDS_BIG_ENDIAN, r);
3177 decode_ieee_extended (fmt, r, buf)
3178 const struct real_format *fmt;
3182 unsigned long image_hi, sig_hi, sig_lo;
3186 if (FLOAT_WORDS_BIG_ENDIAN)
3187 image_hi = buf[0] >> 16, sig_hi = buf[1], sig_lo = buf[2];
3189 sig_lo = buf[0], sig_hi = buf[1], image_hi = buf[2];
3190 sig_lo &= 0xffffffff;
3191 sig_hi &= 0xffffffff;
3192 image_hi &= 0xffffffff;
3194 sign = (image_hi >> 15) & 1;
3195 exp = image_hi & 0x7fff;
3197 memset (r, 0, sizeof (*r));
3201 if ((sig_hi || sig_lo) && fmt->has_denorm)
3203 r->class = rvc_normal;
3206 /* When the IEEE format contains a hidden bit, we know that
3207 it's zero at this point, and so shift up the significand
3208 and decrease the exponent to match. In this case, Motorola
3209 defines the explicit integer bit to be valid, so we don't
3210 know whether the msb is set or not. */
3212 if (HOST_BITS_PER_LONG == 32)
3214 r->sig[SIGSZ-1] = sig_hi;
3215 r->sig[SIGSZ-2] = sig_lo;
3218 r->sig[SIGSZ-1] = (sig_hi << 31 << 1) | sig_lo;
3222 else if (fmt->has_signed_zero)
3225 else if (exp == 32767 && (fmt->has_nans || fmt->has_inf))
3227 /* See above re "pseudo-infinities" and "pseudo-nans".
3228 Short summary is that the MSB will likely always be
3229 set, and that we don't care about it. */
3230 sig_hi &= 0x7fffffff;
3232 if (sig_hi || sig_lo)
3236 r->signalling = ((sig_hi >> 30) & 1) ^ fmt->qnan_msb_set;
3237 if (HOST_BITS_PER_LONG == 32)
3239 r->sig[SIGSZ-1] = sig_hi;
3240 r->sig[SIGSZ-2] = sig_lo;
3243 r->sig[SIGSZ-1] = (sig_hi << 31 << 1) | sig_lo;
3253 r->class = rvc_normal;
3255 r->exp = exp - 16383 + 1;
3256 if (HOST_BITS_PER_LONG == 32)
3258 r->sig[SIGSZ-1] = sig_hi;
3259 r->sig[SIGSZ-2] = sig_lo;
3262 r->sig[SIGSZ-1] = (sig_hi << 31 << 1) | sig_lo;
3267 decode_ieee_extended_128 (fmt, r, buf)
3268 const struct real_format *fmt;
3272 decode_ieee_extended (fmt, r, buf+!!FLOAT_WORDS_BIG_ENDIAN);
3275 const struct real_format ieee_extended_motorola_format =
3277 encode_ieee_extended,
3278 decode_ieee_extended,
3293 const struct real_format ieee_extended_intel_96_format =
3295 encode_ieee_extended,
3296 decode_ieee_extended,
3311 const struct real_format ieee_extended_intel_128_format =
3313 encode_ieee_extended_128,
3314 decode_ieee_extended_128,
3329 /* The following caters to i386 systems that set the rounding precision
3330 to 53 bits instead of 64, e.g. FreeBSD. */
3331 const struct real_format ieee_extended_intel_96_round_53_format =
3333 encode_ieee_extended,
3334 decode_ieee_extended,
3349 /* IBM 128-bit extended precision format: a pair of IEEE double precision
3350 numbers whose sum is equal to the extended precision value. The number
3351 with greater magnitude is first. This format has the same magnitude
3352 range as an IEEE double precision value, but effectively 106 bits of
3353 significand precision. Infinity and NaN are represented by their IEEE
3354 double precision value stored in the first number, the second number is
3355 ignored. Zeroes, Infinities, and NaNs are set in both doubles
3356 due to precedent. */
3358 static void encode_ibm_extended PARAMS ((const struct real_format *fmt,
3359 long *, const REAL_VALUE_TYPE *));
3360 static void decode_ibm_extended PARAMS ((const struct real_format *,
3361 REAL_VALUE_TYPE *, const long *));
3364 encode_ibm_extended (fmt, buf, r)
3365 const struct real_format *fmt;
3367 const REAL_VALUE_TYPE *r;
3369 REAL_VALUE_TYPE u, v;
3370 const struct real_format *base_fmt;
3372 base_fmt = fmt->qnan_msb_set ? &ieee_double_format : &mips_double_format;
3377 /* Both doubles have sign bit set. */
3378 buf[0] = FLOAT_WORDS_BIG_ENDIAN ? r->sign << 31 : 0;
3379 buf[1] = FLOAT_WORDS_BIG_ENDIAN ? 0 : r->sign << 31;
3386 /* Both doubles set to Inf / NaN. */
3387 encode_ieee_double (base_fmt, &buf[0], r);
3393 /* u = IEEE double precision portion of significand. */
3395 clear_significand_below (&u, SIGNIFICAND_BITS - 53);
3398 /* If the upper double is zero, we have a denormal double, so
3399 move it to the first double and leave the second as zero. */
3400 if (u.class == rvc_zero)
3408 /* v = remainder containing additional 53 bits of significand. */
3409 do_add (&v, r, &u, 1);
3410 round_for_format (base_fmt, &v);
3413 round_for_format (base_fmt, &u);
3415 encode_ieee_double (base_fmt, &buf[0], &u);
3416 encode_ieee_double (base_fmt, &buf[2], &v);
3425 decode_ibm_extended (fmt, r, buf)
3426 const struct real_format *fmt ATTRIBUTE_UNUSED;
3430 REAL_VALUE_TYPE u, v;
3431 const struct real_format *base_fmt;
3433 base_fmt = fmt->qnan_msb_set ? &ieee_double_format : &mips_double_format;
3434 decode_ieee_double (base_fmt, &u, &buf[0]);
3436 if (u.class != rvc_zero && u.class != rvc_inf && u.class != rvc_nan)
3438 decode_ieee_double (base_fmt, &v, &buf[2]);
3439 do_add (r, &u, &v, 0);
3445 const struct real_format ibm_extended_format =
3447 encode_ibm_extended,
3448 decode_ibm_extended,
3463 const struct real_format mips_extended_format =
3465 encode_ibm_extended,
3466 decode_ibm_extended,
3482 /* IEEE quad precision format. */
3484 static void encode_ieee_quad PARAMS ((const struct real_format *fmt,
3485 long *, const REAL_VALUE_TYPE *));
3486 static void decode_ieee_quad PARAMS ((const struct real_format *,
3487 REAL_VALUE_TYPE *, const long *));
3490 encode_ieee_quad (fmt, buf, r)
3491 const struct real_format *fmt;
3493 const REAL_VALUE_TYPE *r;
3495 unsigned long image3, image2, image1, image0, exp;
3496 bool denormal = (r->sig[SIGSZ-1] & SIG_MSB) == 0;
3499 image3 = r->sign << 31;
3504 rshift_significand (&u, r, SIGNIFICAND_BITS - 113);
3513 image3 |= 32767 << 16;
3516 image3 |= 0x7fffffff;
3517 image2 = 0xffffffff;
3518 image1 = 0xffffffff;
3519 image0 = 0xffffffff;
3526 image3 |= 32767 << 16;
3530 /* Don't use bits from the significand. The
3531 initialization above is right. */
3533 else if (HOST_BITS_PER_LONG == 32)
3538 image3 |= u.sig[3] & 0xffff;
3543 image1 = image0 >> 31 >> 1;
3545 image3 |= (image2 >> 31 >> 1) & 0xffff;
3546 image0 &= 0xffffffff;
3547 image2 &= 0xffffffff;
3549 if (r->signalling == fmt->qnan_msb_set)
3553 /* We overload qnan_msb_set here: it's only clear for
3554 mips_ieee_single, which wants all mantissa bits but the
3555 quiet/signalling one set in canonical NaNs (at least
3557 if (r->canonical && !fmt->qnan_msb_set)
3560 image2 = image1 = image0 = 0xffffffff;
3562 else if (((image3 & 0xffff) | image2 | image1 | image0) == 0)
3567 image3 |= 0x7fffffff;
3568 image2 = 0xffffffff;
3569 image1 = 0xffffffff;
3570 image0 = 0xffffffff;
3575 /* Recall that IEEE numbers are interpreted as 1.F x 2**exp,
3576 whereas the intermediate representation is 0.F x 2**exp.
3577 Which means we're off by one. */
3581 exp = r->exp + 16383 - 1;
3582 image3 |= exp << 16;
3584 if (HOST_BITS_PER_LONG == 32)
3589 image3 |= u.sig[3] & 0xffff;
3594 image1 = image0 >> 31 >> 1;
3596 image3 |= (image2 >> 31 >> 1) & 0xffff;
3597 image0 &= 0xffffffff;
3598 image2 &= 0xffffffff;
3606 if (FLOAT_WORDS_BIG_ENDIAN)
3623 decode_ieee_quad (fmt, r, buf)
3624 const struct real_format *fmt;
3628 unsigned long image3, image2, image1, image0;
3632 if (FLOAT_WORDS_BIG_ENDIAN)
3646 image0 &= 0xffffffff;
3647 image1 &= 0xffffffff;
3648 image2 &= 0xffffffff;
3650 sign = (image3 >> 31) & 1;
3651 exp = (image3 >> 16) & 0x7fff;
3654 memset (r, 0, sizeof (*r));
3658 if ((image3 | image2 | image1 | image0) && fmt->has_denorm)
3660 r->class = rvc_normal;
3663 r->exp = -16382 + (SIGNIFICAND_BITS - 112);
3664 if (HOST_BITS_PER_LONG == 32)
3673 r->sig[0] = (image1 << 31 << 1) | image0;
3674 r->sig[1] = (image3 << 31 << 1) | image2;
3679 else if (fmt->has_signed_zero)
3682 else if (exp == 32767 && (fmt->has_nans || fmt->has_inf))
3684 if (image3 | image2 | image1 | image0)
3688 r->signalling = ((image3 >> 15) & 1) ^ fmt->qnan_msb_set;
3690 if (HOST_BITS_PER_LONG == 32)
3699 r->sig[0] = (image1 << 31 << 1) | image0;
3700 r->sig[1] = (image3 << 31 << 1) | image2;
3702 lshift_significand (r, r, SIGNIFICAND_BITS - 113);
3712 r->class = rvc_normal;
3714 r->exp = exp - 16383 + 1;
3716 if (HOST_BITS_PER_LONG == 32)
3725 r->sig[0] = (image1 << 31 << 1) | image0;
3726 r->sig[1] = (image3 << 31 << 1) | image2;
3728 lshift_significand (r, r, SIGNIFICAND_BITS - 113);
3729 r->sig[SIGSZ-1] |= SIG_MSB;
3733 const struct real_format ieee_quad_format =
3751 const struct real_format mips_quad_format =
3769 /* Descriptions of VAX floating point formats can be found beginning at
3771 http://www.openvms.compaq.com:8000/73final/4515/4515pro_013.html#f_floating_point_format
3773 The thing to remember is that they're almost IEEE, except for word
3774 order, exponent bias, and the lack of infinities, nans, and denormals.
3776 We don't implement the H_floating format here, simply because neither
3777 the VAX or Alpha ports use it. */
3779 static void encode_vax_f PARAMS ((const struct real_format *fmt,
3780 long *, const REAL_VALUE_TYPE *));
3781 static void decode_vax_f PARAMS ((const struct real_format *,
3782 REAL_VALUE_TYPE *, const long *));
3783 static void encode_vax_d PARAMS ((const struct real_format *fmt,
3784 long *, const REAL_VALUE_TYPE *));
3785 static void decode_vax_d PARAMS ((const struct real_format *,
3786 REAL_VALUE_TYPE *, const long *));
3787 static void encode_vax_g PARAMS ((const struct real_format *fmt,
3788 long *, const REAL_VALUE_TYPE *));
3789 static void decode_vax_g PARAMS ((const struct real_format *,
3790 REAL_VALUE_TYPE *, const long *));
3793 encode_vax_f (fmt, buf, r)
3794 const struct real_format *fmt ATTRIBUTE_UNUSED;
3796 const REAL_VALUE_TYPE *r;
3798 unsigned long sign, exp, sig, image;
3800 sign = r->sign << 15;
3810 image = 0xffff7fff | sign;
3814 sig = (r->sig[SIGSZ-1] >> (HOST_BITS_PER_LONG - 24)) & 0x7fffff;
3817 image = (sig << 16) & 0xffff0000;
3831 decode_vax_f (fmt, r, buf)
3832 const struct real_format *fmt ATTRIBUTE_UNUSED;
3836 unsigned long image = buf[0] & 0xffffffff;
3837 int exp = (image >> 7) & 0xff;
3839 memset (r, 0, sizeof (*r));
3843 r->class = rvc_normal;
3844 r->sign = (image >> 15) & 1;
3847 image = ((image & 0x7f) << 16) | ((image >> 16) & 0xffff);
3848 r->sig[SIGSZ-1] = (image << (HOST_BITS_PER_LONG - 24)) | SIG_MSB;
3853 encode_vax_d (fmt, buf, r)
3854 const struct real_format *fmt ATTRIBUTE_UNUSED;
3856 const REAL_VALUE_TYPE *r;
3858 unsigned long image0, image1, sign = r->sign << 15;
3863 image0 = image1 = 0;
3868 image0 = 0xffff7fff | sign;
3869 image1 = 0xffffffff;
3873 /* Extract the significand into straight hi:lo. */
3874 if (HOST_BITS_PER_LONG == 64)
3876 image0 = r->sig[SIGSZ-1];
3877 image1 = (image0 >> (64 - 56)) & 0xffffffff;
3878 image0 = (image0 >> (64 - 56 + 1) >> 31) & 0x7fffff;
3882 image0 = r->sig[SIGSZ-1];
3883 image1 = r->sig[SIGSZ-2];
3884 image1 = (image0 << 24) | (image1 >> 8);
3885 image0 = (image0 >> 8) & 0xffffff;
3888 /* Rearrange the half-words of the significand to match the
3890 image0 = ((image0 << 16) | (image0 >> 16)) & 0xffff007f;
3891 image1 = ((image1 << 16) | (image1 >> 16)) & 0xffffffff;
3893 /* Add the sign and exponent. */
3895 image0 |= (r->exp + 128) << 7;
3902 if (FLOAT_WORDS_BIG_ENDIAN)
3903 buf[0] = image1, buf[1] = image0;
3905 buf[0] = image0, buf[1] = image1;
3909 decode_vax_d (fmt, r, buf)
3910 const struct real_format *fmt ATTRIBUTE_UNUSED;
3914 unsigned long image0, image1;
3917 if (FLOAT_WORDS_BIG_ENDIAN)
3918 image1 = buf[0], image0 = buf[1];
3920 image0 = buf[0], image1 = buf[1];
3921 image0 &= 0xffffffff;
3922 image1 &= 0xffffffff;
3924 exp = (image0 >> 7) & 0x7f;
3926 memset (r, 0, sizeof (*r));
3930 r->class = rvc_normal;
3931 r->sign = (image0 >> 15) & 1;
3934 /* Rearrange the half-words of the external format into
3935 proper ascending order. */
3936 image0 = ((image0 & 0x7f) << 16) | ((image0 >> 16) & 0xffff);
3937 image1 = ((image1 & 0xffff) << 16) | ((image1 >> 16) & 0xffff);
3939 if (HOST_BITS_PER_LONG == 64)
3941 image0 = (image0 << 31 << 1) | image1;
3944 r->sig[SIGSZ-1] = image0;
3948 r->sig[SIGSZ-1] = image0;
3949 r->sig[SIGSZ-2] = image1;
3950 lshift_significand (r, r, 2*HOST_BITS_PER_LONG - 56);
3951 r->sig[SIGSZ-1] |= SIG_MSB;
3957 encode_vax_g (fmt, buf, r)
3958 const struct real_format *fmt ATTRIBUTE_UNUSED;
3960 const REAL_VALUE_TYPE *r;
3962 unsigned long image0, image1, sign = r->sign << 15;
3967 image0 = image1 = 0;
3972 image0 = 0xffff7fff | sign;
3973 image1 = 0xffffffff;
3977 /* Extract the significand into straight hi:lo. */
3978 if (HOST_BITS_PER_LONG == 64)
3980 image0 = r->sig[SIGSZ-1];
3981 image1 = (image0 >> (64 - 53)) & 0xffffffff;
3982 image0 = (image0 >> (64 - 53 + 1) >> 31) & 0xfffff;
3986 image0 = r->sig[SIGSZ-1];
3987 image1 = r->sig[SIGSZ-2];
3988 image1 = (image0 << 21) | (image1 >> 11);
3989 image0 = (image0 >> 11) & 0xfffff;
3992 /* Rearrange the half-words of the significand to match the
3994 image0 = ((image0 << 16) | (image0 >> 16)) & 0xffff000f;
3995 image1 = ((image1 << 16) | (image1 >> 16)) & 0xffffffff;
3997 /* Add the sign and exponent. */
3999 image0 |= (r->exp + 1024) << 4;
4006 if (FLOAT_WORDS_BIG_ENDIAN)
4007 buf[0] = image1, buf[1] = image0;
4009 buf[0] = image0, buf[1] = image1;
4013 decode_vax_g (fmt, r, buf)
4014 const struct real_format *fmt ATTRIBUTE_UNUSED;
4018 unsigned long image0, image1;
4021 if (FLOAT_WORDS_BIG_ENDIAN)
4022 image1 = buf[0], image0 = buf[1];
4024 image0 = buf[0], image1 = buf[1];
4025 image0 &= 0xffffffff;
4026 image1 &= 0xffffffff;
4028 exp = (image0 >> 4) & 0x7ff;
4030 memset (r, 0, sizeof (*r));
4034 r->class = rvc_normal;
4035 r->sign = (image0 >> 15) & 1;
4036 r->exp = exp - 1024;
4038 /* Rearrange the half-words of the external format into
4039 proper ascending order. */
4040 image0 = ((image0 & 0xf) << 16) | ((image0 >> 16) & 0xffff);
4041 image1 = ((image1 & 0xffff) << 16) | ((image1 >> 16) & 0xffff);
4043 if (HOST_BITS_PER_LONG == 64)
4045 image0 = (image0 << 31 << 1) | image1;
4048 r->sig[SIGSZ-1] = image0;
4052 r->sig[SIGSZ-1] = image0;
4053 r->sig[SIGSZ-2] = image1;
4054 lshift_significand (r, r, 64 - 53);
4055 r->sig[SIGSZ-1] |= SIG_MSB;
4060 const struct real_format vax_f_format =
4078 const struct real_format vax_d_format =
4096 const struct real_format vax_g_format =
4114 /* A good reference for these can be found in chapter 9 of
4115 "ESA/390 Principles of Operation", IBM document number SA22-7201-01.
4116 An on-line version can be found here:
4118 http://publibz.boulder.ibm.com/cgi-bin/bookmgr_OS390/BOOKS/DZ9AR001/9.1?DT=19930923083613
4121 static void encode_i370_single PARAMS ((const struct real_format *fmt,
4122 long *, const REAL_VALUE_TYPE *));
4123 static void decode_i370_single PARAMS ((const struct real_format *,
4124 REAL_VALUE_TYPE *, const long *));
4125 static void encode_i370_double PARAMS ((const struct real_format *fmt,
4126 long *, const REAL_VALUE_TYPE *));
4127 static void decode_i370_double PARAMS ((const struct real_format *,
4128 REAL_VALUE_TYPE *, const long *));
4131 encode_i370_single (fmt, buf, r)
4132 const struct real_format *fmt ATTRIBUTE_UNUSED;
4134 const REAL_VALUE_TYPE *r;
4136 unsigned long sign, exp, sig, image;
4138 sign = r->sign << 31;
4148 image = 0x7fffffff | sign;
4152 sig = (r->sig[SIGSZ-1] >> (HOST_BITS_PER_LONG - 24)) & 0xffffff;
4153 exp = ((r->exp / 4) + 64) << 24;
4154 image = sign | exp | sig;
4165 decode_i370_single (fmt, r, buf)
4166 const struct real_format *fmt ATTRIBUTE_UNUSED;
4170 unsigned long sign, sig, image = buf[0];
4173 sign = (image >> 31) & 1;
4174 exp = (image >> 24) & 0x7f;
4175 sig = image & 0xffffff;
4177 memset (r, 0, sizeof (*r));
4181 r->class = rvc_normal;
4183 r->exp = (exp - 64) * 4;
4184 r->sig[SIGSZ-1] = sig << (HOST_BITS_PER_LONG - 24);
4190 encode_i370_double (fmt, buf, r)
4191 const struct real_format *fmt ATTRIBUTE_UNUSED;
4193 const REAL_VALUE_TYPE *r;
4195 unsigned long sign, exp, image_hi, image_lo;
4197 sign = r->sign << 31;
4202 image_hi = image_lo = 0;
4207 image_hi = 0x7fffffff | sign;
4208 image_lo = 0xffffffff;
4212 if (HOST_BITS_PER_LONG == 64)
4214 image_hi = r->sig[SIGSZ-1];
4215 image_lo = (image_hi >> (64 - 56)) & 0xffffffff;
4216 image_hi = (image_hi >> (64 - 56 + 1) >> 31) & 0xffffff;
4220 image_hi = r->sig[SIGSZ-1];
4221 image_lo = r->sig[SIGSZ-2];
4222 image_lo = (image_lo >> 8) | (image_hi << 24);
4226 exp = ((r->exp / 4) + 64) << 24;
4227 image_hi |= sign | exp;
4234 if (FLOAT_WORDS_BIG_ENDIAN)
4235 buf[0] = image_hi, buf[1] = image_lo;
4237 buf[0] = image_lo, buf[1] = image_hi;
4241 decode_i370_double (fmt, r, buf)
4242 const struct real_format *fmt ATTRIBUTE_UNUSED;
4246 unsigned long sign, image_hi, image_lo;
4249 if (FLOAT_WORDS_BIG_ENDIAN)
4250 image_hi = buf[0], image_lo = buf[1];
4252 image_lo = buf[0], image_hi = buf[1];
4254 sign = (image_hi >> 31) & 1;
4255 exp = (image_hi >> 24) & 0x7f;
4256 image_hi &= 0xffffff;
4257 image_lo &= 0xffffffff;
4259 memset (r, 0, sizeof (*r));
4261 if (exp || image_hi || image_lo)
4263 r->class = rvc_normal;
4265 r->exp = (exp - 64) * 4 + (SIGNIFICAND_BITS - 56);
4267 if (HOST_BITS_PER_LONG == 32)
4269 r->sig[0] = image_lo;
4270 r->sig[1] = image_hi;
4273 r->sig[0] = image_lo | (image_hi << 31 << 1);
4279 const struct real_format i370_single_format =
4292 false, /* ??? The encoding does allow for "unnormals". */
4293 false, /* ??? The encoding does allow for "unnormals". */
4297 const struct real_format i370_double_format =
4310 false, /* ??? The encoding does allow for "unnormals". */
4311 false, /* ??? The encoding does allow for "unnormals". */
4315 /* The "twos-complement" c4x format is officially defined as
4319 This is rather misleading. One must remember that F is signed.
4320 A better description would be
4322 x = -1**s * ((s + 1 + .f) * 2**e
4324 So if we have a (4 bit) fraction of .1000 with a sign bit of 1,
4325 that's -1 * (1+1+(-.5)) == -1.5. I think.
4327 The constructions here are taken from Tables 5-1 and 5-2 of the
4328 TMS320C4x User's Guide wherein step-by-step instructions for
4329 conversion from IEEE are presented. That's close enough to our
4330 internal representation so as to make things easy.
4332 See http://www-s.ti.com/sc/psheets/spru063c/spru063c.pdf */
4334 static void encode_c4x_single PARAMS ((const struct real_format *fmt,
4335 long *, const REAL_VALUE_TYPE *));
4336 static void decode_c4x_single PARAMS ((const struct real_format *,
4337 REAL_VALUE_TYPE *, const long *));
4338 static void encode_c4x_extended PARAMS ((const struct real_format *fmt,
4339 long *, const REAL_VALUE_TYPE *));
4340 static void decode_c4x_extended PARAMS ((const struct real_format *,
4341 REAL_VALUE_TYPE *, const long *));
4344 encode_c4x_single (fmt, buf, r)
4345 const struct real_format *fmt ATTRIBUTE_UNUSED;
4347 const REAL_VALUE_TYPE *r;
4349 unsigned long image, exp, sig;
4361 sig = 0x800000 - r->sign;
4366 sig = (r->sig[SIGSZ-1] >> (HOST_BITS_PER_LONG - 24)) & 0x7fffff;
4381 image = ((exp & 0xff) << 24) | (sig & 0xffffff);
4386 decode_c4x_single (fmt, r, buf)
4387 const struct real_format *fmt ATTRIBUTE_UNUSED;
4391 unsigned long image = buf[0];
4395 exp = (((image >> 24) & 0xff) ^ 0x80) - 0x80;
4396 sf = ((image & 0xffffff) ^ 0x800000) - 0x800000;
4398 memset (r, 0, sizeof (*r));
4402 r->class = rvc_normal;
4404 sig = sf & 0x7fffff;
4413 sig = (sig << (HOST_BITS_PER_LONG - 24)) | SIG_MSB;
4416 r->sig[SIGSZ-1] = sig;
4421 encode_c4x_extended (fmt, buf, r)
4422 const struct real_format *fmt ATTRIBUTE_UNUSED;
4424 const REAL_VALUE_TYPE *r;
4426 unsigned long exp, sig;
4438 sig = 0x80000000 - r->sign;
4444 sig = r->sig[SIGSZ-1];
4445 if (HOST_BITS_PER_LONG == 64)
4446 sig = sig >> 1 >> 31;
4463 exp = (exp & 0xff) << 24;
4466 if (FLOAT_WORDS_BIG_ENDIAN)
4467 buf[0] = exp, buf[1] = sig;
4469 buf[0] = sig, buf[0] = exp;
4473 decode_c4x_extended (fmt, r, buf)
4474 const struct real_format *fmt ATTRIBUTE_UNUSED;
4481 if (FLOAT_WORDS_BIG_ENDIAN)
4482 exp = buf[0], sf = buf[1];
4484 sf = buf[0], exp = buf[1];
4486 exp = (((exp >> 24) & 0xff) & 0x80) - 0x80;
4487 sf = ((sf & 0xffffffff) ^ 0x80000000) - 0x80000000;
4489 memset (r, 0, sizeof (*r));
4493 r->class = rvc_normal;
4495 sig = sf & 0x7fffffff;
4504 if (HOST_BITS_PER_LONG == 64)
4505 sig = sig << 1 << 31;
4509 r->sig[SIGSZ-1] = sig;
4513 const struct real_format c4x_single_format =
4531 const struct real_format c4x_extended_format =
4533 encode_c4x_extended,
4534 decode_c4x_extended,
4550 /* A synthetic "format" for internal arithmetic. It's the size of the
4551 internal significand minus the two bits needed for proper rounding.
4552 The encode and decode routines exist only to satisfy our paranoia
4555 static void encode_internal PARAMS ((const struct real_format *fmt,
4556 long *, const REAL_VALUE_TYPE *));
4557 static void decode_internal PARAMS ((const struct real_format *,
4558 REAL_VALUE_TYPE *, const long *));
4561 encode_internal (fmt, buf, r)
4562 const struct real_format *fmt ATTRIBUTE_UNUSED;
4564 const REAL_VALUE_TYPE *r;
4566 memcpy (buf, r, sizeof (*r));
4570 decode_internal (fmt, r, buf)
4571 const struct real_format *fmt ATTRIBUTE_UNUSED;
4575 memcpy (r, buf, sizeof (*r));
4578 const struct real_format real_internal_format =
4584 SIGNIFICAND_BITS - 2,
4585 SIGNIFICAND_BITS - 2,
4596 /* Set up default mode to format mapping for IEEE. Everyone else has
4597 to set these values in OVERRIDE_OPTIONS. */
4599 const struct real_format *real_format_for_mode[TFmode - QFmode + 1] =
4604 &ieee_single_format, /* SFmode */
4605 &ieee_double_format, /* DFmode */
4607 /* We explicitly don't handle XFmode. There are two formats,
4608 pretty much equally common. Choose one in OVERRIDE_OPTIONS. */
4610 &ieee_quad_format /* TFmode */
4614 /* Calculate the square root of X in mode MODE, and store the result
4615 in R. Return TRUE if the operation does not raise an exception.
4616 For details see "High Precision Division and Square Root",
4617 Alan H. Karp and Peter Markstein, HP Lab Report 93-93-42, June
4618 1993. http://www.hpl.hp.com/techreports/93/HPL-93-42.pdf. */
4621 real_sqrt (r, mode, x)
4623 enum machine_mode mode;
4624 const REAL_VALUE_TYPE *x;
4626 static REAL_VALUE_TYPE halfthree;
4627 static bool init = false;
4628 REAL_VALUE_TYPE h, t, i;
4631 /* sqrt(-0.0) is -0.0. */
4632 if (real_isnegzero (x))
4638 /* Negative arguments return NaN. */
4641 /* Mode is ignored for canonical NaN. */
4642 real_nan (r, "", 1, SFmode);
4646 /* Infinity and NaN return themselves. */
4647 if (real_isinf (x) || real_isnan (x))
4655 do_add (&halfthree, &dconst1, &dconsthalf, 0);
4659 /* Initial guess for reciprocal sqrt, i. */
4660 exp = real_exponent (x);
4661 real_ldexp (&i, &dconst1, -exp/2);
4663 /* Newton's iteration for reciprocal sqrt, i. */
4664 for (iter = 0; iter < 16; iter++)
4666 /* i(n+1) = i(n) * (1.5 - 0.5*i(n)*i(n)*x). */
4667 do_multiply (&t, x, &i);
4668 do_multiply (&h, &t, &i);
4669 do_multiply (&t, &h, &dconsthalf);
4670 do_add (&h, &halfthree, &t, 1);
4671 do_multiply (&t, &i, &h);
4673 /* Check for early convergence. */
4674 if (iter >= 6 && real_identical (&i, &t))
4677 /* ??? Unroll loop to avoid copying. */
4681 /* Final iteration: r = i*x + 0.5*i*x*(1.0 - i*(i*x)). */
4682 do_multiply (&t, x, &i);
4683 do_multiply (&h, &t, &i);
4684 do_add (&i, &dconst1, &h, 1);
4685 do_multiply (&h, &t, &i);
4686 do_multiply (&i, &dconsthalf, &h);
4687 do_add (&h, &t, &i, 0);
4689 /* ??? We need a Tuckerman test to get the last bit. */
4691 real_convert (r, mode, &h);
4695 /* Calculate X raised to the integer exponent N in mode MODE and store
4696 the result in R. Return true if the result may be inexact due to
4697 loss of precision. The algorithm is the classic "left-to-right binary
4698 method" described in section 4.6.3 of Donald Knuth's "Seminumerical
4699 Algorithms", "The Art of Computer Programming", Volume 2. */
4702 real_powi (r, mode, x, n)
4704 enum machine_mode mode;
4705 const REAL_VALUE_TYPE *x;
4708 unsigned HOST_WIDE_INT bit;
4710 bool inexact = false;
4722 /* Don't worry about overflow, from now on n is unsigned. */
4730 bit = (unsigned HOST_WIDE_INT) 1 << (HOST_BITS_PER_WIDE_INT - 1);
4731 for (i = 0; i < HOST_BITS_PER_WIDE_INT; i++)
4735 inexact |= do_multiply (&t, &t, &t);
4737 inexact |= do_multiply (&t, &t, x);
4745 inexact |= do_divide (&t, &dconst1, &t);
4747 real_convert (r, mode, &t);