1 /* real.c - software floating point emulation.
2 Copyright (C) 1993, 1994, 1995, 1996, 1997, 1998, 1999,
3 2000, 2002, 2003 Free Software Foundation, Inc.
4 Contributed by Stephen L. Moshier (moshier@world.std.com).
5 Re-written by Richard Henderson <rth@redhat.com>
7 This file is part of GCC.
9 GCC is free software; you can redistribute it and/or modify it under
10 the terms of the GNU General Public License as published by the Free
11 Software Foundation; either version 2, or (at your option) any later
14 GCC is distributed in the hope that it will be useful, but WITHOUT ANY
15 WARRANTY; without even the implied warranty of MERCHANTABILITY or
16 FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19 You should have received a copy of the GNU General Public License
20 along with GCC; see the file COPYING. If not, write to the Free
21 Software Foundation, 59 Temple Place - Suite 330, Boston, MA
26 #include "coretypes.h"
33 /* The floating point model used internally is not exactly IEEE 754
34 compliant, and close to the description in the ISO C standard,
35 section 5.2.4.2.2 Characteristics of floating types.
39 x = s * b^e * \sum_{k=1}^p f_k * b^{-k}
43 b = base or radix, here always 2
45 p = precision (the number of base-b digits in the significand)
46 f_k = the digits of the significand.
48 We differ from typical IEEE 754 encodings in that the entire
49 significand is fractional. Normalized significands are in the
52 A requirement of the model is that P be larger than than the
53 largest supported target floating-point type by at least 2 bits.
54 This gives us proper rounding when we truncate to the target type.
55 In addition, E must be large enough to hold the smallest supported
56 denormal number in a normalized form.
58 Both of these requirements are easily satisfied. The largest target
59 significand is 113 bits; we store at least 160. The smallest
60 denormal number fits in 17 exponent bits; we store 29.
62 Note that the decimal string conversion routines are sensitive to
63 rounding error. Since the raw arithmetic routines do not themselves
64 have guard digits or rounding, the computation of 10**exp can
65 accumulate more than a few digits of error. The previous incarnation
66 of real.c successfully used a 144 bit fraction; given the current
67 layout of REAL_VALUE_TYPE we're forced to expand to at least 160 bits.
69 Target floating point models that use base 16 instead of base 2
70 (i.e. IBM 370), are handled during round_for_format, in which we
71 canonicalize the exponent to be a multiple of 4 (log2(16)), and
72 adjust the significand to match. */
75 /* Used to classify two numbers simultaneously. */
76 #define CLASS2(A, B) ((A) << 2 | (B))
78 #if HOST_BITS_PER_LONG != 64 && HOST_BITS_PER_LONG != 32
79 #error "Some constant folding done by hand to avoid shift count warnings"
82 static void get_zero PARAMS ((REAL_VALUE_TYPE *, int));
83 static void get_canonical_qnan PARAMS ((REAL_VALUE_TYPE *, int));
84 static void get_canonical_snan PARAMS ((REAL_VALUE_TYPE *, int));
85 static void get_inf PARAMS ((REAL_VALUE_TYPE *, int));
86 static bool sticky_rshift_significand PARAMS ((REAL_VALUE_TYPE *,
87 const REAL_VALUE_TYPE *,
89 static void rshift_significand PARAMS ((REAL_VALUE_TYPE *,
90 const REAL_VALUE_TYPE *,
92 static void lshift_significand PARAMS ((REAL_VALUE_TYPE *,
93 const REAL_VALUE_TYPE *,
95 static void lshift_significand_1 PARAMS ((REAL_VALUE_TYPE *,
96 const REAL_VALUE_TYPE *));
97 static bool add_significands PARAMS ((REAL_VALUE_TYPE *r,
98 const REAL_VALUE_TYPE *,
99 const REAL_VALUE_TYPE *));
100 static bool sub_significands PARAMS ((REAL_VALUE_TYPE *,
101 const REAL_VALUE_TYPE *,
102 const REAL_VALUE_TYPE *, int));
103 static void neg_significand PARAMS ((REAL_VALUE_TYPE *,
104 const REAL_VALUE_TYPE *));
105 static int cmp_significands PARAMS ((const REAL_VALUE_TYPE *,
106 const REAL_VALUE_TYPE *));
107 static int cmp_significand_0 PARAMS ((const REAL_VALUE_TYPE *));
108 static void set_significand_bit PARAMS ((REAL_VALUE_TYPE *, unsigned int));
109 static void clear_significand_bit PARAMS ((REAL_VALUE_TYPE *, unsigned int));
110 static bool test_significand_bit PARAMS ((REAL_VALUE_TYPE *, unsigned int));
111 static void clear_significand_below PARAMS ((REAL_VALUE_TYPE *,
113 static bool div_significands PARAMS ((REAL_VALUE_TYPE *,
114 const REAL_VALUE_TYPE *,
115 const REAL_VALUE_TYPE *));
116 static void normalize PARAMS ((REAL_VALUE_TYPE *));
118 static void do_add PARAMS ((REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *,
119 const REAL_VALUE_TYPE *, int));
120 static void do_multiply PARAMS ((REAL_VALUE_TYPE *,
121 const REAL_VALUE_TYPE *,
122 const REAL_VALUE_TYPE *));
123 static void do_divide PARAMS ((REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *,
124 const REAL_VALUE_TYPE *));
125 static int do_compare PARAMS ((const REAL_VALUE_TYPE *,
126 const REAL_VALUE_TYPE *, int));
127 static void do_fix_trunc PARAMS ((REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *));
129 static unsigned long rtd_divmod PARAMS ((REAL_VALUE_TYPE *,
132 static const REAL_VALUE_TYPE * ten_to_ptwo PARAMS ((int));
133 static const REAL_VALUE_TYPE * ten_to_mptwo PARAMS ((int));
134 static const REAL_VALUE_TYPE * real_digit PARAMS ((int));
135 static void times_pten PARAMS ((REAL_VALUE_TYPE *, int));
137 static void round_for_format PARAMS ((const struct real_format *,
140 /* Initialize R with a positive zero. */
147 memset (r, 0, sizeof (*r));
151 /* Initialize R with the canonical quiet NaN. */
154 get_canonical_qnan (r, sign)
158 memset (r, 0, sizeof (*r));
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 /* Return R = A + (SUBTRACT_P ? -B : B). */
569 do_add (r, a, b, subtract_p)
571 const REAL_VALUE_TYPE *a, *b;
576 bool inexact = false;
578 /* Determine if we need to add or subtract. */
580 subtract_p = (sign ^ b->sign) ^ subtract_p;
582 switch (CLASS2 (a->class, b->class))
584 case CLASS2 (rvc_zero, rvc_zero):
585 /* -0 + -0 = -0, -0 - +0 = -0; all other cases yield +0. */
586 get_zero (r, sign & !subtract_p);
589 case CLASS2 (rvc_zero, rvc_normal):
590 case CLASS2 (rvc_zero, rvc_inf):
591 case CLASS2 (rvc_zero, rvc_nan):
593 case CLASS2 (rvc_normal, rvc_nan):
594 case CLASS2 (rvc_inf, rvc_nan):
595 case CLASS2 (rvc_nan, rvc_nan):
596 /* ANY + NaN = NaN. */
597 case CLASS2 (rvc_normal, rvc_inf):
600 r->sign = sign ^ subtract_p;
603 case CLASS2 (rvc_normal, rvc_zero):
604 case CLASS2 (rvc_inf, rvc_zero):
605 case CLASS2 (rvc_nan, rvc_zero):
607 case CLASS2 (rvc_nan, rvc_normal):
608 case CLASS2 (rvc_nan, rvc_inf):
609 /* NaN + ANY = NaN. */
610 case CLASS2 (rvc_inf, rvc_normal):
615 case CLASS2 (rvc_inf, rvc_inf):
617 /* Inf - Inf = NaN. */
618 get_canonical_qnan (r, 0);
620 /* Inf + Inf = Inf. */
624 case CLASS2 (rvc_normal, rvc_normal):
631 /* Swap the arguments such that A has the larger exponent. */
632 dexp = a->exp - b->exp;
635 const REAL_VALUE_TYPE *t;
642 /* If the exponents are not identical, we need to shift the
643 significand of B down. */
646 /* If the exponents are too far apart, the significands
647 do not overlap, which makes the subtraction a noop. */
648 if (dexp >= SIGNIFICAND_BITS)
655 inexact |= sticky_rshift_significand (&t, b, dexp);
661 if (sub_significands (r, a, b, inexact))
663 /* We got a borrow out of the subtraction. That means that
664 A and B had the same exponent, and B had the larger
665 significand. We need to swap the sign and negate the
668 neg_significand (r, r);
673 if (add_significands (r, a, b))
675 /* We got carry out of the addition. This means we need to
676 shift the significand back down one bit and increase the
678 inexact |= sticky_rshift_significand (r, r, 1);
679 r->sig[SIGSZ-1] |= SIG_MSB;
688 r->class = rvc_normal;
692 /* Re-normalize the result. */
695 /* Special case: if the subtraction results in zero, the result
697 if (r->class == rvc_zero)
700 r->sig[0] |= inexact;
703 /* Return R = A * B. */
706 do_multiply (r, a, b)
708 const REAL_VALUE_TYPE *a, *b;
710 REAL_VALUE_TYPE u, t, *rr;
711 unsigned int i, j, k;
712 int sign = a->sign ^ b->sign;
714 switch (CLASS2 (a->class, b->class))
716 case CLASS2 (rvc_zero, rvc_zero):
717 case CLASS2 (rvc_zero, rvc_normal):
718 case CLASS2 (rvc_normal, rvc_zero):
719 /* +-0 * ANY = 0 with appropriate sign. */
723 case CLASS2 (rvc_zero, rvc_nan):
724 case CLASS2 (rvc_normal, rvc_nan):
725 case CLASS2 (rvc_inf, rvc_nan):
726 case CLASS2 (rvc_nan, rvc_nan):
727 /* ANY * NaN = NaN. */
732 case CLASS2 (rvc_nan, rvc_zero):
733 case CLASS2 (rvc_nan, rvc_normal):
734 case CLASS2 (rvc_nan, rvc_inf):
735 /* NaN * ANY = NaN. */
740 case CLASS2 (rvc_zero, rvc_inf):
741 case CLASS2 (rvc_inf, rvc_zero):
743 get_canonical_qnan (r, sign);
746 case CLASS2 (rvc_inf, rvc_inf):
747 case CLASS2 (rvc_normal, rvc_inf):
748 case CLASS2 (rvc_inf, rvc_normal):
749 /* Inf * Inf = Inf, R * Inf = Inf */
754 case CLASS2 (rvc_normal, rvc_normal):
761 if (r == a || r == b)
767 /* Collect all the partial products. Since we don't have sure access
768 to a widening multiply, we split each long into two half-words.
770 Consider the long-hand form of a four half-word multiplication:
780 We construct partial products of the widened half-word products
781 that are known to not overlap, e.g. DF+DH. Each such partial
782 product is given its proper exponent, which allows us to sum them
783 and obtain the finished product. */
785 for (i = 0; i < SIGSZ * 2; ++i)
787 unsigned long ai = a->sig[i / 2];
789 ai >>= HOST_BITS_PER_LONG / 2;
791 ai &= ((unsigned long)1 << (HOST_BITS_PER_LONG / 2)) - 1;
796 for (j = 0; j < 2; ++j)
798 int exp = (a->exp - (2*SIGSZ-1-i)*(HOST_BITS_PER_LONG/2)
799 + (b->exp - (1-j)*(HOST_BITS_PER_LONG/2)));
804 /* Would underflow to zero, which we shouldn't bother adding. */
807 u.class = rvc_normal;
811 for (k = j; k < SIGSZ * 2; k += 2)
813 unsigned long bi = b->sig[k / 2];
815 bi >>= HOST_BITS_PER_LONG / 2;
817 bi &= ((unsigned long)1 << (HOST_BITS_PER_LONG / 2)) - 1;
819 u.sig[k / 2] = ai * bi;
823 do_add (rr, rr, &u, 0);
832 /* Return R = A / B. */
837 const REAL_VALUE_TYPE *a, *b;
839 int exp, sign = a->sign ^ b->sign;
840 REAL_VALUE_TYPE t, *rr;
843 switch (CLASS2 (a->class, b->class))
845 case CLASS2 (rvc_zero, rvc_zero):
847 case CLASS2 (rvc_inf, rvc_inf):
848 /* Inf / Inf = NaN. */
849 get_canonical_qnan (r, sign);
852 case CLASS2 (rvc_zero, rvc_normal):
853 case CLASS2 (rvc_zero, rvc_inf):
855 case CLASS2 (rvc_normal, rvc_inf):
861 case CLASS2 (rvc_normal, rvc_zero):
863 case CLASS2 (rvc_inf, rvc_zero):
868 case CLASS2 (rvc_zero, rvc_nan):
869 case CLASS2 (rvc_normal, rvc_nan):
870 case CLASS2 (rvc_inf, rvc_nan):
871 case CLASS2 (rvc_nan, rvc_nan):
872 /* ANY / NaN = NaN. */
877 case CLASS2 (rvc_nan, rvc_zero):
878 case CLASS2 (rvc_nan, rvc_normal):
879 case CLASS2 (rvc_nan, rvc_inf):
880 /* NaN / ANY = NaN. */
885 case CLASS2 (rvc_inf, rvc_normal):
891 case CLASS2 (rvc_normal, rvc_normal):
898 if (r == a || r == b)
903 rr->class = rvc_normal;
906 exp = a->exp - b->exp + 1;
913 inexact = div_significands (rr, a, b);
915 /* Re-normalize the result. */
917 rr->sig[0] |= inexact;
923 /* Return a tri-state comparison of A vs B. Return NAN_RESULT if
924 one of the two operands is a NaN. */
927 do_compare (a, b, nan_result)
928 const REAL_VALUE_TYPE *a, *b;
933 switch (CLASS2 (a->class, b->class))
935 case CLASS2 (rvc_zero, rvc_zero):
936 /* Sign of zero doesn't matter for compares. */
939 case CLASS2 (rvc_inf, rvc_zero):
940 case CLASS2 (rvc_inf, rvc_normal):
941 case CLASS2 (rvc_normal, rvc_zero):
942 return (a->sign ? -1 : 1);
944 case CLASS2 (rvc_inf, rvc_inf):
945 return -a->sign - -b->sign;
947 case CLASS2 (rvc_zero, rvc_normal):
948 case CLASS2 (rvc_zero, rvc_inf):
949 case CLASS2 (rvc_normal, rvc_inf):
950 return (b->sign ? 1 : -1);
952 case CLASS2 (rvc_zero, rvc_nan):
953 case CLASS2 (rvc_normal, rvc_nan):
954 case CLASS2 (rvc_inf, rvc_nan):
955 case CLASS2 (rvc_nan, rvc_nan):
956 case CLASS2 (rvc_nan, rvc_zero):
957 case CLASS2 (rvc_nan, rvc_normal):
958 case CLASS2 (rvc_nan, rvc_inf):
961 case CLASS2 (rvc_normal, rvc_normal):
968 if (a->sign != b->sign)
969 return -a->sign - -b->sign;
973 else if (a->exp < b->exp)
976 ret = cmp_significands (a, b);
978 return (a->sign ? -ret : ret);
981 /* Return A truncated to an integral value toward zero. */
986 const REAL_VALUE_TYPE *a;
999 get_zero (r, r->sign);
1000 else if (r->exp < SIGNIFICAND_BITS)
1001 clear_significand_below (r, SIGNIFICAND_BITS - r->exp);
1009 /* Perform the binary or unary operation described by CODE.
1010 For a unary operation, leave OP1 NULL. */
1013 real_arithmetic (r, icode, op0, op1)
1016 const REAL_VALUE_TYPE *op0, *op1;
1018 enum tree_code code = icode;
1023 do_add (r, op0, op1, 0);
1027 do_add (r, op0, op1, 1);
1031 do_multiply (r, op0, op1);
1035 do_divide (r, op0, op1);
1039 if (op1->class == rvc_nan)
1041 else if (do_compare (op0, op1, -1) < 0)
1048 if (op1->class == rvc_nan)
1050 else if (do_compare (op0, op1, 1) < 0)
1066 case FIX_TRUNC_EXPR:
1067 do_fix_trunc (r, op0);
1075 /* Legacy. Similar, but return the result directly. */
1078 real_arithmetic2 (icode, op0, op1)
1080 const REAL_VALUE_TYPE *op0, *op1;
1083 real_arithmetic (&r, icode, op0, op1);
1088 real_compare (icode, op0, op1)
1090 const REAL_VALUE_TYPE *op0, *op1;
1092 enum tree_code code = icode;
1097 return do_compare (op0, op1, 1) < 0;
1099 return do_compare (op0, op1, 1) <= 0;
1101 return do_compare (op0, op1, -1) > 0;
1103 return do_compare (op0, op1, -1) >= 0;
1105 return do_compare (op0, op1, -1) == 0;
1107 return do_compare (op0, op1, -1) != 0;
1108 case UNORDERED_EXPR:
1109 return op0->class == rvc_nan || op1->class == rvc_nan;
1111 return op0->class != rvc_nan && op1->class != rvc_nan;
1113 return do_compare (op0, op1, -1) < 0;
1115 return do_compare (op0, op1, -1) <= 0;
1117 return do_compare (op0, op1, 1) > 0;
1119 return do_compare (op0, op1, 1) >= 0;
1121 return do_compare (op0, op1, 0) == 0;
1128 /* Return floor log2(R). */
1132 const REAL_VALUE_TYPE *r;
1140 return (unsigned int)-1 >> 1;
1148 /* R = OP0 * 2**EXP. */
1151 real_ldexp (r, op0, exp)
1153 const REAL_VALUE_TYPE *op0;
1167 get_inf (r, r->sign);
1168 else if (exp < -MAX_EXP)
1169 get_zero (r, r->sign);
1179 /* Determine whether a floating-point value X is infinite. */
1183 const REAL_VALUE_TYPE *r;
1185 return (r->class == rvc_inf);
1188 /* Determine whether a floating-point value X is a NaN. */
1192 const REAL_VALUE_TYPE *r;
1194 return (r->class == rvc_nan);
1197 /* Determine whether a floating-point value X is negative. */
1201 const REAL_VALUE_TYPE *r;
1206 /* Determine whether a floating-point value X is minus zero. */
1210 const REAL_VALUE_TYPE *r;
1212 return r->sign && r->class == rvc_zero;
1215 /* Compare two floating-point objects for bitwise identity. */
1218 real_identical (a, b)
1219 const REAL_VALUE_TYPE *a, *b;
1223 if (a->class != b->class)
1225 if (a->sign != b->sign)
1235 if (a->exp != b->exp)
1240 if (a->signalling != b->signalling)
1242 /* The significand is ignored for canonical NaNs. */
1243 if (a->canonical || b->canonical)
1244 return a->canonical == b->canonical;
1251 for (i = 0; i < SIGSZ; ++i)
1252 if (a->sig[i] != b->sig[i])
1258 /* Try to change R into its exact multiplicative inverse in machine
1259 mode MODE. Return true if successful. */
1262 exact_real_inverse (mode, r)
1263 enum machine_mode mode;
1266 const REAL_VALUE_TYPE *one = real_digit (1);
1270 if (r->class != rvc_normal)
1273 /* Check for a power of two: all significand bits zero except the MSB. */
1274 for (i = 0; i < SIGSZ-1; ++i)
1277 if (r->sig[SIGSZ-1] != SIG_MSB)
1280 /* Find the inverse and truncate to the required mode. */
1281 do_divide (&u, one, r);
1282 real_convert (&u, mode, &u);
1284 /* The rounding may have overflowed. */
1285 if (u.class != rvc_normal)
1287 for (i = 0; i < SIGSZ-1; ++i)
1290 if (u.sig[SIGSZ-1] != SIG_MSB)
1297 /* Render R as an integer. */
1301 const REAL_VALUE_TYPE *r;
1303 unsigned HOST_WIDE_INT i;
1314 i = (unsigned HOST_WIDE_INT) 1 << (HOST_BITS_PER_WIDE_INT - 1);
1322 /* Only force overflow for unsigned overflow. Signed overflow is
1323 undefined, so it doesn't matter what we return, and some callers
1324 expect to be able to use this routine for both signed and
1325 unsigned conversions. */
1326 if (r->exp > HOST_BITS_PER_WIDE_INT)
1329 if (HOST_BITS_PER_WIDE_INT == HOST_BITS_PER_LONG)
1330 i = r->sig[SIGSZ-1];
1331 else if (HOST_BITS_PER_WIDE_INT == 2*HOST_BITS_PER_LONG)
1333 i = r->sig[SIGSZ-1];
1334 i = i << (HOST_BITS_PER_LONG - 1) << 1;
1335 i |= r->sig[SIGSZ-2];
1340 i >>= HOST_BITS_PER_WIDE_INT - r->exp;
1351 /* Likewise, but to an integer pair, HI+LOW. */
1354 real_to_integer2 (plow, phigh, r)
1355 HOST_WIDE_INT *plow, *phigh;
1356 const REAL_VALUE_TYPE *r;
1359 HOST_WIDE_INT low, high;
1372 high = (unsigned HOST_WIDE_INT) 1 << (HOST_BITS_PER_WIDE_INT - 1);
1386 /* Only force overflow for unsigned overflow. Signed overflow is
1387 undefined, so it doesn't matter what we return, and some callers
1388 expect to be able to use this routine for both signed and
1389 unsigned conversions. */
1390 if (exp > 2*HOST_BITS_PER_WIDE_INT)
1393 rshift_significand (&t, r, 2*HOST_BITS_PER_WIDE_INT - exp);
1394 if (HOST_BITS_PER_WIDE_INT == HOST_BITS_PER_LONG)
1396 high = t.sig[SIGSZ-1];
1397 low = t.sig[SIGSZ-2];
1399 else if (HOST_BITS_PER_WIDE_INT == 2*HOST_BITS_PER_LONG)
1401 high = t.sig[SIGSZ-1];
1402 high = high << (HOST_BITS_PER_LONG - 1) << 1;
1403 high |= t.sig[SIGSZ-2];
1405 low = t.sig[SIGSZ-3];
1406 low = low << (HOST_BITS_PER_LONG - 1) << 1;
1407 low |= t.sig[SIGSZ-4];
1417 low = -low, high = ~high;
1429 /* A subroutine of real_to_decimal. Compute the quotient and remainder
1430 of NUM / DEN. Return the quotient and place the remainder in NUM.
1431 It is expected that NUM / DEN are close enough that the quotient is
1434 static unsigned long
1435 rtd_divmod (num, den)
1436 REAL_VALUE_TYPE *num, *den;
1438 unsigned long q, msb;
1439 int expn = num->exp, expd = den->exp;
1448 msb = num->sig[SIGSZ-1] & SIG_MSB;
1450 lshift_significand_1 (num, num);
1452 if (msb || cmp_significands (num, den) >= 0)
1454 sub_significands (num, num, den, 0);
1458 while (--expn >= expd);
1466 /* Render R as a decimal floating point constant. Emit DIGITS significant
1467 digits in the result, bounded by BUF_SIZE. If DIGITS is 0, choose the
1468 maximum for the representation. If CROP_TRAILING_ZEROS, strip trailing
1471 #define M_LOG10_2 0.30102999566398119521
1474 real_to_decimal (str, r_orig, buf_size, digits, crop_trailing_zeros)
1476 const REAL_VALUE_TYPE *r_orig;
1477 size_t buf_size, digits;
1478 int crop_trailing_zeros;
1480 const REAL_VALUE_TYPE *one, *ten;
1481 REAL_VALUE_TYPE r, pten, u, v;
1482 int dec_exp, cmp_one, digit;
1484 char *p, *first, *last;
1491 strcpy (str, (r.sign ? "-0.0" : "0.0"));
1496 strcpy (str, (r.sign ? "-Inf" : "+Inf"));
1499 /* ??? Print the significand as well, if not canonical? */
1500 strcpy (str, (r.sign ? "-NaN" : "+NaN"));
1506 /* Bound the number of digits printed by the size of the representation. */
1507 max_digits = SIGNIFICAND_BITS * M_LOG10_2;
1508 if (digits == 0 || digits > max_digits)
1509 digits = max_digits;
1511 /* Estimate the decimal exponent, and compute the length of the string it
1512 will print as. Be conservative and add one to account for possible
1513 overflow or rounding error. */
1514 dec_exp = r.exp * M_LOG10_2;
1515 for (max_digits = 1; dec_exp ; max_digits++)
1518 /* Bound the number of digits printed by the size of the output buffer. */
1519 max_digits = buf_size - 1 - 1 - 2 - max_digits - 1;
1520 if (max_digits > buf_size)
1522 if (digits > max_digits)
1523 digits = max_digits;
1525 one = real_digit (1);
1526 ten = ten_to_ptwo (0);
1534 cmp_one = do_compare (&r, one, 0);
1539 /* Number is greater than one. Convert significand to an integer
1540 and strip trailing decimal zeros. */
1543 u.exp = SIGNIFICAND_BITS - 1;
1545 /* Largest M, such that 10**2**M fits within SIGNIFICAND_BITS. */
1546 m = floor_log2 (max_digits);
1548 /* Iterate over the bits of the possible powers of 10 that might
1549 be present in U and eliminate them. That is, if we find that
1550 10**2**M divides U evenly, keep the division and increase
1556 do_divide (&t, &u, ten_to_ptwo (m));
1557 do_fix_trunc (&v, &t);
1558 if (cmp_significands (&v, &t) == 0)
1566 /* Revert the scaling to integer that we performed earlier. */
1567 u.exp += r.exp - (SIGNIFICAND_BITS - 1);
1570 /* Find power of 10. Do this by dividing out 10**2**M when
1571 this is larger than the current remainder. Fill PTEN with
1572 the power of 10 that we compute. */
1575 m = floor_log2 ((int)(r.exp * M_LOG10_2)) + 1;
1578 const REAL_VALUE_TYPE *ptentwo = ten_to_ptwo (m);
1579 if (do_compare (&u, ptentwo, 0) >= 0)
1581 do_divide (&u, &u, ptentwo);
1582 do_multiply (&pten, &pten, ptentwo);
1589 /* We managed to divide off enough tens in the above reduction
1590 loop that we've now got a negative exponent. Fall into the
1591 less-than-one code to compute the proper value for PTEN. */
1598 /* Number is less than one. Pad significand with leading
1604 /* Stop if we'd shift bits off the bottom. */
1608 do_multiply (&u, &v, ten);
1610 /* Stop if we're now >= 1. */
1619 /* Find power of 10. Do this by multiplying in P=10**2**M when
1620 the current remainder is smaller than 1/P. Fill PTEN with the
1621 power of 10 that we compute. */
1622 m = floor_log2 ((int)(-r.exp * M_LOG10_2)) + 1;
1625 const REAL_VALUE_TYPE *ptentwo = ten_to_ptwo (m);
1626 const REAL_VALUE_TYPE *ptenmtwo = ten_to_mptwo (m);
1628 if (do_compare (&v, ptenmtwo, 0) <= 0)
1630 do_multiply (&v, &v, ptentwo);
1631 do_multiply (&pten, &pten, ptentwo);
1637 /* Invert the positive power of 10 that we've collected so far. */
1638 do_divide (&pten, one, &pten);
1646 /* At this point, PTEN should contain the nearest power of 10 smaller
1647 than R, such that this division produces the first digit.
1649 Using a divide-step primitive that returns the complete integral
1650 remainder avoids the rounding error that would be produced if
1651 we were to use do_divide here and then simply multiply by 10 for
1652 each subsequent digit. */
1654 digit = rtd_divmod (&r, &pten);
1656 /* Be prepared for error in that division via underflow ... */
1657 if (digit == 0 && cmp_significand_0 (&r))
1659 /* Multiply by 10 and try again. */
1660 do_multiply (&r, &r, ten);
1661 digit = rtd_divmod (&r, &pten);
1667 /* ... or overflow. */
1675 else if (digit > 10)
1680 /* Generate subsequent digits. */
1681 while (--digits > 0)
1683 do_multiply (&r, &r, ten);
1684 digit = rtd_divmod (&r, &pten);
1689 /* Generate one more digit with which to do rounding. */
1690 do_multiply (&r, &r, ten);
1691 digit = rtd_divmod (&r, &pten);
1693 /* Round the result. */
1696 /* Round to nearest. If R is nonzero there are additional
1697 nonzero digits to be extracted. */
1698 if (cmp_significand_0 (&r))
1700 /* Round to even. */
1701 else if ((p[-1] - '0') & 1)
1718 /* Carry out of the first digit. This means we had all 9's and
1719 now have all 0's. "Prepend" a 1 by overwriting the first 0. */
1727 /* Insert the decimal point. */
1728 first[0] = first[1];
1731 /* If requested, drop trailing zeros. Never crop past "1.0". */
1732 if (crop_trailing_zeros)
1733 while (last > first + 3 && last[-1] == '0')
1736 /* Append the exponent. */
1737 sprintf (last, "e%+d", dec_exp);
1740 /* Render R as a hexadecimal floating point constant. Emit DIGITS
1741 significant digits in the result, bounded by BUF_SIZE. If DIGITS is 0,
1742 choose the maximum for the representation. If CROP_TRAILING_ZEROS,
1743 strip trailing zeros. */
1746 real_to_hexadecimal (str, r, buf_size, digits, crop_trailing_zeros)
1748 const REAL_VALUE_TYPE *r;
1749 size_t buf_size, digits;
1750 int crop_trailing_zeros;
1752 int i, j, exp = r->exp;
1765 strcpy (str, (r->sign ? "-Inf" : "+Inf"));
1768 /* ??? Print the significand as well, if not canonical? */
1769 strcpy (str, (r->sign ? "-NaN" : "+NaN"));
1776 digits = SIGNIFICAND_BITS / 4;
1778 /* Bound the number of digits printed by the size of the output buffer. */
1780 sprintf (exp_buf, "p%+d", exp);
1781 max_digits = buf_size - strlen (exp_buf) - r->sign - 4 - 1;
1782 if (max_digits > buf_size)
1784 if (digits > max_digits)
1785 digits = max_digits;
1796 for (i = SIGSZ - 1; i >= 0; --i)
1797 for (j = HOST_BITS_PER_LONG - 4; j >= 0; j -= 4)
1799 *p++ = "0123456789abcdef"[(r->sig[i] >> j) & 15];
1805 if (crop_trailing_zeros)
1806 while (p > first + 1 && p[-1] == '0')
1809 sprintf (p, "p%+d", exp);
1812 /* Initialize R from a decimal or hexadecimal string. The string is
1813 assumed to have been syntax checked already. */
1816 real_from_string (r, str)
1830 else if (*str == '+')
1833 if (str[0] == '0' && str[1] == 'x')
1835 /* Hexadecimal floating point. */
1836 int pos = SIGNIFICAND_BITS - 4, d;
1844 d = hex_value (*str);
1849 r->sig[pos / HOST_BITS_PER_LONG]
1850 |= (unsigned long) d << (pos % HOST_BITS_PER_LONG);
1859 if (pos == SIGNIFICAND_BITS - 4)
1866 d = hex_value (*str);
1871 r->sig[pos / HOST_BITS_PER_LONG]
1872 |= (unsigned long) d << (pos % HOST_BITS_PER_LONG);
1878 if (*str == 'p' || *str == 'P')
1880 bool exp_neg = false;
1888 else if (*str == '+')
1892 while (ISDIGIT (*str))
1898 /* Overflowed the exponent. */
1912 r->class = rvc_normal;
1919 /* Decimal floating point. */
1920 const REAL_VALUE_TYPE *ten = ten_to_ptwo (0);
1925 while (ISDIGIT (*str))
1928 do_multiply (r, r, ten);
1930 do_add (r, r, real_digit (d), 0);
1935 if (r->class == rvc_zero)
1940 while (ISDIGIT (*str))
1943 do_multiply (r, r, ten);
1945 do_add (r, r, real_digit (d), 0);
1950 if (*str == 'e' || *str == 'E')
1952 bool exp_neg = false;
1960 else if (*str == '+')
1964 while (ISDIGIT (*str))
1970 /* Overflowed the exponent. */
1984 times_pten (r, exp);
1999 /* Legacy. Similar, but return the result directly. */
2002 real_from_string2 (s, mode)
2004 enum machine_mode mode;
2008 real_from_string (&r, s);
2009 if (mode != VOIDmode)
2010 real_convert (&r, mode, &r);
2015 /* Initialize R from the integer pair HIGH+LOW. */
2018 real_from_integer (r, mode, low, high, unsigned_p)
2020 enum machine_mode mode;
2021 unsigned HOST_WIDE_INT low;
2025 if (low == 0 && high == 0)
2029 r->class = rvc_normal;
2030 r->sign = high < 0 && !unsigned_p;
2031 r->exp = 2 * HOST_BITS_PER_WIDE_INT;
2042 if (HOST_BITS_PER_LONG == HOST_BITS_PER_WIDE_INT)
2044 r->sig[SIGSZ-1] = high;
2045 r->sig[SIGSZ-2] = low;
2046 memset (r->sig, 0, sizeof(long)*(SIGSZ-2));
2048 else if (HOST_BITS_PER_LONG*2 == HOST_BITS_PER_WIDE_INT)
2050 r->sig[SIGSZ-1] = high >> (HOST_BITS_PER_LONG - 1) >> 1;
2051 r->sig[SIGSZ-2] = high;
2052 r->sig[SIGSZ-3] = low >> (HOST_BITS_PER_LONG - 1) >> 1;
2053 r->sig[SIGSZ-4] = low;
2055 memset (r->sig, 0, sizeof(long)*(SIGSZ-4));
2063 if (mode != VOIDmode)
2064 real_convert (r, mode, r);
2067 /* Returns 10**2**N. */
2069 static const REAL_VALUE_TYPE *
2073 static REAL_VALUE_TYPE tens[EXP_BITS];
2075 if (n < 0 || n >= EXP_BITS)
2078 if (tens[n].class == rvc_zero)
2080 if (n < (HOST_BITS_PER_WIDE_INT == 64 ? 5 : 4))
2082 HOST_WIDE_INT t = 10;
2085 for (i = 0; i < n; ++i)
2088 real_from_integer (&tens[n], VOIDmode, t, 0, 1);
2092 const REAL_VALUE_TYPE *t = ten_to_ptwo (n - 1);
2093 do_multiply (&tens[n], t, t);
2100 /* Returns 10**(-2**N). */
2102 static const REAL_VALUE_TYPE *
2106 static REAL_VALUE_TYPE tens[EXP_BITS];
2108 if (n < 0 || n >= EXP_BITS)
2111 if (tens[n].class == rvc_zero)
2112 do_divide (&tens[n], real_digit (1), ten_to_ptwo (n));
2119 static const REAL_VALUE_TYPE *
2123 static REAL_VALUE_TYPE num[10];
2128 if (n > 0 && num[n].class == rvc_zero)
2129 real_from_integer (&num[n], VOIDmode, n, 0, 1);
2134 /* Multiply R by 10**EXP. */
2141 REAL_VALUE_TYPE pten, *rr;
2142 bool negative = (exp < 0);
2148 pten = *real_digit (1);
2154 for (i = 0; exp > 0; ++i, exp >>= 1)
2156 do_multiply (rr, rr, ten_to_ptwo (i));
2159 do_divide (r, r, &pten);
2162 /* Fills R with +Inf. */
2171 /* Fills R with a NaN whose significand is described by STR. If QUIET,
2172 we force a QNaN, else we force an SNaN. The string, if not empty,
2173 is parsed as a number and placed in the significand. Return true
2174 if the string was successfully parsed. */
2177 real_nan (r, str, quiet, mode)
2181 enum machine_mode mode;
2183 const struct real_format *fmt;
2185 fmt = real_format_for_mode[mode - QFmode];
2192 get_canonical_qnan (r, 0);
2194 get_canonical_snan (r, 0);
2201 memset (r, 0, sizeof (*r));
2204 /* Parse akin to strtol into the significand of R. */
2206 while (ISSPACE (*str))
2210 else if (*str == '+')
2220 while ((d = hex_value (*str)) < base)
2227 lshift_significand (r, r, 3);
2230 lshift_significand (r, r, 4);
2233 lshift_significand_1 (&u, r);
2234 lshift_significand (r, r, 3);
2235 add_significands (r, r, &u);
2243 add_significands (r, r, &u);
2248 /* Must have consumed the entire string for success. */
2252 /* Shift the significand into place such that the bits
2253 are in the most significant bits for the format. */
2254 lshift_significand (r, r, SIGNIFICAND_BITS - fmt->pnan);
2256 /* Our MSB is always unset for NaNs. */
2257 r->sig[SIGSZ-1] &= ~SIG_MSB;
2259 /* Force quiet or signalling NaN. */
2260 r->signalling = !quiet;
2266 /* Fills R with 2**N. */
2273 memset (r, 0, sizeof (*r));
2278 else if (n < -MAX_EXP)
2282 r->class = rvc_normal;
2284 r->sig[SIGSZ-1] = SIG_MSB;
2290 round_for_format (fmt, r)
2291 const struct real_format *fmt;
2295 unsigned long sticky;
2299 p2 = fmt->p * fmt->log2_b;
2300 emin2m1 = (fmt->emin - 1) * fmt->log2_b;
2301 emax2 = fmt->emax * fmt->log2_b;
2303 np2 = SIGNIFICAND_BITS - p2;
2307 get_zero (r, r->sign);
2309 if (!fmt->has_signed_zero)
2314 get_inf (r, r->sign);
2319 clear_significand_below (r, np2);
2329 /* If we're not base2, normalize the exponent to a multiple of
2331 if (fmt->log2_b != 1)
2333 int shift = r->exp & (fmt->log2_b - 1);
2336 shift = fmt->log2_b - shift;
2337 r->sig[0] |= sticky_rshift_significand (r, r, shift);
2342 /* Check the range of the exponent. If we're out of range,
2343 either underflow or overflow. */
2346 else if (r->exp <= emin2m1)
2350 if (!fmt->has_denorm)
2352 /* Don't underflow completely until we've had a chance to round. */
2353 if (r->exp < emin2m1)
2358 diff = emin2m1 - r->exp + 1;
2362 /* De-normalize the significand. */
2363 r->sig[0] |= sticky_rshift_significand (r, r, diff);
2368 /* There are P2 true significand bits, followed by one guard bit,
2369 followed by one sticky bit, followed by stuff. Fold nonzero
2370 stuff into the sticky bit. */
2373 for (i = 0, w = (np2 - 1) / HOST_BITS_PER_LONG; i < w; ++i)
2374 sticky |= r->sig[i];
2376 r->sig[w] & (((unsigned long)1 << ((np2 - 1) % HOST_BITS_PER_LONG)) - 1);
2378 guard = test_significand_bit (r, np2 - 1);
2379 lsb = test_significand_bit (r, np2);
2381 /* Round to even. */
2382 if (guard && (sticky || lsb))
2386 set_significand_bit (&u, np2);
2388 if (add_significands (r, r, &u))
2390 /* Overflow. Means the significand had been all ones, and
2391 is now all zeros. Need to increase the exponent, and
2392 possibly re-normalize it. */
2393 if (++r->exp > emax2)
2395 r->sig[SIGSZ-1] = SIG_MSB;
2397 if (fmt->log2_b != 1)
2399 int shift = r->exp & (fmt->log2_b - 1);
2402 shift = fmt->log2_b - shift;
2403 rshift_significand (r, r, shift);
2412 /* Catch underflow that we deferred until after rounding. */
2413 if (r->exp <= emin2m1)
2416 /* Clear out trailing garbage. */
2417 clear_significand_below (r, np2);
2420 /* Extend or truncate to a new mode. */
2423 real_convert (r, mode, a)
2425 enum machine_mode mode;
2426 const REAL_VALUE_TYPE *a;
2428 const struct real_format *fmt;
2430 fmt = real_format_for_mode[mode - QFmode];
2435 round_for_format (fmt, r);
2437 /* round_for_format de-normalizes denormals. Undo just that part. */
2438 if (r->class == rvc_normal)
2442 /* Legacy. Likewise, except return the struct directly. */
2445 real_value_truncate (mode, a)
2446 enum machine_mode mode;
2450 real_convert (&r, mode, &a);
2454 /* Return true if truncating to MODE is exact. */
2457 exact_real_truncate (mode, a)
2458 enum machine_mode mode;
2459 const REAL_VALUE_TYPE *a;
2462 real_convert (&t, mode, a);
2463 return real_identical (&t, a);
2466 /* Write R to the given target format. Place the words of the result
2467 in target word order in BUF. There are always 32 bits in each
2468 long, no matter the size of the host long.
2470 Legacy: return word 0 for implementing REAL_VALUE_TO_TARGET_SINGLE. */
2473 real_to_target_fmt (buf, r_orig, fmt)
2475 const REAL_VALUE_TYPE *r_orig;
2476 const struct real_format *fmt;
2482 round_for_format (fmt, &r);
2486 (*fmt->encode) (fmt, buf, &r);
2491 /* Similar, but look up the format from MODE. */
2494 real_to_target (buf, r, mode)
2496 const REAL_VALUE_TYPE *r;
2497 enum machine_mode mode;
2499 const struct real_format *fmt;
2501 fmt = real_format_for_mode[mode - QFmode];
2505 return real_to_target_fmt (buf, r, fmt);
2508 /* Read R from the given target format. Read the words of the result
2509 in target word order in BUF. There are always 32 bits in each
2510 long, no matter the size of the host long. */
2513 real_from_target_fmt (r, buf, fmt)
2516 const struct real_format *fmt;
2518 (*fmt->decode) (fmt, r, buf);
2521 /* Similar, but look up the format from MODE. */
2524 real_from_target (r, buf, mode)
2527 enum machine_mode mode;
2529 const struct real_format *fmt;
2531 fmt = real_format_for_mode[mode - QFmode];
2535 (*fmt->decode) (fmt, r, buf);
2538 /* Return the number of bits in the significand for MODE. */
2539 /* ??? Legacy. Should get access to real_format directly. */
2542 significand_size (mode)
2543 enum machine_mode mode;
2545 const struct real_format *fmt;
2547 fmt = real_format_for_mode[mode - QFmode];
2551 return fmt->p * fmt->log2_b;
2554 /* Return a hash value for the given real value. */
2555 /* ??? The "unsigned int" return value is intended to be hashval_t,
2556 but I didn't want to pull hashtab.h into real.h. */
2560 const REAL_VALUE_TYPE *r;
2565 h = r->class | (r->sign << 2);
2578 h ^= (unsigned int)-1;
2587 if (sizeof(unsigned long) > sizeof(unsigned int))
2588 for (i = 0; i < SIGSZ; ++i)
2590 unsigned long s = r->sig[i];
2591 h ^= s ^ (s >> (HOST_BITS_PER_LONG / 2));
2594 for (i = 0; i < SIGSZ; ++i)
2600 /* IEEE single-precision format. */
2602 static void encode_ieee_single PARAMS ((const struct real_format *fmt,
2603 long *, const REAL_VALUE_TYPE *));
2604 static void decode_ieee_single PARAMS ((const struct real_format *,
2605 REAL_VALUE_TYPE *, const long *));
2608 encode_ieee_single (fmt, buf, r)
2609 const struct real_format *fmt;
2611 const REAL_VALUE_TYPE *r;
2613 unsigned long image, sig, exp;
2614 bool denormal = (r->sig[SIGSZ-1] & SIG_MSB) == 0;
2616 image = r->sign << 31;
2617 sig = (r->sig[SIGSZ-1] >> (HOST_BITS_PER_LONG - 24)) & 0x7fffff;
2628 image |= 0x7fffffff;
2636 if (r->signalling == fmt->qnan_msb_set)
2640 /* We overload qnan_msb_set here: it's only clear for
2641 mips_ieee_single, which wants all mantissa bits but the
2642 quiet/signalling one set in canonical NaNs (at least
2644 if (r->canonical && !fmt->qnan_msb_set)
2645 sig |= (1 << 22) - 1;
2653 image |= 0x7fffffff;
2657 /* Recall that IEEE numbers are interpreted as 1.F x 2**exp,
2658 whereas the intermediate representation is 0.F x 2**exp.
2659 Which means we're off by one. */
2663 exp = r->exp + 127 - 1;
2676 decode_ieee_single (fmt, r, buf)
2677 const struct real_format *fmt;
2681 unsigned long image = buf[0] & 0xffffffff;
2682 bool sign = (image >> 31) & 1;
2683 int exp = (image >> 23) & 0xff;
2685 memset (r, 0, sizeof (*r));
2686 image <<= HOST_BITS_PER_LONG - 24;
2691 if (image && fmt->has_denorm)
2693 r->class = rvc_normal;
2696 r->sig[SIGSZ-1] = image << 1;
2699 else if (fmt->has_signed_zero)
2702 else if (exp == 255 && (fmt->has_nans || fmt->has_inf))
2708 r->signalling = ((image >> 22) & 1) ^ fmt->qnan_msb_set;
2709 r->sig[SIGSZ-1] = image;
2719 r->class = rvc_normal;
2721 r->exp = exp - 127 + 1;
2722 r->sig[SIGSZ-1] = image | SIG_MSB;
2726 const struct real_format ieee_single_format =
2744 const struct real_format mips_single_format =
2763 /* IEEE double-precision format. */
2765 static void encode_ieee_double PARAMS ((const struct real_format *fmt,
2766 long *, const REAL_VALUE_TYPE *));
2767 static void decode_ieee_double PARAMS ((const struct real_format *,
2768 REAL_VALUE_TYPE *, const long *));
2771 encode_ieee_double (fmt, buf, r)
2772 const struct real_format *fmt;
2774 const REAL_VALUE_TYPE *r;
2776 unsigned long image_lo, image_hi, sig_lo, sig_hi, exp;
2777 bool denormal = (r->sig[SIGSZ-1] & SIG_MSB) == 0;
2779 image_hi = r->sign << 31;
2782 if (HOST_BITS_PER_LONG == 64)
2784 sig_hi = r->sig[SIGSZ-1];
2785 sig_lo = (sig_hi >> (64 - 53)) & 0xffffffff;
2786 sig_hi = (sig_hi >> (64 - 53 + 1) >> 31) & 0xfffff;
2790 sig_hi = r->sig[SIGSZ-1];
2791 sig_lo = r->sig[SIGSZ-2];
2792 sig_lo = (sig_hi << 21) | (sig_lo >> 11);
2793 sig_hi = (sig_hi >> 11) & 0xfffff;
2803 image_hi |= 2047 << 20;
2806 image_hi |= 0x7fffffff;
2807 image_lo = 0xffffffff;
2815 sig_hi = sig_lo = 0;
2816 if (r->signalling == fmt->qnan_msb_set)
2817 sig_hi &= ~(1 << 19);
2820 /* We overload qnan_msb_set here: it's only clear for
2821 mips_ieee_single, which wants all mantissa bits but the
2822 quiet/signalling one set in canonical NaNs (at least
2824 if (r->canonical && !fmt->qnan_msb_set)
2826 sig_hi |= (1 << 19) - 1;
2827 sig_lo = 0xffffffff;
2829 else if (sig_hi == 0 && sig_lo == 0)
2832 image_hi |= 2047 << 20;
2838 image_hi |= 0x7fffffff;
2839 image_lo = 0xffffffff;
2844 /* Recall that IEEE numbers are interpreted as 1.F x 2**exp,
2845 whereas the intermediate representation is 0.F x 2**exp.
2846 Which means we're off by one. */
2850 exp = r->exp + 1023 - 1;
2851 image_hi |= exp << 20;
2860 if (FLOAT_WORDS_BIG_ENDIAN)
2861 buf[0] = image_hi, buf[1] = image_lo;
2863 buf[0] = image_lo, buf[1] = image_hi;
2867 decode_ieee_double (fmt, r, buf)
2868 const struct real_format *fmt;
2872 unsigned long image_hi, image_lo;
2876 if (FLOAT_WORDS_BIG_ENDIAN)
2877 image_hi = buf[0], image_lo = buf[1];
2879 image_lo = buf[0], image_hi = buf[1];
2880 image_lo &= 0xffffffff;
2881 image_hi &= 0xffffffff;
2883 sign = (image_hi >> 31) & 1;
2884 exp = (image_hi >> 20) & 0x7ff;
2886 memset (r, 0, sizeof (*r));
2888 image_hi <<= 32 - 21;
2889 image_hi |= image_lo >> 21;
2890 image_hi &= 0x7fffffff;
2891 image_lo <<= 32 - 21;
2895 if ((image_hi || image_lo) && fmt->has_denorm)
2897 r->class = rvc_normal;
2900 if (HOST_BITS_PER_LONG == 32)
2902 image_hi = (image_hi << 1) | (image_lo >> 31);
2904 r->sig[SIGSZ-1] = image_hi;
2905 r->sig[SIGSZ-2] = image_lo;
2909 image_hi = (image_hi << 31 << 2) | (image_lo << 1);
2910 r->sig[SIGSZ-1] = image_hi;
2914 else if (fmt->has_signed_zero)
2917 else if (exp == 2047 && (fmt->has_nans || fmt->has_inf))
2919 if (image_hi || image_lo)
2923 r->signalling = ((image_hi >> 30) & 1) ^ fmt->qnan_msb_set;
2924 if (HOST_BITS_PER_LONG == 32)
2926 r->sig[SIGSZ-1] = image_hi;
2927 r->sig[SIGSZ-2] = image_lo;
2930 r->sig[SIGSZ-1] = (image_hi << 31 << 1) | image_lo;
2940 r->class = rvc_normal;
2942 r->exp = exp - 1023 + 1;
2943 if (HOST_BITS_PER_LONG == 32)
2945 r->sig[SIGSZ-1] = image_hi | SIG_MSB;
2946 r->sig[SIGSZ-2] = image_lo;
2949 r->sig[SIGSZ-1] = (image_hi << 31 << 1) | image_lo | SIG_MSB;
2953 const struct real_format ieee_double_format =
2971 const struct real_format mips_double_format =
2990 /* IEEE extended double precision format. This comes in three
2991 flavours: Intel's as a 12 byte image, Intel's as a 16 byte image,
2994 static void encode_ieee_extended PARAMS ((const struct real_format *fmt,
2995 long *, const REAL_VALUE_TYPE *));
2996 static void decode_ieee_extended PARAMS ((const struct real_format *,
2997 REAL_VALUE_TYPE *, const long *));
2999 static void encode_ieee_extended_128 PARAMS ((const struct real_format *fmt,
3001 const REAL_VALUE_TYPE *));
3002 static void decode_ieee_extended_128 PARAMS ((const struct real_format *,
3007 encode_ieee_extended (fmt, buf, r)
3008 const struct real_format *fmt;
3010 const REAL_VALUE_TYPE *r;
3012 unsigned long image_hi, sig_hi, sig_lo;
3013 bool denormal = (r->sig[SIGSZ-1] & SIG_MSB) == 0;
3015 image_hi = r->sign << 15;
3016 sig_hi = sig_lo = 0;
3028 /* Intel requires the explicit integer bit to be set, otherwise
3029 it considers the value a "pseudo-infinity". Motorola docs
3030 say it doesn't care. */
3031 sig_hi = 0x80000000;
3036 sig_lo = sig_hi = 0xffffffff;
3044 if (HOST_BITS_PER_LONG == 32)
3046 sig_hi = r->sig[SIGSZ-1];
3047 sig_lo = r->sig[SIGSZ-2];
3051 sig_lo = r->sig[SIGSZ-1];
3052 sig_hi = sig_lo >> 31 >> 1;
3053 sig_lo &= 0xffffffff;
3055 if (r->signalling == fmt->qnan_msb_set)
3056 sig_hi &= ~(1 << 30);
3059 if ((sig_hi & 0x7fffffff) == 0 && sig_lo == 0)
3062 /* Intel requires the explicit integer bit to be set, otherwise
3063 it considers the value a "pseudo-nan". Motorola docs say it
3065 sig_hi |= 0x80000000;
3070 sig_lo = sig_hi = 0xffffffff;
3078 /* Recall that IEEE numbers are interpreted as 1.F x 2**exp,
3079 whereas the intermediate representation is 0.F x 2**exp.
3080 Which means we're off by one.
3082 Except for Motorola, which consider exp=0 and explicit
3083 integer bit set to continue to be normalized. In theory
3084 this discrepancy has been taken care of by the difference
3085 in fmt->emin in round_for_format. */
3097 if (HOST_BITS_PER_LONG == 32)
3099 sig_hi = r->sig[SIGSZ-1];
3100 sig_lo = r->sig[SIGSZ-2];
3104 sig_lo = r->sig[SIGSZ-1];
3105 sig_hi = sig_lo >> 31 >> 1;
3106 sig_lo &= 0xffffffff;
3115 if (FLOAT_WORDS_BIG_ENDIAN)
3116 buf[0] = image_hi << 16, buf[1] = sig_hi, buf[2] = sig_lo;
3118 buf[0] = sig_lo, buf[1] = sig_hi, buf[2] = image_hi;
3122 encode_ieee_extended_128 (fmt, buf, r)
3123 const struct real_format *fmt;
3125 const REAL_VALUE_TYPE *r;
3127 buf[3 * !FLOAT_WORDS_BIG_ENDIAN] = 0;
3128 encode_ieee_extended (fmt, buf+!!FLOAT_WORDS_BIG_ENDIAN, r);
3132 decode_ieee_extended (fmt, r, buf)
3133 const struct real_format *fmt;
3137 unsigned long image_hi, sig_hi, sig_lo;
3141 if (FLOAT_WORDS_BIG_ENDIAN)
3142 image_hi = buf[0] >> 16, sig_hi = buf[1], sig_lo = buf[2];
3144 sig_lo = buf[0], sig_hi = buf[1], image_hi = buf[2];
3145 sig_lo &= 0xffffffff;
3146 sig_hi &= 0xffffffff;
3147 image_hi &= 0xffffffff;
3149 sign = (image_hi >> 15) & 1;
3150 exp = image_hi & 0x7fff;
3152 memset (r, 0, sizeof (*r));
3156 if ((sig_hi || sig_lo) && fmt->has_denorm)
3158 r->class = rvc_normal;
3161 /* When the IEEE format contains a hidden bit, we know that
3162 it's zero at this point, and so shift up the significand
3163 and decrease the exponent to match. In this case, Motorola
3164 defines the explicit integer bit to be valid, so we don't
3165 know whether the msb is set or not. */
3167 if (HOST_BITS_PER_LONG == 32)
3169 r->sig[SIGSZ-1] = sig_hi;
3170 r->sig[SIGSZ-2] = sig_lo;
3173 r->sig[SIGSZ-1] = (sig_hi << 31 << 1) | sig_lo;
3177 else if (fmt->has_signed_zero)
3180 else if (exp == 32767 && (fmt->has_nans || fmt->has_inf))
3182 /* See above re "pseudo-infinities" and "pseudo-nans".
3183 Short summary is that the MSB will likely always be
3184 set, and that we don't care about it. */
3185 sig_hi &= 0x7fffffff;
3187 if (sig_hi || sig_lo)
3191 r->signalling = ((sig_hi >> 30) & 1) ^ fmt->qnan_msb_set;
3192 if (HOST_BITS_PER_LONG == 32)
3194 r->sig[SIGSZ-1] = sig_hi;
3195 r->sig[SIGSZ-2] = sig_lo;
3198 r->sig[SIGSZ-1] = (sig_hi << 31 << 1) | sig_lo;
3208 r->class = rvc_normal;
3210 r->exp = exp - 16383 + 1;
3211 if (HOST_BITS_PER_LONG == 32)
3213 r->sig[SIGSZ-1] = sig_hi;
3214 r->sig[SIGSZ-2] = sig_lo;
3217 r->sig[SIGSZ-1] = (sig_hi << 31 << 1) | sig_lo;
3222 decode_ieee_extended_128 (fmt, r, buf)
3223 const struct real_format *fmt;
3227 decode_ieee_extended (fmt, r, buf+!!FLOAT_WORDS_BIG_ENDIAN);
3230 const struct real_format ieee_extended_motorola_format =
3232 encode_ieee_extended,
3233 decode_ieee_extended,
3248 const struct real_format ieee_extended_intel_96_format =
3250 encode_ieee_extended,
3251 decode_ieee_extended,
3266 const struct real_format ieee_extended_intel_128_format =
3268 encode_ieee_extended_128,
3269 decode_ieee_extended_128,
3285 /* IBM 128-bit extended precision format: a pair of IEEE double precision
3286 numbers whose sum is equal to the extended precision value. The number
3287 with greater magnitude is first. This format has the same magnitude
3288 range as an IEEE double precision value, but effectively 106 bits of
3289 significand precision. Infinity and NaN are represented by their IEEE
3290 double precision value stored in the first number, the second number is
3291 ignored. Zeroes, Infinities, and NaNs are set in both doubles
3292 due to precedent. */
3294 static void encode_ibm_extended PARAMS ((const struct real_format *fmt,
3295 long *, const REAL_VALUE_TYPE *));
3296 static void decode_ibm_extended PARAMS ((const struct real_format *,
3297 REAL_VALUE_TYPE *, const long *));
3300 encode_ibm_extended (fmt, buf, r)
3301 const struct real_format *fmt;
3303 const REAL_VALUE_TYPE *r;
3305 REAL_VALUE_TYPE u, v;
3306 const struct real_format *base_fmt;
3308 base_fmt = fmt->qnan_msb_set ? &ieee_double_format : &mips_double_format;
3313 /* Both doubles have sign bit set. */
3314 buf[0] = FLOAT_WORDS_BIG_ENDIAN ? r->sign << 31 : 0;
3315 buf[1] = FLOAT_WORDS_BIG_ENDIAN ? 0 : r->sign << 31;
3322 /* Both doubles set to Inf / NaN. */
3323 encode_ieee_double (base_fmt, &buf[0], r);
3329 /* u = IEEE double precision portion of significand. */
3331 clear_significand_below (&u, SIGNIFICAND_BITS - 53);
3334 /* If the upper double is zero, we have a denormal double, so
3335 move it to the first double and leave the second as zero. */
3336 if (u.class == rvc_zero)
3344 /* v = remainder containing additional 53 bits of significand. */
3345 do_add (&v, r, &u, 1);
3346 round_for_format (base_fmt, &v);
3349 round_for_format (base_fmt, &u);
3351 encode_ieee_double (base_fmt, &buf[0], &u);
3352 encode_ieee_double (base_fmt, &buf[2], &v);
3361 decode_ibm_extended (fmt, r, buf)
3362 const struct real_format *fmt ATTRIBUTE_UNUSED;
3366 REAL_VALUE_TYPE u, v;
3367 const struct real_format *base_fmt;
3369 base_fmt = fmt->qnan_msb_set ? &ieee_double_format : &mips_double_format;
3370 decode_ieee_double (base_fmt, &u, &buf[0]);
3372 if (u.class != rvc_zero && u.class != rvc_inf && u.class != rvc_nan)
3374 decode_ieee_double (base_fmt, &v, &buf[2]);
3375 do_add (r, &u, &v, 0);
3381 const struct real_format ibm_extended_format =
3383 encode_ibm_extended,
3384 decode_ibm_extended,
3399 const struct real_format mips_extended_format =
3401 encode_ibm_extended,
3402 decode_ibm_extended,
3418 /* IEEE quad precision format. */
3420 static void encode_ieee_quad PARAMS ((const struct real_format *fmt,
3421 long *, const REAL_VALUE_TYPE *));
3422 static void decode_ieee_quad PARAMS ((const struct real_format *,
3423 REAL_VALUE_TYPE *, const long *));
3426 encode_ieee_quad (fmt, buf, r)
3427 const struct real_format *fmt;
3429 const REAL_VALUE_TYPE *r;
3431 unsigned long image3, image2, image1, image0, exp;
3432 bool denormal = (r->sig[SIGSZ-1] & SIG_MSB) == 0;
3435 image3 = r->sign << 31;
3440 rshift_significand (&u, r, SIGNIFICAND_BITS - 113);
3449 image3 |= 32767 << 16;
3452 image3 |= 0x7fffffff;
3453 image2 = 0xffffffff;
3454 image1 = 0xffffffff;
3455 image0 = 0xffffffff;
3462 image3 |= 32767 << 16;
3466 /* Don't use bits from the significand. The
3467 initialization above is right. */
3469 else if (HOST_BITS_PER_LONG == 32)
3474 image3 |= u.sig[3] & 0xffff;
3479 image1 = image0 >> 31 >> 1;
3481 image3 |= (image2 >> 31 >> 1) & 0xffff;
3482 image0 &= 0xffffffff;
3483 image2 &= 0xffffffff;
3485 if (r->signalling == fmt->qnan_msb_set)
3489 /* We overload qnan_msb_set here: it's only clear for
3490 mips_ieee_single, which wants all mantissa bits but the
3491 quiet/signalling one set in canonical NaNs (at least
3493 if (r->canonical && !fmt->qnan_msb_set)
3496 image2 = image1 = image0 = 0xffffffff;
3498 else if (((image3 & 0xffff) | image2 | image1 | image0) == 0)
3503 image3 |= 0x7fffffff;
3504 image2 = 0xffffffff;
3505 image1 = 0xffffffff;
3506 image0 = 0xffffffff;
3511 /* Recall that IEEE numbers are interpreted as 1.F x 2**exp,
3512 whereas the intermediate representation is 0.F x 2**exp.
3513 Which means we're off by one. */
3517 exp = r->exp + 16383 - 1;
3518 image3 |= exp << 16;
3520 if (HOST_BITS_PER_LONG == 32)
3525 image3 |= u.sig[3] & 0xffff;
3530 image1 = image0 >> 31 >> 1;
3532 image3 |= (image2 >> 31 >> 1) & 0xffff;
3533 image0 &= 0xffffffff;
3534 image2 &= 0xffffffff;
3542 if (FLOAT_WORDS_BIG_ENDIAN)
3559 decode_ieee_quad (fmt, r, buf)
3560 const struct real_format *fmt;
3564 unsigned long image3, image2, image1, image0;
3568 if (FLOAT_WORDS_BIG_ENDIAN)
3582 image0 &= 0xffffffff;
3583 image1 &= 0xffffffff;
3584 image2 &= 0xffffffff;
3586 sign = (image3 >> 31) & 1;
3587 exp = (image3 >> 16) & 0x7fff;
3590 memset (r, 0, sizeof (*r));
3594 if ((image3 | image2 | image1 | image0) && fmt->has_denorm)
3596 r->class = rvc_normal;
3599 r->exp = -16382 + (SIGNIFICAND_BITS - 112);
3600 if (HOST_BITS_PER_LONG == 32)
3609 r->sig[0] = (image1 << 31 << 1) | image0;
3610 r->sig[1] = (image3 << 31 << 1) | image2;
3615 else if (fmt->has_signed_zero)
3618 else if (exp == 32767 && (fmt->has_nans || fmt->has_inf))
3620 if (image3 | image2 | image1 | image0)
3624 r->signalling = ((image3 >> 15) & 1) ^ fmt->qnan_msb_set;
3626 if (HOST_BITS_PER_LONG == 32)
3635 r->sig[0] = (image1 << 31 << 1) | image0;
3636 r->sig[1] = (image3 << 31 << 1) | image2;
3638 lshift_significand (r, r, SIGNIFICAND_BITS - 113);
3648 r->class = rvc_normal;
3650 r->exp = exp - 16383 + 1;
3652 if (HOST_BITS_PER_LONG == 32)
3661 r->sig[0] = (image1 << 31 << 1) | image0;
3662 r->sig[1] = (image3 << 31 << 1) | image2;
3664 lshift_significand (r, r, SIGNIFICAND_BITS - 113);
3665 r->sig[SIGSZ-1] |= SIG_MSB;
3669 const struct real_format ieee_quad_format =
3687 const struct real_format mips_quad_format =
3705 /* Descriptions of VAX floating point formats can be found beginning at
3707 http://www.openvms.compaq.com:8000/73final/4515/4515pro_013.html#f_floating_point_format
3709 The thing to remember is that they're almost IEEE, except for word
3710 order, exponent bias, and the lack of infinities, nans, and denormals.
3712 We don't implement the H_floating format here, simply because neither
3713 the VAX or Alpha ports use it. */
3715 static void encode_vax_f PARAMS ((const struct real_format *fmt,
3716 long *, const REAL_VALUE_TYPE *));
3717 static void decode_vax_f PARAMS ((const struct real_format *,
3718 REAL_VALUE_TYPE *, const long *));
3719 static void encode_vax_d PARAMS ((const struct real_format *fmt,
3720 long *, const REAL_VALUE_TYPE *));
3721 static void decode_vax_d PARAMS ((const struct real_format *,
3722 REAL_VALUE_TYPE *, const long *));
3723 static void encode_vax_g PARAMS ((const struct real_format *fmt,
3724 long *, const REAL_VALUE_TYPE *));
3725 static void decode_vax_g PARAMS ((const struct real_format *,
3726 REAL_VALUE_TYPE *, const long *));
3729 encode_vax_f (fmt, buf, r)
3730 const struct real_format *fmt ATTRIBUTE_UNUSED;
3732 const REAL_VALUE_TYPE *r;
3734 unsigned long sign, exp, sig, image;
3736 sign = r->sign << 15;
3746 image = 0xffff7fff | sign;
3750 sig = (r->sig[SIGSZ-1] >> (HOST_BITS_PER_LONG - 24)) & 0x7fffff;
3753 image = (sig << 16) & 0xffff0000;
3767 decode_vax_f (fmt, r, buf)
3768 const struct real_format *fmt ATTRIBUTE_UNUSED;
3772 unsigned long image = buf[0] & 0xffffffff;
3773 int exp = (image >> 7) & 0xff;
3775 memset (r, 0, sizeof (*r));
3779 r->class = rvc_normal;
3780 r->sign = (image >> 15) & 1;
3783 image = ((image & 0x7f) << 16) | ((image >> 16) & 0xffff);
3784 r->sig[SIGSZ-1] = (image << (HOST_BITS_PER_LONG - 24)) | SIG_MSB;
3789 encode_vax_d (fmt, buf, r)
3790 const struct real_format *fmt ATTRIBUTE_UNUSED;
3792 const REAL_VALUE_TYPE *r;
3794 unsigned long image0, image1, sign = r->sign << 15;
3799 image0 = image1 = 0;
3804 image0 = 0xffff7fff | sign;
3805 image1 = 0xffffffff;
3809 /* Extract the significand into straight hi:lo. */
3810 if (HOST_BITS_PER_LONG == 64)
3812 image0 = r->sig[SIGSZ-1];
3813 image1 = (image0 >> (64 - 56)) & 0xffffffff;
3814 image0 = (image0 >> (64 - 56 + 1) >> 31) & 0x7fffff;
3818 image0 = r->sig[SIGSZ-1];
3819 image1 = r->sig[SIGSZ-2];
3820 image1 = (image0 << 24) | (image1 >> 8);
3821 image0 = (image0 >> 8) & 0xffffff;
3824 /* Rearrange the half-words of the significand to match the
3826 image0 = ((image0 << 16) | (image0 >> 16)) & 0xffff007f;
3827 image1 = ((image1 << 16) | (image1 >> 16)) & 0xffffffff;
3829 /* Add the sign and exponent. */
3831 image0 |= (r->exp + 128) << 7;
3838 if (FLOAT_WORDS_BIG_ENDIAN)
3839 buf[0] = image1, buf[1] = image0;
3841 buf[0] = image0, buf[1] = image1;
3845 decode_vax_d (fmt, r, buf)
3846 const struct real_format *fmt ATTRIBUTE_UNUSED;
3850 unsigned long image0, image1;
3853 if (FLOAT_WORDS_BIG_ENDIAN)
3854 image1 = buf[0], image0 = buf[1];
3856 image0 = buf[0], image1 = buf[1];
3857 image0 &= 0xffffffff;
3858 image1 &= 0xffffffff;
3860 exp = (image0 >> 7) & 0x7f;
3862 memset (r, 0, sizeof (*r));
3866 r->class = rvc_normal;
3867 r->sign = (image0 >> 15) & 1;
3870 /* Rearrange the half-words of the external format into
3871 proper ascending order. */
3872 image0 = ((image0 & 0x7f) << 16) | ((image0 >> 16) & 0xffff);
3873 image1 = ((image1 & 0xffff) << 16) | ((image1 >> 16) & 0xffff);
3875 if (HOST_BITS_PER_LONG == 64)
3877 image0 = (image0 << 31 << 1) | image1;
3880 r->sig[SIGSZ-1] = image0;
3884 r->sig[SIGSZ-1] = image0;
3885 r->sig[SIGSZ-2] = image1;
3886 lshift_significand (r, r, 2*HOST_BITS_PER_LONG - 56);
3887 r->sig[SIGSZ-1] |= SIG_MSB;
3893 encode_vax_g (fmt, buf, r)
3894 const struct real_format *fmt ATTRIBUTE_UNUSED;
3896 const REAL_VALUE_TYPE *r;
3898 unsigned long image0, image1, sign = r->sign << 15;
3903 image0 = image1 = 0;
3908 image0 = 0xffff7fff | sign;
3909 image1 = 0xffffffff;
3913 /* Extract the significand into straight hi:lo. */
3914 if (HOST_BITS_PER_LONG == 64)
3916 image0 = r->sig[SIGSZ-1];
3917 image1 = (image0 >> (64 - 53)) & 0xffffffff;
3918 image0 = (image0 >> (64 - 53 + 1) >> 31) & 0xfffff;
3922 image0 = r->sig[SIGSZ-1];
3923 image1 = r->sig[SIGSZ-2];
3924 image1 = (image0 << 21) | (image1 >> 11);
3925 image0 = (image0 >> 11) & 0xfffff;
3928 /* Rearrange the half-words of the significand to match the
3930 image0 = ((image0 << 16) | (image0 >> 16)) & 0xffff000f;
3931 image1 = ((image1 << 16) | (image1 >> 16)) & 0xffffffff;
3933 /* Add the sign and exponent. */
3935 image0 |= (r->exp + 1024) << 4;
3942 if (FLOAT_WORDS_BIG_ENDIAN)
3943 buf[0] = image1, buf[1] = image0;
3945 buf[0] = image0, buf[1] = image1;
3949 decode_vax_g (fmt, r, buf)
3950 const struct real_format *fmt ATTRIBUTE_UNUSED;
3954 unsigned long image0, image1;
3957 if (FLOAT_WORDS_BIG_ENDIAN)
3958 image1 = buf[0], image0 = buf[1];
3960 image0 = buf[0], image1 = buf[1];
3961 image0 &= 0xffffffff;
3962 image1 &= 0xffffffff;
3964 exp = (image0 >> 4) & 0x7ff;
3966 memset (r, 0, sizeof (*r));
3970 r->class = rvc_normal;
3971 r->sign = (image0 >> 15) & 1;
3972 r->exp = exp - 1024;
3974 /* Rearrange the half-words of the external format into
3975 proper ascending order. */
3976 image0 = ((image0 & 0xf) << 16) | ((image0 >> 16) & 0xffff);
3977 image1 = ((image1 & 0xffff) << 16) | ((image1 >> 16) & 0xffff);
3979 if (HOST_BITS_PER_LONG == 64)
3981 image0 = (image0 << 31 << 1) | image1;
3984 r->sig[SIGSZ-1] = image0;
3988 r->sig[SIGSZ-1] = image0;
3989 r->sig[SIGSZ-2] = image1;
3990 lshift_significand (r, r, 64 - 53);
3991 r->sig[SIGSZ-1] |= SIG_MSB;
3996 const struct real_format vax_f_format =
4014 const struct real_format vax_d_format =
4032 const struct real_format vax_g_format =
4050 /* A good reference for these can be found in chapter 9 of
4051 "ESA/390 Principles of Operation", IBM document number SA22-7201-01.
4052 An on-line version can be found here:
4054 http://publibz.boulder.ibm.com/cgi-bin/bookmgr_OS390/BOOKS/DZ9AR001/9.1?DT=19930923083613
4057 static void encode_i370_single PARAMS ((const struct real_format *fmt,
4058 long *, const REAL_VALUE_TYPE *));
4059 static void decode_i370_single PARAMS ((const struct real_format *,
4060 REAL_VALUE_TYPE *, const long *));
4061 static void encode_i370_double PARAMS ((const struct real_format *fmt,
4062 long *, const REAL_VALUE_TYPE *));
4063 static void decode_i370_double PARAMS ((const struct real_format *,
4064 REAL_VALUE_TYPE *, const long *));
4067 encode_i370_single (fmt, buf, r)
4068 const struct real_format *fmt ATTRIBUTE_UNUSED;
4070 const REAL_VALUE_TYPE *r;
4072 unsigned long sign, exp, sig, image;
4074 sign = r->sign << 31;
4084 image = 0x7fffffff | sign;
4088 sig = (r->sig[SIGSZ-1] >> (HOST_BITS_PER_LONG - 24)) & 0xffffff;
4089 exp = ((r->exp / 4) + 64) << 24;
4090 image = sign | exp | sig;
4101 decode_i370_single (fmt, r, buf)
4102 const struct real_format *fmt ATTRIBUTE_UNUSED;
4106 unsigned long sign, sig, image = buf[0];
4109 sign = (image >> 31) & 1;
4110 exp = (image >> 24) & 0x7f;
4111 sig = image & 0xffffff;
4113 memset (r, 0, sizeof (*r));
4117 r->class = rvc_normal;
4119 r->exp = (exp - 64) * 4;
4120 r->sig[SIGSZ-1] = sig << (HOST_BITS_PER_LONG - 24);
4126 encode_i370_double (fmt, buf, r)
4127 const struct real_format *fmt ATTRIBUTE_UNUSED;
4129 const REAL_VALUE_TYPE *r;
4131 unsigned long sign, exp, image_hi, image_lo;
4133 sign = r->sign << 31;
4138 image_hi = image_lo = 0;
4143 image_hi = 0x7fffffff | sign;
4144 image_lo = 0xffffffff;
4148 if (HOST_BITS_PER_LONG == 64)
4150 image_hi = r->sig[SIGSZ-1];
4151 image_lo = (image_hi >> (64 - 56)) & 0xffffffff;
4152 image_hi = (image_hi >> (64 - 56 + 1) >> 31) & 0xffffff;
4156 image_hi = r->sig[SIGSZ-1];
4157 image_lo = r->sig[SIGSZ-2];
4158 image_lo = (image_lo >> 8) | (image_hi << 24);
4162 exp = ((r->exp / 4) + 64) << 24;
4163 image_hi |= sign | exp;
4170 if (FLOAT_WORDS_BIG_ENDIAN)
4171 buf[0] = image_hi, buf[1] = image_lo;
4173 buf[0] = image_lo, buf[1] = image_hi;
4177 decode_i370_double (fmt, r, buf)
4178 const struct real_format *fmt ATTRIBUTE_UNUSED;
4182 unsigned long sign, image_hi, image_lo;
4185 if (FLOAT_WORDS_BIG_ENDIAN)
4186 image_hi = buf[0], image_lo = buf[1];
4188 image_lo = buf[0], image_hi = buf[1];
4190 sign = (image_hi >> 31) & 1;
4191 exp = (image_hi >> 24) & 0x7f;
4192 image_hi &= 0xffffff;
4193 image_lo &= 0xffffffff;
4195 memset (r, 0, sizeof (*r));
4197 if (exp || image_hi || image_lo)
4199 r->class = rvc_normal;
4201 r->exp = (exp - 64) * 4 + (SIGNIFICAND_BITS - 56);
4203 if (HOST_BITS_PER_LONG == 32)
4205 r->sig[0] = image_lo;
4206 r->sig[1] = image_hi;
4209 r->sig[0] = image_lo | (image_hi << 31 << 1);
4215 const struct real_format i370_single_format =
4228 false, /* ??? The encoding does allow for "unnormals". */
4229 false, /* ??? The encoding does allow for "unnormals". */
4233 const struct real_format i370_double_format =
4246 false, /* ??? The encoding does allow for "unnormals". */
4247 false, /* ??? The encoding does allow for "unnormals". */
4251 /* The "twos-complement" c4x format is officially defined as
4255 This is rather misleading. One must remember that F is signed.
4256 A better description would be
4258 x = -1**s * ((s + 1 + .f) * 2**e
4260 So if we have a (4 bit) fraction of .1000 with a sign bit of 1,
4261 that's -1 * (1+1+(-.5)) == -1.5. I think.
4263 The constructions here are taken from Tables 5-1 and 5-2 of the
4264 TMS320C4x User's Guide wherein step-by-step instructions for
4265 conversion from IEEE are presented. That's close enough to our
4266 internal representation so as to make things easy.
4268 See http://www-s.ti.com/sc/psheets/spru063c/spru063c.pdf */
4270 static void encode_c4x_single PARAMS ((const struct real_format *fmt,
4271 long *, const REAL_VALUE_TYPE *));
4272 static void decode_c4x_single PARAMS ((const struct real_format *,
4273 REAL_VALUE_TYPE *, const long *));
4274 static void encode_c4x_extended PARAMS ((const struct real_format *fmt,
4275 long *, const REAL_VALUE_TYPE *));
4276 static void decode_c4x_extended PARAMS ((const struct real_format *,
4277 REAL_VALUE_TYPE *, const long *));
4280 encode_c4x_single (fmt, buf, r)
4281 const struct real_format *fmt ATTRIBUTE_UNUSED;
4283 const REAL_VALUE_TYPE *r;
4285 unsigned long image, exp, sig;
4297 sig = 0x800000 - r->sign;
4302 sig = (r->sig[SIGSZ-1] >> (HOST_BITS_PER_LONG - 24)) & 0x7fffff;
4317 image = ((exp & 0xff) << 24) | (sig & 0xffffff);
4322 decode_c4x_single (fmt, r, buf)
4323 const struct real_format *fmt ATTRIBUTE_UNUSED;
4327 unsigned long image = buf[0];
4331 exp = (((image >> 24) & 0xff) ^ 0x80) - 0x80;
4332 sf = ((image & 0xffffff) ^ 0x800000) - 0x800000;
4334 memset (r, 0, sizeof (*r));
4338 r->class = rvc_normal;
4340 sig = sf & 0x7fffff;
4349 sig = (sig << (HOST_BITS_PER_LONG - 24)) | SIG_MSB;
4352 r->sig[SIGSZ-1] = sig;
4357 encode_c4x_extended (fmt, buf, r)
4358 const struct real_format *fmt ATTRIBUTE_UNUSED;
4360 const REAL_VALUE_TYPE *r;
4362 unsigned long exp, sig;
4374 sig = 0x80000000 - r->sign;
4380 sig = r->sig[SIGSZ-1];
4381 if (HOST_BITS_PER_LONG == 64)
4382 sig = sig >> 1 >> 31;
4399 exp = (exp & 0xff) << 24;
4402 if (FLOAT_WORDS_BIG_ENDIAN)
4403 buf[0] = exp, buf[1] = sig;
4405 buf[0] = sig, buf[0] = exp;
4409 decode_c4x_extended (fmt, r, buf)
4410 const struct real_format *fmt ATTRIBUTE_UNUSED;
4417 if (FLOAT_WORDS_BIG_ENDIAN)
4418 exp = buf[0], sf = buf[1];
4420 sf = buf[0], exp = buf[1];
4422 exp = (((exp >> 24) & 0xff) & 0x80) - 0x80;
4423 sf = ((sf & 0xffffffff) ^ 0x80000000) - 0x80000000;
4425 memset (r, 0, sizeof (*r));
4429 r->class = rvc_normal;
4431 sig = sf & 0x7fffffff;
4440 if (HOST_BITS_PER_LONG == 64)
4441 sig = sig << 1 << 31;
4445 r->sig[SIGSZ-1] = sig;
4449 const struct real_format c4x_single_format =
4467 const struct real_format c4x_extended_format =
4469 encode_c4x_extended,
4470 decode_c4x_extended,
4486 /* A synthetic "format" for internal arithmetic. It's the size of the
4487 internal significand minus the two bits needed for proper rounding.
4488 The encode and decode routines exist only to satisfy our paranoia
4491 static void encode_internal PARAMS ((const struct real_format *fmt,
4492 long *, const REAL_VALUE_TYPE *));
4493 static void decode_internal PARAMS ((const struct real_format *,
4494 REAL_VALUE_TYPE *, const long *));
4497 encode_internal (fmt, buf, r)
4498 const struct real_format *fmt ATTRIBUTE_UNUSED;
4500 const REAL_VALUE_TYPE *r;
4502 memcpy (buf, r, sizeof (*r));
4506 decode_internal (fmt, r, buf)
4507 const struct real_format *fmt ATTRIBUTE_UNUSED;
4511 memcpy (r, buf, sizeof (*r));
4514 const struct real_format real_internal_format =
4520 SIGNIFICAND_BITS - 2,
4521 SIGNIFICAND_BITS - 2,
4532 /* Set up default mode to format mapping for IEEE. Everyone else has
4533 to set these values in OVERRIDE_OPTIONS. */
4535 const struct real_format *real_format_for_mode[TFmode - QFmode + 1] =
4540 &ieee_single_format, /* SFmode */
4541 &ieee_double_format, /* DFmode */
4543 /* We explicitly don't handle XFmode. There are two formats,
4544 pretty much equally common. Choose one in OVERRIDE_OPTIONS. */
4546 &ieee_quad_format /* TFmode */
4550 /* Calculate the square root of X in mode MODE, and store the result
4551 in R. Return TRUE if the operation does not raise an exception.
4552 For details see "High Precision Division and Square Root",
4553 Alan H. Karp and Peter Markstein, HP Lab Report 93-93-42, June
4554 1993. http://www.hpl.hp.com/techreports/93/HPL-93-42.pdf. */
4557 real_sqrt (r, mode, x)
4559 enum machine_mode mode;
4560 const REAL_VALUE_TYPE *x;
4562 static REAL_VALUE_TYPE halfthree;
4563 static bool init = false;
4564 REAL_VALUE_TYPE h, t, i;
4567 /* sqrt(-0.0) is -0.0. */
4568 if (real_isnegzero (x))
4574 /* Negative arguments return NaN. */
4577 /* Mode is ignored for canonical NaN. */
4578 real_nan (r, "", 1, SFmode);
4582 /* Infinity and NaN return themselves. */
4583 if (real_isinf (x) || real_isnan (x))
4591 real_arithmetic (&halfthree, PLUS_EXPR, &dconst1, &dconsthalf);
4595 /* Initial guess for reciprocal sqrt, i. */
4596 exp = real_exponent (x);
4597 real_ldexp (&i, &dconst1, -exp/2);
4599 /* Newton's iteration for reciprocal sqrt, i. */
4600 for (iter = 0; iter < 16; iter++)
4602 /* i(n+1) = i(n) * (1.5 - 0.5*i(n)*i(n)*x). */
4603 real_arithmetic (&t, MULT_EXPR, x, &i);
4604 real_arithmetic (&h, MULT_EXPR, &t, &i);
4605 real_arithmetic (&t, MULT_EXPR, &h, &dconsthalf);
4606 real_arithmetic (&h, MINUS_EXPR, &halfthree, &t);
4607 real_arithmetic (&t, MULT_EXPR, &i, &h);
4609 /* Check for early convergence. */
4610 if (iter >= 6 && real_identical (&i, &t))
4613 /* ??? Unroll loop to avoid copying. */
4617 /* Final iteration: r = i*x + 0.5*i*x*(1.0 - i*(i*x)). */
4618 real_arithmetic (&t, MULT_EXPR, x, &i);
4619 real_arithmetic (&h, MULT_EXPR, &t, &i);
4620 real_arithmetic (&i, MINUS_EXPR, &dconst1, &h);
4621 real_arithmetic (&h, MULT_EXPR, &t, &i);
4622 real_arithmetic (&i, MULT_EXPR, &dconsthalf, &h);
4623 real_arithmetic (&h, PLUS_EXPR, &t, &i);
4625 /* ??? We need a Tuckerman test to get the last bit. */
4627 real_convert (r, mode, &h);