OSDN Git Service

* MAINTAINERS (c4x port): Remove.
[pf3gnuchains/gcc-fork.git] / gcc / real.c
1 /* real.c - software floating point emulation.
2    Copyright (C) 1993, 1994, 1995, 1996, 1997, 1998, 1999,
3    2000, 2002, 2003, 2004, 2005, 2007, 2008 Free Software Foundation, Inc.
4    Contributed by Stephen L. Moshier (moshier@world.std.com).
5    Re-written by Richard Henderson <rth@redhat.com>
6
7    This file is part of GCC.
8
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 3, or (at your option) any later
12    version.
13
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
17    for more details.
18
19    You should have received a copy of the GNU General Public License
20    along with GCC; see the file COPYING3.  If not see
21    <http://www.gnu.org/licenses/>.  */
22
23 #include "config.h"
24 #include "system.h"
25 #include "coretypes.h"
26 #include "tm.h"
27 #include "tree.h"
28 #include "toplev.h"
29 #include "real.h"
30 #include "tm_p.h"
31 #include "dfp.h"
32
33 /* The floating point model used internally is not exactly IEEE 754
34    compliant, and close to the description in the ISO C99 standard,
35    section 5.2.4.2.2 Characteristics of floating types.
36
37    Specifically
38
39         x = s * b^e * \sum_{k=1}^p f_k * b^{-k}
40
41         where
42                 s = sign (+- 1)
43                 b = base or radix, here always 2
44                 e = exponent
45                 p = precision (the number of base-b digits in the significand)
46                 f_k = the digits of the significand.
47
48    We differ from typical IEEE 754 encodings in that the entire
49    significand is fractional.  Normalized significands are in the
50    range [0.5, 1.0).
51
52    A requirement of the model is that P be larger than the largest
53    supported target floating-point type by at least 2 bits.  This gives
54    us proper rounding when we truncate to the target type.  In addition,
55    E must be large enough to hold the smallest supported denormal number
56    in a normalized form.
57
58    Both of these requirements are easily satisfied.  The largest target
59    significand is 113 bits; we store at least 160.  The smallest
60    denormal number fits in 17 exponent bits; we store 27.
61
62    Note that the decimal string conversion routines are sensitive to
63    rounding errors.  Since the raw arithmetic routines do not themselves
64    have guard digits or rounding, the computation of 10**exp can
65    accumulate more than a few digits of error.  The previous incarnation
66    of real.c successfully used a 144-bit fraction; given the current
67    layout of REAL_VALUE_TYPE we're forced to expand to at least 160 bits.  */
68
69
70 /* Used to classify two numbers simultaneously.  */
71 #define CLASS2(A, B)  ((A) << 2 | (B))
72
73 #if HOST_BITS_PER_LONG != 64 && HOST_BITS_PER_LONG != 32
74  #error "Some constant folding done by hand to avoid shift count warnings"
75 #endif
76
77 static void get_zero (REAL_VALUE_TYPE *, int);
78 static void get_canonical_qnan (REAL_VALUE_TYPE *, int);
79 static void get_canonical_snan (REAL_VALUE_TYPE *, int);
80 static void get_inf (REAL_VALUE_TYPE *, int);
81 static bool sticky_rshift_significand (REAL_VALUE_TYPE *,
82                                        const REAL_VALUE_TYPE *, unsigned int);
83 static void rshift_significand (REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *,
84                                 unsigned int);
85 static void lshift_significand (REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *,
86                                 unsigned int);
87 static void lshift_significand_1 (REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *);
88 static bool add_significands (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *,
89                               const REAL_VALUE_TYPE *);
90 static bool sub_significands (REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *,
91                               const REAL_VALUE_TYPE *, int);
92 static void neg_significand (REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *);
93 static int cmp_significands (const REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *);
94 static int cmp_significand_0 (const REAL_VALUE_TYPE *);
95 static void set_significand_bit (REAL_VALUE_TYPE *, unsigned int);
96 static void clear_significand_bit (REAL_VALUE_TYPE *, unsigned int);
97 static bool test_significand_bit (REAL_VALUE_TYPE *, unsigned int);
98 static void clear_significand_below (REAL_VALUE_TYPE *, unsigned int);
99 static bool div_significands (REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *,
100                               const REAL_VALUE_TYPE *);
101 static void normalize (REAL_VALUE_TYPE *);
102
103 static bool do_add (REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *,
104                     const REAL_VALUE_TYPE *, int);
105 static bool do_multiply (REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *,
106                          const REAL_VALUE_TYPE *);
107 static bool do_divide (REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *,
108                        const REAL_VALUE_TYPE *);
109 static int do_compare (const REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *, int);
110 static void do_fix_trunc (REAL_VALUE_TYPE *, const REAL_VALUE_TYPE *);
111
112 static unsigned long rtd_divmod (REAL_VALUE_TYPE *, REAL_VALUE_TYPE *);
113
114 static const REAL_VALUE_TYPE * ten_to_ptwo (int);
115 static const REAL_VALUE_TYPE * ten_to_mptwo (int);
116 static const REAL_VALUE_TYPE * real_digit (int);
117 static void times_pten (REAL_VALUE_TYPE *, int);
118
119 static void round_for_format (const struct real_format *, REAL_VALUE_TYPE *);
120 \f
121 /* Initialize R with a positive zero.  */
122
123 static inline void
124 get_zero (REAL_VALUE_TYPE *r, int sign)
125 {
126   memset (r, 0, sizeof (*r));
127   r->sign = sign;
128 }
129
130 /* Initialize R with the canonical quiet NaN.  */
131
132 static inline void
133 get_canonical_qnan (REAL_VALUE_TYPE *r, int sign)
134 {
135   memset (r, 0, sizeof (*r));
136   r->cl = rvc_nan;
137   r->sign = sign;
138   r->canonical = 1;
139 }
140
141 static inline void
142 get_canonical_snan (REAL_VALUE_TYPE *r, int sign)
143 {
144   memset (r, 0, sizeof (*r));
145   r->cl = rvc_nan;
146   r->sign = sign;
147   r->signalling = 1;
148   r->canonical = 1;
149 }
150
151 static inline void
152 get_inf (REAL_VALUE_TYPE *r, int sign)
153 {
154   memset (r, 0, sizeof (*r));
155   r->cl = rvc_inf;
156   r->sign = sign;
157 }
158
159 \f
160 /* Right-shift the significand of A by N bits; put the result in the
161    significand of R.  If any one bits are shifted out, return true.  */
162
163 static bool
164 sticky_rshift_significand (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *a,
165                            unsigned int n)
166 {
167   unsigned long sticky = 0;
168   unsigned int i, ofs = 0;
169
170   if (n >= HOST_BITS_PER_LONG)
171     {
172       for (i = 0, ofs = n / HOST_BITS_PER_LONG; i < ofs; ++i)
173         sticky |= a->sig[i];
174       n &= HOST_BITS_PER_LONG - 1;
175     }
176
177   if (n != 0)
178     {
179       sticky |= a->sig[ofs] & (((unsigned long)1 << n) - 1);
180       for (i = 0; i < SIGSZ; ++i)
181         {
182           r->sig[i]
183             = (((ofs + i >= SIGSZ ? 0 : a->sig[ofs + i]) >> n)
184                | ((ofs + i + 1 >= SIGSZ ? 0 : a->sig[ofs + i + 1])
185                   << (HOST_BITS_PER_LONG - n)));
186         }
187     }
188   else
189     {
190       for (i = 0; ofs + i < SIGSZ; ++i)
191         r->sig[i] = a->sig[ofs + i];
192       for (; i < SIGSZ; ++i)
193         r->sig[i] = 0;
194     }
195
196   return sticky != 0;
197 }
198
199 /* Right-shift the significand of A by N bits; put the result in the
200    significand of R.  */
201
202 static void
203 rshift_significand (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *a,
204                     unsigned int n)
205 {
206   unsigned int i, ofs = n / HOST_BITS_PER_LONG;
207
208   n &= HOST_BITS_PER_LONG - 1;
209   if (n != 0)
210     {
211       for (i = 0; i < SIGSZ; ++i)
212         {
213           r->sig[i]
214             = (((ofs + i >= SIGSZ ? 0 : a->sig[ofs + i]) >> n)
215                | ((ofs + i + 1 >= SIGSZ ? 0 : a->sig[ofs + i + 1])
216                   << (HOST_BITS_PER_LONG - n)));
217         }
218     }
219   else
220     {
221       for (i = 0; ofs + i < SIGSZ; ++i)
222         r->sig[i] = a->sig[ofs + i];
223       for (; i < SIGSZ; ++i)
224         r->sig[i] = 0;
225     }
226 }
227
228 /* Left-shift the significand of A by N bits; put the result in the
229    significand of R.  */
230
231 static void
232 lshift_significand (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *a,
233                     unsigned int n)
234 {
235   unsigned int i, ofs = n / HOST_BITS_PER_LONG;
236
237   n &= HOST_BITS_PER_LONG - 1;
238   if (n == 0)
239     {
240       for (i = 0; ofs + i < SIGSZ; ++i)
241         r->sig[SIGSZ-1-i] = a->sig[SIGSZ-1-i-ofs];
242       for (; i < SIGSZ; ++i)
243         r->sig[SIGSZ-1-i] = 0;
244     }
245   else
246     for (i = 0; i < SIGSZ; ++i)
247       {
248         r->sig[SIGSZ-1-i]
249           = (((ofs + i >= SIGSZ ? 0 : a->sig[SIGSZ-1-i-ofs]) << n)
250              | ((ofs + i + 1 >= SIGSZ ? 0 : a->sig[SIGSZ-1-i-ofs-1])
251                 >> (HOST_BITS_PER_LONG - n)));
252       }
253 }
254
255 /* Likewise, but N is specialized to 1.  */
256
257 static inline void
258 lshift_significand_1 (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *a)
259 {
260   unsigned int i;
261
262   for (i = SIGSZ - 1; i > 0; --i)
263     r->sig[i] = (a->sig[i] << 1) | (a->sig[i-1] >> (HOST_BITS_PER_LONG - 1));
264   r->sig[0] = a->sig[0] << 1;
265 }
266
267 /* Add the significands of A and B, placing the result in R.  Return
268    true if there was carry out of the most significant word.  */
269
270 static inline bool
271 add_significands (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *a,
272                   const REAL_VALUE_TYPE *b)
273 {
274   bool carry = false;
275   int i;
276
277   for (i = 0; i < SIGSZ; ++i)
278     {
279       unsigned long ai = a->sig[i];
280       unsigned long ri = ai + b->sig[i];
281
282       if (carry)
283         {
284           carry = ri < ai;
285           carry |= ++ri == 0;
286         }
287       else
288         carry = ri < ai;
289
290       r->sig[i] = ri;
291     }
292
293   return carry;
294 }
295
296 /* Subtract the significands of A and B, placing the result in R.  CARRY is
297    true if there's a borrow incoming to the least significant word.
298    Return true if there was borrow out of the most significant word.  */
299
300 static inline bool
301 sub_significands (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *a,
302                   const REAL_VALUE_TYPE *b, int carry)
303 {
304   int i;
305
306   for (i = 0; i < SIGSZ; ++i)
307     {
308       unsigned long ai = a->sig[i];
309       unsigned long ri = ai - b->sig[i];
310
311       if (carry)
312         {
313           carry = ri > ai;
314           carry |= ~--ri == 0;
315         }
316       else
317         carry = ri > ai;
318
319       r->sig[i] = ri;
320     }
321
322   return carry;
323 }
324
325 /* Negate the significand A, placing the result in R.  */
326
327 static inline void
328 neg_significand (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *a)
329 {
330   bool carry = true;
331   int i;
332
333   for (i = 0; i < SIGSZ; ++i)
334     {
335       unsigned long ri, ai = a->sig[i];
336
337       if (carry)
338         {
339           if (ai)
340             {
341               ri = -ai;
342               carry = false;
343             }
344           else
345             ri = ai;
346         }
347       else
348         ri = ~ai;
349
350       r->sig[i] = ri;
351     }
352 }
353
354 /* Compare significands.  Return tri-state vs zero.  */
355
356 static inline int
357 cmp_significands (const REAL_VALUE_TYPE *a, const REAL_VALUE_TYPE *b)
358 {
359   int i;
360
361   for (i = SIGSZ - 1; i >= 0; --i)
362     {
363       unsigned long ai = a->sig[i];
364       unsigned long bi = b->sig[i];
365
366       if (ai > bi)
367         return 1;
368       if (ai < bi)
369         return -1;
370     }
371
372   return 0;
373 }
374
375 /* Return true if A is nonzero.  */
376
377 static inline int
378 cmp_significand_0 (const REAL_VALUE_TYPE *a)
379 {
380   int i;
381
382   for (i = SIGSZ - 1; i >= 0; --i)
383     if (a->sig[i])
384       return 1;
385
386   return 0;
387 }
388
389 /* Set bit N of the significand of R.  */
390
391 static inline void
392 set_significand_bit (REAL_VALUE_TYPE *r, unsigned int n)
393 {
394   r->sig[n / HOST_BITS_PER_LONG]
395     |= (unsigned long)1 << (n % HOST_BITS_PER_LONG);
396 }
397
398 /* Clear bit N of the significand of R.  */
399
400 static inline void
401 clear_significand_bit (REAL_VALUE_TYPE *r, unsigned int n)
402 {
403   r->sig[n / HOST_BITS_PER_LONG]
404     &= ~((unsigned long)1 << (n % HOST_BITS_PER_LONG));
405 }
406
407 /* Test bit N of the significand of R.  */
408
409 static inline bool
410 test_significand_bit (REAL_VALUE_TYPE *r, unsigned int n)
411 {
412   /* ??? Compiler bug here if we return this expression directly.
413      The conversion to bool strips the "&1" and we wind up testing
414      e.g. 2 != 0 -> true.  Seen in gcc version 3.2 20020520.  */
415   int t = (r->sig[n / HOST_BITS_PER_LONG] >> (n % HOST_BITS_PER_LONG)) & 1;
416   return t;
417 }
418
419 /* Clear bits 0..N-1 of the significand of R.  */
420
421 static void
422 clear_significand_below (REAL_VALUE_TYPE *r, unsigned int n)
423 {
424   int i, w = n / HOST_BITS_PER_LONG;
425
426   for (i = 0; i < w; ++i)
427     r->sig[i] = 0;
428
429   r->sig[w] &= ~(((unsigned long)1 << (n % HOST_BITS_PER_LONG)) - 1);
430 }
431
432 /* Divide the significands of A and B, placing the result in R.  Return
433    true if the division was inexact.  */
434
435 static inline bool
436 div_significands (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *a,
437                   const REAL_VALUE_TYPE *b)
438 {
439   REAL_VALUE_TYPE u;
440   int i, bit = SIGNIFICAND_BITS - 1;
441   unsigned long msb, inexact;
442
443   u = *a;
444   memset (r->sig, 0, sizeof (r->sig));
445
446   msb = 0;
447   goto start;
448   do
449     {
450       msb = u.sig[SIGSZ-1] & SIG_MSB;
451       lshift_significand_1 (&u, &u);
452     start:
453       if (msb || cmp_significands (&u, b) >= 0)
454         {
455           sub_significands (&u, &u, b, 0);
456           set_significand_bit (r, bit);
457         }
458     }
459   while (--bit >= 0);
460
461   for (i = 0, inexact = 0; i < SIGSZ; i++)
462     inexact |= u.sig[i];
463
464   return inexact != 0;
465 }
466
467 /* Adjust the exponent and significand of R such that the most
468    significant bit is set.  We underflow to zero and overflow to
469    infinity here, without denormals.  (The intermediate representation
470    exponent is large enough to handle target denormals normalized.)  */
471
472 static void
473 normalize (REAL_VALUE_TYPE *r)
474 {
475   int shift = 0, exp;
476   int i, j;
477
478   if (r->decimal)
479     return;
480
481   /* Find the first word that is nonzero.  */
482   for (i = SIGSZ - 1; i >= 0; i--)
483     if (r->sig[i] == 0)
484       shift += HOST_BITS_PER_LONG;
485     else
486       break;
487
488   /* Zero significand flushes to zero.  */
489   if (i < 0)
490     {
491       r->cl = rvc_zero;
492       SET_REAL_EXP (r, 0);
493       return;
494     }
495
496   /* Find the first bit that is nonzero.  */
497   for (j = 0; ; j++)
498     if (r->sig[i] & ((unsigned long)1 << (HOST_BITS_PER_LONG - 1 - j)))
499       break;
500   shift += j;
501
502   if (shift > 0)
503     {
504       exp = REAL_EXP (r) - shift;
505       if (exp > MAX_EXP)
506         get_inf (r, r->sign);
507       else if (exp < -MAX_EXP)
508         get_zero (r, r->sign);
509       else
510         {
511           SET_REAL_EXP (r, exp);
512           lshift_significand (r, r, shift);
513         }
514     }
515 }
516 \f
517 /* Calculate R = A + (SUBTRACT_P ? -B : B).  Return true if the
518    result may be inexact due to a loss of precision.  */
519
520 static bool
521 do_add (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *a,
522         const REAL_VALUE_TYPE *b, int subtract_p)
523 {
524   int dexp, sign, exp;
525   REAL_VALUE_TYPE t;
526   bool inexact = false;
527
528   /* Determine if we need to add or subtract.  */
529   sign = a->sign;
530   subtract_p = (sign ^ b->sign) ^ subtract_p;
531
532   switch (CLASS2 (a->cl, b->cl))
533     {
534     case CLASS2 (rvc_zero, rvc_zero):
535       /* -0 + -0 = -0, -0 - +0 = -0; all other cases yield +0.  */
536       get_zero (r, sign & !subtract_p);
537       return false;
538
539     case CLASS2 (rvc_zero, rvc_normal):
540     case CLASS2 (rvc_zero, rvc_inf):
541     case CLASS2 (rvc_zero, rvc_nan):
542       /* 0 + ANY = ANY.  */
543     case CLASS2 (rvc_normal, rvc_nan):
544     case CLASS2 (rvc_inf, rvc_nan):
545     case CLASS2 (rvc_nan, rvc_nan):
546       /* ANY + NaN = NaN.  */
547     case CLASS2 (rvc_normal, rvc_inf):
548       /* R + Inf = Inf.  */
549       *r = *b;
550       r->sign = sign ^ subtract_p;
551       return false;
552
553     case CLASS2 (rvc_normal, rvc_zero):
554     case CLASS2 (rvc_inf, rvc_zero):
555     case CLASS2 (rvc_nan, rvc_zero):
556       /* ANY + 0 = ANY.  */
557     case CLASS2 (rvc_nan, rvc_normal):
558     case CLASS2 (rvc_nan, rvc_inf):
559       /* NaN + ANY = NaN.  */
560     case CLASS2 (rvc_inf, rvc_normal):
561       /* Inf + R = Inf.  */
562       *r = *a;
563       return false;
564
565     case CLASS2 (rvc_inf, rvc_inf):
566       if (subtract_p)
567         /* Inf - Inf = NaN.  */
568         get_canonical_qnan (r, 0);
569       else
570         /* Inf + Inf = Inf.  */
571         *r = *a;
572       return false;
573
574     case CLASS2 (rvc_normal, rvc_normal):
575       break;
576
577     default:
578       gcc_unreachable ();
579     }
580
581   /* Swap the arguments such that A has the larger exponent.  */
582   dexp = REAL_EXP (a) - REAL_EXP (b);
583   if (dexp < 0)
584     {
585       const REAL_VALUE_TYPE *t;
586       t = a, a = b, b = t;
587       dexp = -dexp;
588       sign ^= subtract_p;
589     }
590   exp = REAL_EXP (a);
591
592   /* If the exponents are not identical, we need to shift the
593      significand of B down.  */
594   if (dexp > 0)
595     {
596       /* If the exponents are too far apart, the significands
597          do not overlap, which makes the subtraction a noop.  */
598       if (dexp >= SIGNIFICAND_BITS)
599         {
600           *r = *a;
601           r->sign = sign;
602           return true;
603         }
604
605       inexact |= sticky_rshift_significand (&t, b, dexp);
606       b = &t;
607     }
608
609   if (subtract_p)
610     {
611       if (sub_significands (r, a, b, inexact))
612         {
613           /* We got a borrow out of the subtraction.  That means that
614              A and B had the same exponent, and B had the larger
615              significand.  We need to swap the sign and negate the
616              significand.  */
617           sign ^= 1;
618           neg_significand (r, r);
619         }
620     }
621   else
622     {
623       if (add_significands (r, a, b))
624         {
625           /* We got carry out of the addition.  This means we need to
626              shift the significand back down one bit and increase the
627              exponent.  */
628           inexact |= sticky_rshift_significand (r, r, 1);
629           r->sig[SIGSZ-1] |= SIG_MSB;
630           if (++exp > MAX_EXP)
631             {
632               get_inf (r, sign);
633               return true;
634             }
635         }
636     }
637
638   r->cl = rvc_normal;
639   r->sign = sign;
640   SET_REAL_EXP (r, exp);
641   /* Zero out the remaining fields.  */
642   r->signalling = 0;
643   r->canonical = 0;
644   r->decimal = 0;
645
646   /* Re-normalize the result.  */
647   normalize (r);
648
649   /* Special case: if the subtraction results in zero, the result
650      is positive.  */
651   if (r->cl == rvc_zero)
652     r->sign = 0;
653   else
654     r->sig[0] |= inexact;
655
656   return inexact;
657 }
658
659 /* Calculate R = A * B.  Return true if the result may be inexact.  */
660
661 static bool
662 do_multiply (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *a,
663              const REAL_VALUE_TYPE *b)
664 {
665   REAL_VALUE_TYPE u, t, *rr;
666   unsigned int i, j, k;
667   int sign = a->sign ^ b->sign;
668   bool inexact = false;
669
670   switch (CLASS2 (a->cl, b->cl))
671     {
672     case CLASS2 (rvc_zero, rvc_zero):
673     case CLASS2 (rvc_zero, rvc_normal):
674     case CLASS2 (rvc_normal, rvc_zero):
675       /* +-0 * ANY = 0 with appropriate sign.  */
676       get_zero (r, sign);
677       return false;
678
679     case CLASS2 (rvc_zero, rvc_nan):
680     case CLASS2 (rvc_normal, rvc_nan):
681     case CLASS2 (rvc_inf, rvc_nan):
682     case CLASS2 (rvc_nan, rvc_nan):
683       /* ANY * NaN = NaN.  */
684       *r = *b;
685       r->sign = sign;
686       return false;
687
688     case CLASS2 (rvc_nan, rvc_zero):
689     case CLASS2 (rvc_nan, rvc_normal):
690     case CLASS2 (rvc_nan, rvc_inf):
691       /* NaN * ANY = NaN.  */
692       *r = *a;
693       r->sign = sign;
694       return false;
695
696     case CLASS2 (rvc_zero, rvc_inf):
697     case CLASS2 (rvc_inf, rvc_zero):
698       /* 0 * Inf = NaN */
699       get_canonical_qnan (r, sign);
700       return false;
701
702     case CLASS2 (rvc_inf, rvc_inf):
703     case CLASS2 (rvc_normal, rvc_inf):
704     case CLASS2 (rvc_inf, rvc_normal):
705       /* Inf * Inf = Inf, R * Inf = Inf */
706       get_inf (r, sign);
707       return false;
708
709     case CLASS2 (rvc_normal, rvc_normal):
710       break;
711
712     default:
713       gcc_unreachable ();
714     }
715
716   if (r == a || r == b)
717     rr = &t;
718   else
719     rr = r;
720   get_zero (rr, 0);
721
722   /* Collect all the partial products.  Since we don't have sure access
723      to a widening multiply, we split each long into two half-words.
724
725      Consider the long-hand form of a four half-word multiplication:
726
727                  A  B  C  D
728               *  E  F  G  H
729              --------------
730                 DE DF DG DH
731              CE CF CG CH
732           BE BF BG BH
733        AE AF AG AH
734
735      We construct partial products of the widened half-word products
736      that are known to not overlap, e.g. DF+DH.  Each such partial
737      product is given its proper exponent, which allows us to sum them
738      and obtain the finished product.  */
739
740   for (i = 0; i < SIGSZ * 2; ++i)
741     {
742       unsigned long ai = a->sig[i / 2];
743       if (i & 1)
744         ai >>= HOST_BITS_PER_LONG / 2;
745       else
746         ai &= ((unsigned long)1 << (HOST_BITS_PER_LONG / 2)) - 1;
747
748       if (ai == 0)
749         continue;
750
751       for (j = 0; j < 2; ++j)
752         {
753           int exp = (REAL_EXP (a) - (2*SIGSZ-1-i)*(HOST_BITS_PER_LONG/2)
754                      + (REAL_EXP (b) - (1-j)*(HOST_BITS_PER_LONG/2)));
755
756           if (exp > MAX_EXP)
757             {
758               get_inf (r, sign);
759               return true;
760             }
761           if (exp < -MAX_EXP)
762             {
763               /* Would underflow to zero, which we shouldn't bother adding.  */
764               inexact = true;
765               continue;
766             }
767
768           memset (&u, 0, sizeof (u));
769           u.cl = rvc_normal;
770           SET_REAL_EXP (&u, exp);
771
772           for (k = j; k < SIGSZ * 2; k += 2)
773             {
774               unsigned long bi = b->sig[k / 2];
775               if (k & 1)
776                 bi >>= HOST_BITS_PER_LONG / 2;
777               else
778                 bi &= ((unsigned long)1 << (HOST_BITS_PER_LONG / 2)) - 1;
779
780               u.sig[k / 2] = ai * bi;
781             }
782
783           normalize (&u);
784           inexact |= do_add (rr, rr, &u, 0);
785         }
786     }
787
788   rr->sign = sign;
789   if (rr != r)
790     *r = t;
791
792   return inexact;
793 }
794
795 /* Calculate R = A / B.  Return true if the result may be inexact.  */
796
797 static bool
798 do_divide (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *a,
799            const REAL_VALUE_TYPE *b)
800 {
801   int exp, sign = a->sign ^ b->sign;
802   REAL_VALUE_TYPE t, *rr;
803   bool inexact;
804
805   switch (CLASS2 (a->cl, b->cl))
806     {
807     case CLASS2 (rvc_zero, rvc_zero):
808       /* 0 / 0 = NaN.  */
809     case CLASS2 (rvc_inf, rvc_inf):
810       /* Inf / Inf = NaN.  */
811       get_canonical_qnan (r, sign);
812       return false;
813
814     case CLASS2 (rvc_zero, rvc_normal):
815     case CLASS2 (rvc_zero, rvc_inf):
816       /* 0 / ANY = 0.  */
817     case CLASS2 (rvc_normal, rvc_inf):
818       /* R / Inf = 0.  */
819       get_zero (r, sign);
820       return false;
821
822     case CLASS2 (rvc_normal, rvc_zero):
823       /* R / 0 = Inf.  */
824     case CLASS2 (rvc_inf, rvc_zero):
825       /* Inf / 0 = Inf.  */
826       get_inf (r, sign);
827       return false;
828
829     case CLASS2 (rvc_zero, rvc_nan):
830     case CLASS2 (rvc_normal, rvc_nan):
831     case CLASS2 (rvc_inf, rvc_nan):
832     case CLASS2 (rvc_nan, rvc_nan):
833       /* ANY / NaN = NaN.  */
834       *r = *b;
835       r->sign = sign;
836       return false;
837
838     case CLASS2 (rvc_nan, rvc_zero):
839     case CLASS2 (rvc_nan, rvc_normal):
840     case CLASS2 (rvc_nan, rvc_inf):
841       /* NaN / ANY = NaN.  */
842       *r = *a;
843       r->sign = sign;
844       return false;
845
846     case CLASS2 (rvc_inf, rvc_normal):
847       /* Inf / R = Inf.  */
848       get_inf (r, sign);
849       return false;
850
851     case CLASS2 (rvc_normal, rvc_normal):
852       break;
853
854     default:
855       gcc_unreachable ();
856     }
857
858   if (r == a || r == b)
859     rr = &t;
860   else
861     rr = r;
862
863   /* Make sure all fields in the result are initialized.  */
864   get_zero (rr, 0);
865   rr->cl = rvc_normal;
866   rr->sign = sign;
867
868   exp = REAL_EXP (a) - REAL_EXP (b) + 1;
869   if (exp > MAX_EXP)
870     {
871       get_inf (r, sign);
872       return true;
873     }
874   if (exp < -MAX_EXP)
875     {
876       get_zero (r, sign);
877       return true;
878     }
879   SET_REAL_EXP (rr, exp);
880
881   inexact = div_significands (rr, a, b);
882
883   /* Re-normalize the result.  */
884   normalize (rr);
885   rr->sig[0] |= inexact;
886
887   if (rr != r)
888     *r = t;
889
890   return inexact;
891 }
892
893 /* Return a tri-state comparison of A vs B.  Return NAN_RESULT if
894    one of the two operands is a NaN.  */
895
896 static int
897 do_compare (const REAL_VALUE_TYPE *a, const REAL_VALUE_TYPE *b,
898             int nan_result)
899 {
900   int ret;
901
902   switch (CLASS2 (a->cl, b->cl))
903     {
904     case CLASS2 (rvc_zero, rvc_zero):
905       /* Sign of zero doesn't matter for compares.  */
906       return 0;
907
908     case CLASS2 (rvc_inf, rvc_zero):
909     case CLASS2 (rvc_inf, rvc_normal):
910     case CLASS2 (rvc_normal, rvc_zero):
911       return (a->sign ? -1 : 1);
912
913     case CLASS2 (rvc_inf, rvc_inf):
914       return -a->sign - -b->sign;
915
916     case CLASS2 (rvc_zero, rvc_normal):
917     case CLASS2 (rvc_zero, rvc_inf):
918     case CLASS2 (rvc_normal, rvc_inf):
919       return (b->sign ? 1 : -1);
920
921     case CLASS2 (rvc_zero, rvc_nan):
922     case CLASS2 (rvc_normal, rvc_nan):
923     case CLASS2 (rvc_inf, rvc_nan):
924     case CLASS2 (rvc_nan, rvc_nan):
925     case CLASS2 (rvc_nan, rvc_zero):
926     case CLASS2 (rvc_nan, rvc_normal):
927     case CLASS2 (rvc_nan, rvc_inf):
928       return nan_result;
929
930     case CLASS2 (rvc_normal, rvc_normal):
931       break;
932
933     default:
934       gcc_unreachable ();
935     }
936
937   if (a->sign != b->sign)
938     return -a->sign - -b->sign;
939
940   if (a->decimal || b->decimal)
941     return decimal_do_compare (a, b, nan_result);
942
943   if (REAL_EXP (a) > REAL_EXP (b))
944     ret = 1;
945   else if (REAL_EXP (a) < REAL_EXP (b))
946     ret = -1;
947   else
948     ret = cmp_significands (a, b);
949
950   return (a->sign ? -ret : ret);
951 }
952
953 /* Return A truncated to an integral value toward zero.  */
954
955 static void
956 do_fix_trunc (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *a)
957 {
958   *r = *a;
959
960   switch (r->cl)
961     {
962     case rvc_zero:
963     case rvc_inf:
964     case rvc_nan:
965       break;
966
967     case rvc_normal:
968       if (r->decimal)
969         {
970           decimal_do_fix_trunc (r, a);
971           return;
972         }
973       if (REAL_EXP (r) <= 0)
974         get_zero (r, r->sign);
975       else if (REAL_EXP (r) < SIGNIFICAND_BITS)
976         clear_significand_below (r, SIGNIFICAND_BITS - REAL_EXP (r));
977       break;
978
979     default:
980       gcc_unreachable ();
981     }
982 }
983
984 /* Perform the binary or unary operation described by CODE.
985    For a unary operation, leave OP1 NULL.  This function returns
986    true if the result may be inexact due to loss of precision.  */
987
988 bool
989 real_arithmetic (REAL_VALUE_TYPE *r, int icode, const REAL_VALUE_TYPE *op0,
990                  const REAL_VALUE_TYPE *op1)
991 {
992   enum tree_code code = icode;
993
994   if (op0->decimal || (op1 && op1->decimal))
995     return decimal_real_arithmetic (r, icode, op0, op1);
996
997   switch (code)
998     {
999     case PLUS_EXPR:
1000       return do_add (r, op0, op1, 0);
1001
1002     case MINUS_EXPR:
1003       return do_add (r, op0, op1, 1);
1004
1005     case MULT_EXPR:
1006       return do_multiply (r, op0, op1);
1007
1008     case RDIV_EXPR:
1009       return do_divide (r, op0, op1);
1010
1011     case MIN_EXPR:
1012       if (op1->cl == rvc_nan)
1013         *r = *op1;
1014       else if (do_compare (op0, op1, -1) < 0)
1015         *r = *op0;
1016       else
1017         *r = *op1;
1018       break;
1019
1020     case MAX_EXPR:
1021       if (op1->cl == rvc_nan)
1022         *r = *op1;
1023       else if (do_compare (op0, op1, 1) < 0)
1024         *r = *op1;
1025       else
1026         *r = *op0;
1027       break;
1028
1029     case NEGATE_EXPR:
1030       *r = *op0;
1031       r->sign ^= 1;
1032       break;
1033
1034     case ABS_EXPR:
1035       *r = *op0;
1036       r->sign = 0;
1037       break;
1038
1039     case FIX_TRUNC_EXPR:
1040       do_fix_trunc (r, op0);
1041       break;
1042
1043     default:
1044       gcc_unreachable ();
1045     }
1046   return false;
1047 }
1048
1049 /* Legacy.  Similar, but return the result directly.  */
1050
1051 REAL_VALUE_TYPE
1052 real_arithmetic2 (int icode, const REAL_VALUE_TYPE *op0,
1053                   const REAL_VALUE_TYPE *op1)
1054 {
1055   REAL_VALUE_TYPE r;
1056   real_arithmetic (&r, icode, op0, op1);
1057   return r;
1058 }
1059
1060 bool
1061 real_compare (int icode, const REAL_VALUE_TYPE *op0,
1062               const REAL_VALUE_TYPE *op1)
1063 {
1064   enum tree_code code = icode;
1065
1066   switch (code)
1067     {
1068     case LT_EXPR:
1069       return do_compare (op0, op1, 1) < 0;
1070     case LE_EXPR:
1071       return do_compare (op0, op1, 1) <= 0;
1072     case GT_EXPR:
1073       return do_compare (op0, op1, -1) > 0;
1074     case GE_EXPR:
1075       return do_compare (op0, op1, -1) >= 0;
1076     case EQ_EXPR:
1077       return do_compare (op0, op1, -1) == 0;
1078     case NE_EXPR:
1079       return do_compare (op0, op1, -1) != 0;
1080     case UNORDERED_EXPR:
1081       return op0->cl == rvc_nan || op1->cl == rvc_nan;
1082     case ORDERED_EXPR:
1083       return op0->cl != rvc_nan && op1->cl != rvc_nan;
1084     case UNLT_EXPR:
1085       return do_compare (op0, op1, -1) < 0;
1086     case UNLE_EXPR:
1087       return do_compare (op0, op1, -1) <= 0;
1088     case UNGT_EXPR:
1089       return do_compare (op0, op1, 1) > 0;
1090     case UNGE_EXPR:
1091       return do_compare (op0, op1, 1) >= 0;
1092     case UNEQ_EXPR:
1093       return do_compare (op0, op1, 0) == 0;
1094     case LTGT_EXPR:
1095       return do_compare (op0, op1, 0) != 0;
1096
1097     default:
1098       gcc_unreachable ();
1099     }
1100 }
1101
1102 /* Return floor log2(R).  */
1103
1104 int
1105 real_exponent (const REAL_VALUE_TYPE *r)
1106 {
1107   switch (r->cl)
1108     {
1109     case rvc_zero:
1110       return 0;
1111     case rvc_inf:
1112     case rvc_nan:
1113       return (unsigned int)-1 >> 1;
1114     case rvc_normal:
1115       return REAL_EXP (r);
1116     default:
1117       gcc_unreachable ();
1118     }
1119 }
1120
1121 /* R = OP0 * 2**EXP.  */
1122
1123 void
1124 real_ldexp (REAL_VALUE_TYPE *r, const REAL_VALUE_TYPE *op0, int exp)
1125 {
1126   *r = *op0;
1127   switch (r->cl)
1128     {
1129     case rvc_zero:
1130     case rvc_inf:
1131     case rvc_nan:
1132       break;
1133
1134     case rvc_normal:
1135       exp += REAL_EXP (op0);
1136       if (exp > MAX_EXP)
1137         get_inf (r, r->sign);
1138       else if (exp < -MAX_EXP)
1139         get_zero (r, r->sign);
1140       else
1141         SET_REAL_EXP (r, exp);
1142       break;
1143
1144     default:
1145       gcc_unreachable ();
1146     }
1147 }
1148
1149 /* Determine whether a floating-point value X is infinite.  */
1150
1151 bool
1152 real_isinf (const REAL_VALUE_TYPE *r)
1153 {
1154   return (r->cl == rvc_inf);
1155 }
1156
1157 /* Determine whether a floating-point value X is a NaN.  */
1158
1159 bool
1160 real_isnan (const REAL_VALUE_TYPE *r)
1161 {
1162   return (r->cl == rvc_nan);
1163 }
1164
1165 /* Determine whether a floating-point value X is finite.  */
1166
1167 bool
1168 real_isfinite (const REAL_VALUE_TYPE *r)
1169 {
1170   return (r->cl != rvc_nan) && (r->cl != rvc_inf);
1171 }
1172
1173 /* Determine whether a floating-point value X is negative.  */
1174
1175 bool
1176 real_isneg (const REAL_VALUE_TYPE *r)
1177 {
1178   return r->sign;
1179 }
1180
1181 /* Determine whether a floating-point value X is minus zero.  */
1182
1183 bool
1184 real_isnegzero (const REAL_VALUE_TYPE *r)
1185 {
1186   return r->sign && r->cl == rvc_zero;
1187 }
1188
1189 /* Compare two floating-point objects for bitwise identity.  */
1190
1191 bool
1192 real_identical (const REAL_VALUE_TYPE *a, const REAL_VALUE_TYPE *b)
1193 {
1194   int i;
1195
1196   if (a->cl != b->cl)
1197     return false;
1198   if (a->sign != b->sign)
1199     return false;
1200
1201   switch (a->cl)
1202     {
1203     case rvc_zero:
1204     case rvc_inf:
1205       return true;
1206
1207     case rvc_normal:
1208       if (a->decimal != b->decimal)
1209         return false;
1210       if (REAL_EXP (a) != REAL_EXP (b))
1211         return false;
1212       break;
1213
1214     case rvc_nan:
1215       if (a->signalling != b->signalling)
1216         return false;
1217       /* The significand is ignored for canonical NaNs.  */
1218       if (a->canonical || b->canonical)
1219         return a->canonical == b->canonical;
1220       break;
1221
1222     default:
1223       gcc_unreachable ();
1224     }
1225
1226   for (i = 0; i < SIGSZ; ++i)
1227     if (a->sig[i] != b->sig[i])
1228       return false;
1229
1230   return true;
1231 }
1232
1233 /* Try to change R into its exact multiplicative inverse in machine
1234    mode MODE.  Return true if successful.  */
1235
1236 bool
1237 exact_real_inverse (enum machine_mode mode, REAL_VALUE_TYPE *r)
1238 {
1239   const REAL_VALUE_TYPE *one = real_digit (1);
1240   REAL_VALUE_TYPE u;
1241   int i;
1242
1243   if (r->cl != rvc_normal)
1244     return false;
1245
1246   /* Check for a power of two: all significand bits zero except the MSB.  */
1247   for (i = 0; i < SIGSZ-1; ++i)
1248     if (r->sig[i] != 0)
1249       return false;
1250   if (r->sig[SIGSZ-1] != SIG_MSB)
1251     return false;
1252
1253   /* Find the inverse and truncate to the required mode.  */
1254   do_divide (&u, one, r);
1255   real_convert (&u, mode, &u);
1256
1257   /* The rounding may have overflowed.  */
1258   if (u.cl != rvc_normal)
1259     return false;
1260   for (i = 0; i < SIGSZ-1; ++i)
1261     if (u.sig[i] != 0)
1262       return false;
1263   if (u.sig[SIGSZ-1] != SIG_MSB)
1264     return false;
1265
1266   *r = u;
1267   return true;
1268 }
1269 \f
1270 /* Render R as an integer.  */
1271
1272 HOST_WIDE_INT
1273 real_to_integer (const REAL_VALUE_TYPE *r)
1274 {
1275   unsigned HOST_WIDE_INT i;
1276
1277   switch (r->cl)
1278     {
1279     case rvc_zero:
1280     underflow:
1281       return 0;
1282
1283     case rvc_inf:
1284     case rvc_nan:
1285     overflow:
1286       i = (unsigned HOST_WIDE_INT) 1 << (HOST_BITS_PER_WIDE_INT - 1);
1287       if (!r->sign)
1288         i--;
1289       return i;
1290
1291     case rvc_normal:
1292       if (r->decimal)
1293         return decimal_real_to_integer (r);
1294
1295       if (REAL_EXP (r) <= 0)
1296         goto underflow;
1297       /* Only force overflow for unsigned overflow.  Signed overflow is
1298          undefined, so it doesn't matter what we return, and some callers
1299          expect to be able to use this routine for both signed and
1300          unsigned conversions.  */
1301       if (REAL_EXP (r) > HOST_BITS_PER_WIDE_INT)
1302         goto overflow;
1303
1304       if (HOST_BITS_PER_WIDE_INT == HOST_BITS_PER_LONG)
1305         i = r->sig[SIGSZ-1];
1306       else 
1307         {
1308           gcc_assert (HOST_BITS_PER_WIDE_INT == 2 * HOST_BITS_PER_LONG);
1309           i = r->sig[SIGSZ-1];
1310           i = i << (HOST_BITS_PER_LONG - 1) << 1;
1311           i |= r->sig[SIGSZ-2];
1312         }
1313
1314       i >>= HOST_BITS_PER_WIDE_INT - REAL_EXP (r);
1315
1316       if (r->sign)
1317         i = -i;
1318       return i;
1319
1320     default:
1321       gcc_unreachable ();
1322     }
1323 }
1324
1325 /* Likewise, but to an integer pair, HI+LOW.  */
1326
1327 void
1328 real_to_integer2 (HOST_WIDE_INT *plow, HOST_WIDE_INT *phigh,
1329                   const REAL_VALUE_TYPE *r)
1330 {
1331   REAL_VALUE_TYPE t;
1332   HOST_WIDE_INT low, high;
1333   int exp;
1334
1335   switch (r->cl)
1336     {
1337     case rvc_zero:
1338     underflow:
1339       low = high = 0;
1340       break;
1341
1342     case rvc_inf:
1343     case rvc_nan:
1344     overflow:
1345       high = (unsigned HOST_WIDE_INT) 1 << (HOST_BITS_PER_WIDE_INT - 1);
1346       if (r->sign)
1347         low = 0;
1348       else
1349         {
1350           high--;
1351           low = -1;
1352         }
1353       break;
1354
1355     case rvc_normal:
1356       if (r->decimal)
1357         { 
1358           decimal_real_to_integer2 (plow, phigh, r);
1359           return;
1360         }
1361         
1362       exp = REAL_EXP (r);
1363       if (exp <= 0)
1364         goto underflow;
1365       /* Only force overflow for unsigned overflow.  Signed overflow is
1366          undefined, so it doesn't matter what we return, and some callers
1367          expect to be able to use this routine for both signed and
1368          unsigned conversions.  */
1369       if (exp > 2*HOST_BITS_PER_WIDE_INT)
1370         goto overflow;
1371
1372       rshift_significand (&t, r, 2*HOST_BITS_PER_WIDE_INT - exp);
1373       if (HOST_BITS_PER_WIDE_INT == HOST_BITS_PER_LONG)
1374         {
1375           high = t.sig[SIGSZ-1];
1376           low = t.sig[SIGSZ-2];
1377         }
1378       else 
1379         {
1380           gcc_assert (HOST_BITS_PER_WIDE_INT == 2*HOST_BITS_PER_LONG);
1381           high = t.sig[SIGSZ-1];
1382           high = high << (HOST_BITS_PER_LONG - 1) << 1;
1383           high |= t.sig[SIGSZ-2];
1384
1385           low = t.sig[SIGSZ-3];
1386           low = low << (HOST_BITS_PER_LONG - 1) << 1;
1387           low |= t.sig[SIGSZ-4];
1388         }
1389
1390       if (r->sign)
1391         {
1392           if (low == 0)
1393             high = -high;
1394           else
1395             low = -low, high = ~high;
1396         }
1397       break;
1398
1399     default:
1400       gcc_unreachable ();
1401     }
1402
1403   *plow = low;
1404   *phigh = high;
1405 }
1406
1407 /* A subroutine of real_to_decimal.  Compute the quotient and remainder
1408    of NUM / DEN.  Return the quotient and place the remainder in NUM.
1409    It is expected that NUM / DEN are close enough that the quotient is
1410    small.  */
1411
1412 static unsigned long
1413 rtd_divmod (REAL_VALUE_TYPE *num, REAL_VALUE_TYPE *den)
1414 {
1415   unsigned long q, msb;
1416   int expn = REAL_EXP (num), expd = REAL_EXP (den);
1417
1418   if (expn < expd)
1419     return 0;
1420
1421   q = msb = 0;
1422   goto start;
1423   do
1424     {
1425       msb = num->sig[SIGSZ-1] & SIG_MSB;
1426       q <<= 1;
1427       lshift_significand_1 (num, num);
1428     start:
1429       if (msb || cmp_significands (num, den) >= 0)
1430         {
1431           sub_significands (num, num, den, 0);
1432           q |= 1;
1433         }
1434     }
1435   while (--expn >= expd);
1436
1437   SET_REAL_EXP (num, expd);
1438   normalize (num);
1439
1440   return q;
1441 }
1442
1443 /* Render R as a decimal floating point constant.  Emit DIGITS significant
1444    digits in the result, bounded by BUF_SIZE.  If DIGITS is 0, choose the
1445    maximum for the representation.  If CROP_TRAILING_ZEROS, strip trailing
1446    zeros.  */
1447
1448 #define M_LOG10_2       0.30102999566398119521
1449
1450 void
1451 real_to_decimal (char *str, const REAL_VALUE_TYPE *r_orig, size_t buf_size,
1452                  size_t digits, int crop_trailing_zeros)
1453 {
1454   const REAL_VALUE_TYPE *one, *ten;
1455   REAL_VALUE_TYPE r, pten, u, v;
1456   int dec_exp, cmp_one, digit;
1457   size_t max_digits;
1458   char *p, *first, *last;
1459   bool sign;
1460
1461   r = *r_orig;
1462   switch (r.cl)
1463     {
1464     case rvc_zero:
1465       strcpy (str, (r.sign ? "-0.0" : "0.0"));
1466       return;
1467     case rvc_normal:
1468       break;
1469     case rvc_inf:
1470       strcpy (str, (r.sign ? "-Inf" : "+Inf"));
1471       return;
1472     case rvc_nan:
1473       /* ??? Print the significand as well, if not canonical?  */
1474       strcpy (str, (r.sign ? "-NaN" : "+NaN"));
1475       return;
1476     default:
1477       gcc_unreachable ();
1478     }
1479
1480   if (r.decimal)
1481     {
1482       decimal_real_to_decimal (str, &r, buf_size, digits, crop_trailing_zeros);
1483       return;
1484     }
1485
1486   /* Bound the number of digits printed by the size of the representation.  */
1487   max_digits = SIGNIFICAND_BITS * M_LOG10_2;
1488   if (digits == 0 || digits > max_digits)
1489     digits = max_digits;
1490
1491   /* Estimate the decimal exponent, and compute the length of the string it
1492      will print as.  Be conservative and add one to account for possible
1493      overflow or rounding error.  */
1494   dec_exp = REAL_EXP (&r) * M_LOG10_2;
1495   for (max_digits = 1; dec_exp ; max_digits++)
1496     dec_exp /= 10;
1497
1498   /* Bound the number of digits printed by the size of the output buffer.  */
1499   max_digits = buf_size - 1 - 1 - 2 - max_digits - 1;
1500   gcc_assert (max_digits <= buf_size);
1501   if (digits > max_digits)
1502     digits = max_digits;
1503
1504   one = real_digit (1);
1505   ten = ten_to_ptwo (0);
1506
1507   sign = r.sign;
1508   r.sign = 0;
1509
1510   dec_exp = 0;
1511   pten = *one;
1512
1513   cmp_one = do_compare (&r, one, 0);
1514   if (cmp_one > 0)
1515     {
1516       int m;
1517
1518       /* Number is greater than one.  Convert significand to an integer
1519          and strip trailing decimal zeros.  */
1520
1521       u = r;
1522       SET_REAL_EXP (&u, SIGNIFICAND_BITS - 1);
1523
1524       /* Largest M, such that 10**2**M fits within SIGNIFICAND_BITS.  */
1525       m = floor_log2 (max_digits);
1526
1527       /* Iterate over the bits of the possible powers of 10 that might
1528          be present in U and eliminate them.  That is, if we find that
1529          10**2**M divides U evenly, keep the division and increase
1530          DEC_EXP by 2**M.  */
1531       do
1532         {
1533           REAL_VALUE_TYPE t;
1534
1535           do_divide (&t, &u, ten_to_ptwo (m));
1536           do_fix_trunc (&v, &t);
1537           if (cmp_significands (&v, &t) == 0)
1538             {
1539               u = t;
1540               dec_exp += 1 << m;
1541             }
1542         }
1543       while (--m >= 0);
1544
1545       /* Revert the scaling to integer that we performed earlier.  */
1546       SET_REAL_EXP (&u, REAL_EXP (&u) + REAL_EXP (&r)
1547                     - (SIGNIFICAND_BITS - 1));
1548       r = u;
1549
1550       /* Find power of 10.  Do this by dividing out 10**2**M when
1551          this is larger than the current remainder.  Fill PTEN with
1552          the power of 10 that we compute.  */
1553       if (REAL_EXP (&r) > 0)
1554         {
1555           m = floor_log2 ((int)(REAL_EXP (&r) * M_LOG10_2)) + 1;
1556           do
1557             {
1558               const REAL_VALUE_TYPE *ptentwo = ten_to_ptwo (m);
1559               if (do_compare (&u, ptentwo, 0) >= 0)
1560                 {
1561                   do_divide (&u, &u, ptentwo);
1562                   do_multiply (&pten, &pten, ptentwo);
1563                   dec_exp += 1 << m;
1564                 }
1565             }
1566           while (--m >= 0);
1567         }
1568       else
1569         /* We managed to divide off enough tens in the above reduction
1570            loop that we've now got a negative exponent.  Fall into the
1571            less-than-one code to compute the proper value for PTEN.  */
1572         cmp_one = -1;
1573     }
1574   if (cmp_one < 0)
1575     {
1576       int m;
1577
1578       /* Number is less than one.  Pad significand with leading
1579          decimal zeros.  */
1580
1581       v = r;
1582       while (1)
1583         {
1584           /* Stop if we'd shift bits off the bottom.  */
1585           if (v.sig[0] & 7)
1586             break;
1587
1588           do_multiply (&u, &v, ten);
1589
1590           /* Stop if we're now >= 1.  */
1591           if (REAL_EXP (&u) > 0)
1592             break;
1593
1594           v = u;
1595           dec_exp -= 1;
1596         }
1597       r = v;
1598
1599       /* Find power of 10.  Do this by multiplying in P=10**2**M when
1600          the current remainder is smaller than 1/P.  Fill PTEN with the
1601          power of 10 that we compute.  */
1602       m = floor_log2 ((int)(-REAL_EXP (&r) * M_LOG10_2)) + 1;
1603       do
1604         {
1605           const REAL_VALUE_TYPE *ptentwo = ten_to_ptwo (m);
1606           const REAL_VALUE_TYPE *ptenmtwo = ten_to_mptwo (m);
1607
1608           if (do_compare (&v, ptenmtwo, 0) <= 0)
1609             {
1610               do_multiply (&v, &v, ptentwo);
1611               do_multiply (&pten, &pten, ptentwo);
1612               dec_exp -= 1 << m;
1613             }
1614         }
1615       while (--m >= 0);
1616
1617       /* Invert the positive power of 10 that we've collected so far.  */
1618       do_divide (&pten, one, &pten);
1619     }
1620
1621   p = str;
1622   if (sign)
1623     *p++ = '-';
1624   first = p++;
1625
1626   /* At this point, PTEN should contain the nearest power of 10 smaller
1627      than R, such that this division produces the first digit.
1628
1629      Using a divide-step primitive that returns the complete integral
1630      remainder avoids the rounding error that would be produced if
1631      we were to use do_divide here and then simply multiply by 10 for
1632      each subsequent digit.  */
1633
1634   digit = rtd_divmod (&r, &pten);
1635
1636   /* Be prepared for error in that division via underflow ...  */
1637   if (digit == 0 && cmp_significand_0 (&r))
1638     {
1639       /* Multiply by 10 and try again.  */
1640       do_multiply (&r, &r, ten);
1641       digit = rtd_divmod (&r, &pten);
1642       dec_exp -= 1;
1643       gcc_assert (digit != 0);
1644     }
1645
1646   /* ... or overflow.  */
1647   if (digit == 10)
1648     {
1649       *p++ = '1';
1650       if (--digits > 0)
1651         *p++ = '0';
1652       dec_exp += 1;
1653     }
1654   else
1655     {
1656       gcc_assert (digit <= 10);
1657       *p++ = digit + '0';
1658     }
1659
1660   /* Generate subsequent digits.  */
1661   while (--digits > 0)
1662     {
1663       do_multiply (&r, &r, ten);
1664       digit = rtd_divmod (&r, &pten);
1665       *p++ = digit + '0';
1666     }
1667   last = p;
1668
1669   /* Generate one more digit with which to do rounding.  */
1670   do_multiply (&r, &r, ten);
1671   digit = rtd_divmod (&r, &pten);
1672
1673   /* Round the result.  */
1674   if (digit == 5)
1675     {
1676       /* Round to nearest.  If R is nonzero there are additional
1677          nonzero digits to be extracted.  */
1678       if (cmp_significand_0 (&r))
1679         digit++;
1680       /* Round to even.  */
1681       else if ((p[-1] - '0') & 1)
1682         digit++;
1683     }
1684   if (digit > 5)
1685     {
1686       while (p > first)
1687         {
1688           digit = *--p;
1689           if (digit == '9')
1690             *p = '0';
1691           else
1692             {
1693               *p = digit + 1;
1694               break;
1695             }
1696         }
1697
1698       /* Carry out of the first digit.  This means we had all 9's and
1699          now have all 0's.  "Prepend" a 1 by overwriting the first 0.  */
1700       if (p == first)
1701         {
1702           first[1] = '1';
1703           dec_exp++;
1704         }
1705     }
1706
1707   /* Insert the decimal point.  */
1708   first[0] = first[1];
1709   first[1] = '.';
1710
1711   /* If requested, drop trailing zeros.  Never crop past "1.0".  */
1712   if (crop_trailing_zeros)
1713     while (last > first + 3 && last[-1] == '0')
1714       last--;
1715
1716   /* Append the exponent.  */
1717   sprintf (last, "e%+d", dec_exp);
1718 }
1719
1720 /* Render R as a hexadecimal floating point constant.  Emit DIGITS
1721    significant digits in the result, bounded by BUF_SIZE.  If DIGITS is 0,
1722    choose the maximum for the representation.  If CROP_TRAILING_ZEROS,
1723    strip trailing zeros.  */
1724
1725 void
1726 real_to_hexadecimal (char *str, const REAL_VALUE_TYPE *r, size_t buf_size,
1727                      size_t digits, int crop_trailing_zeros)
1728 {
1729   int i, j, exp = REAL_EXP (r);
1730   char *p, *first;
1731   char exp_buf[16];
1732   size_t max_digits;
1733
1734   switch (r->cl)
1735     {
1736     case rvc_zero:
1737       exp = 0;
1738       break;
1739     case rvc_normal:
1740       break;
1741     case rvc_inf:
1742       strcpy (str, (r->sign ? "-Inf" : "+Inf"));
1743       return;
1744     case rvc_nan:
1745       /* ??? Print the significand as well, if not canonical?  */
1746       strcpy (str, (r->sign ? "-NaN" : "+NaN"));
1747       return;
1748     default:
1749       gcc_unreachable ();
1750     }
1751
1752   if (r->decimal)
1753     {
1754       /* Hexadecimal format for decimal floats is not interesting. */
1755       strcpy (str, "N/A");
1756       return;
1757     }
1758
1759   if (digits == 0)
1760     digits = SIGNIFICAND_BITS / 4;
1761
1762   /* Bound the number of digits printed by the size of the output buffer.  */
1763
1764   sprintf (exp_buf, "p%+d", exp);
1765   max_digits = buf_size - strlen (exp_buf) - r->sign - 4 - 1;
1766   gcc_assert (max_digits <= buf_size);
1767   if (digits > max_digits)
1768     digits = max_digits;
1769
1770   p = str;
1771   if (r->sign)
1772     *p++ = '-';
1773   *p++ = '0';
1774   *p++ = 'x';
1775   *p++ = '0';
1776   *p++ = '.';
1777   first = p;
1778
1779   for (i = SIGSZ - 1; i >= 0; --i)
1780     for (j = HOST_BITS_PER_LONG - 4; j >= 0; j -= 4)
1781       {
1782         *p++ = "0123456789abcdef"[(r->sig[i] >> j) & 15];
1783         if (--digits == 0)
1784           goto out;
1785       }
1786
1787  out:
1788   if (crop_trailing_zeros)
1789     while (p > first + 1 && p[-1] == '0')
1790       p--;
1791
1792   sprintf (p, "p%+d", exp);
1793 }
1794
1795 /* Initialize R from a decimal or hexadecimal string.  The string is
1796    assumed to have been syntax checked already.  Return -1 if the
1797    value underflows, +1 if overflows, and 0 otherwise. */
1798
1799 int
1800 real_from_string (REAL_VALUE_TYPE *r, const char *str)
1801 {
1802   int exp = 0;
1803   bool sign = false;
1804
1805   get_zero (r, 0);
1806
1807   if (*str == '-')
1808     {
1809       sign = true;
1810       str++;
1811     }
1812   else if (*str == '+')
1813     str++;
1814
1815   if (str[0] == '0' && (str[1] == 'x' || str[1] == 'X'))
1816     {
1817       /* Hexadecimal floating point.  */
1818       int pos = SIGNIFICAND_BITS - 4, d;
1819
1820       str += 2;
1821
1822       while (*str == '0')
1823         str++;
1824       while (1)
1825         {
1826           d = hex_value (*str);
1827           if (d == _hex_bad)
1828             break;
1829           if (pos >= 0)
1830             {
1831               r->sig[pos / HOST_BITS_PER_LONG]
1832                 |= (unsigned long) d << (pos % HOST_BITS_PER_LONG);
1833               pos -= 4;
1834             }
1835           else if (d)
1836             /* Ensure correct rounding by setting last bit if there is
1837                a subsequent nonzero digit.  */
1838             r->sig[0] |= 1;
1839           exp += 4;
1840           str++;
1841         }
1842       if (*str == '.')
1843         {
1844           str++;
1845           if (pos == SIGNIFICAND_BITS - 4)
1846             {
1847               while (*str == '0')
1848                 str++, exp -= 4;
1849             }
1850           while (1)
1851             {
1852               d = hex_value (*str);
1853               if (d == _hex_bad)
1854                 break;
1855               if (pos >= 0)
1856                 {
1857                   r->sig[pos / HOST_BITS_PER_LONG]
1858                     |= (unsigned long) d << (pos % HOST_BITS_PER_LONG);
1859                   pos -= 4;
1860                 }
1861               else if (d)
1862                 /* Ensure correct rounding by setting last bit if there is
1863                    a subsequent nonzero digit.  */
1864                 r->sig[0] |= 1;
1865               str++;
1866             }
1867         }
1868
1869       /* If the mantissa is zero, ignore the exponent.  */
1870       if (!cmp_significand_0 (r))
1871         goto is_a_zero;
1872
1873       if (*str == 'p' || *str == 'P')
1874         {
1875           bool exp_neg = false;
1876
1877           str++;
1878           if (*str == '-')
1879             {
1880               exp_neg = true;
1881               str++;
1882             }
1883           else if (*str == '+')
1884             str++;
1885
1886           d = 0;
1887           while (ISDIGIT (*str))
1888             {
1889               d *= 10;
1890               d += *str - '0';
1891               if (d > MAX_EXP)
1892                 {
1893                   /* Overflowed the exponent.  */
1894                   if (exp_neg)
1895                     goto underflow;
1896                   else
1897                     goto overflow;
1898                 }
1899               str++;
1900             }
1901           if (exp_neg)
1902             d = -d;
1903
1904           exp += d;
1905         }
1906
1907       r->cl = rvc_normal;
1908       SET_REAL_EXP (r, exp);
1909
1910       normalize (r);
1911     }
1912   else
1913     {
1914       /* Decimal floating point.  */
1915       const REAL_VALUE_TYPE *ten = ten_to_ptwo (0);
1916       int d;
1917
1918       while (*str == '0')
1919         str++;
1920       while (ISDIGIT (*str))
1921         {
1922           d = *str++ - '0';
1923           do_multiply (r, r, ten);
1924           if (d)
1925             do_add (r, r, real_digit (d), 0);
1926         }
1927       if (*str == '.')
1928         {
1929           str++;
1930           if (r->cl == rvc_zero)
1931             {
1932               while (*str == '0')
1933                 str++, exp--;
1934             }
1935           while (ISDIGIT (*str))
1936             {
1937               d = *str++ - '0';
1938               do_multiply (r, r, ten);
1939               if (d)
1940                 do_add (r, r, real_digit (d), 0);
1941               exp--;
1942             }
1943         }
1944
1945       /* If the mantissa is zero, ignore the exponent.  */
1946       if (r->cl == rvc_zero)
1947         goto is_a_zero;
1948
1949       if (*str == 'e' || *str == 'E')
1950         {
1951           bool exp_neg = false;
1952
1953           str++;
1954           if (*str == '-')
1955             {
1956               exp_neg = true;
1957               str++;
1958             }
1959           else if (*str == '+')
1960             str++;
1961
1962           d = 0;
1963           while (ISDIGIT (*str))
1964             {
1965               d *= 10;
1966               d += *str - '0';
1967               if (d > MAX_EXP)
1968                 {
1969                   /* Overflowed the exponent.  */
1970                   if (exp_neg)
1971                     goto underflow;
1972                   else
1973                     goto overflow;
1974                 }
1975               str++;
1976             }
1977           if (exp_neg)
1978             d = -d;
1979           exp += d;
1980         }
1981
1982       if (exp)
1983         times_pten (r, exp);
1984     }
1985
1986   r->sign = sign;
1987   return 0;
1988
1989  is_a_zero:
1990   get_zero (r, sign);
1991   return 0;
1992
1993  underflow:
1994   get_zero (r, sign);
1995   return -1;
1996
1997  overflow:
1998   get_inf (r, sign);
1999   return 1;
2000 }
2001
2002 /* Legacy.  Similar, but return the result directly.  */
2003
2004 REAL_VALUE_TYPE
2005 real_from_string2 (const char *s, enum machine_mode mode)
2006 {
2007   REAL_VALUE_TYPE r;
2008
2009   real_from_string (&r, s);
2010   if (mode != VOIDmode)
2011     real_convert (&r, mode, &r);
2012
2013   return r;
2014 }
2015
2016 /* Initialize R from string S and desired MODE. */
2017
2018 void 
2019 real_from_string3 (REAL_VALUE_TYPE *r, const char *s, enum machine_mode mode)
2020 {
2021   if (DECIMAL_FLOAT_MODE_P (mode))
2022     decimal_real_from_string (r, s);
2023   else
2024     real_from_string (r, s);
2025
2026   if (mode != VOIDmode)
2027     real_convert (r, mode, r);  
2028
2029
2030 /* Initialize R from the integer pair HIGH+LOW.  */
2031
2032 void
2033 real_from_integer (REAL_VALUE_TYPE *r, enum machine_mode mode,
2034                    unsigned HOST_WIDE_INT low, HOST_WIDE_INT high,
2035                    int unsigned_p)
2036 {
2037   if (low == 0 && high == 0)
2038     get_zero (r, 0);
2039   else
2040     {
2041       memset (r, 0, sizeof (*r));
2042       r->cl = rvc_normal;
2043       r->sign = high < 0 && !unsigned_p;
2044       SET_REAL_EXP (r, 2 * HOST_BITS_PER_WIDE_INT);
2045
2046       if (r->sign)
2047         {
2048           high = ~high;
2049           if (low == 0)
2050             high += 1;
2051           else
2052             low = -low;
2053         }
2054
2055       if (HOST_BITS_PER_LONG == HOST_BITS_PER_WIDE_INT)
2056         {
2057           r->sig[SIGSZ-1] = high;
2058           r->sig[SIGSZ-2] = low;
2059         }
2060       else
2061         {
2062           gcc_assert (HOST_BITS_PER_LONG*2 == HOST_BITS_PER_WIDE_INT);
2063           r->sig[SIGSZ-1] = high >> (HOST_BITS_PER_LONG - 1) >> 1;
2064           r->sig[SIGSZ-2] = high;
2065           r->sig[SIGSZ-3] = low >> (HOST_BITS_PER_LONG - 1) >> 1;
2066           r->sig[SIGSZ-4] = low;
2067         }
2068
2069       normalize (r);
2070     }
2071
2072   if (mode != VOIDmode)
2073     real_convert (r, mode, r);
2074 }
2075
2076 /* Returns 10**2**N.  */
2077
2078 static const REAL_VALUE_TYPE *
2079 ten_to_ptwo (int n)
2080 {
2081   static REAL_VALUE_TYPE tens[EXP_BITS];
2082
2083   gcc_assert (n >= 0);
2084   gcc_assert (n < EXP_BITS);
2085
2086   if (tens[n].cl == rvc_zero)
2087     {
2088       if (n < (HOST_BITS_PER_WIDE_INT == 64 ? 5 : 4))
2089         {
2090           HOST_WIDE_INT t = 10;
2091           int i;
2092
2093           for (i = 0; i < n; ++i)
2094             t *= t;
2095
2096           real_from_integer (&tens[n], VOIDmode, t, 0, 1);
2097         }
2098       else
2099         {
2100           const REAL_VALUE_TYPE *t = ten_to_ptwo (n - 1);
2101           do_multiply (&tens[n], t, t);
2102         }
2103     }
2104
2105   return &tens[n];
2106 }
2107
2108 /* Returns 10**(-2**N).  */
2109
2110 static const REAL_VALUE_TYPE *
2111 ten_to_mptwo (int n)
2112 {
2113   static REAL_VALUE_TYPE tens[EXP_BITS];
2114
2115   gcc_assert (n >= 0);
2116   gcc_assert (n < EXP_BITS);
2117
2118   if (tens[n].cl == rvc_zero)
2119     do_divide (&tens[n], real_digit (1), ten_to_ptwo (n));
2120
2121   return &tens[n];
2122 }
2123
2124 /* Returns N.  */
2125
2126 static const REAL_VALUE_TYPE *
2127 real_digit (int n)
2128 {
2129   static REAL_VALUE_TYPE num[10];
2130
2131   gcc_assert (n >= 0);
2132   gcc_assert (n <= 9);
2133
2134   if (n > 0 && num[n].cl == rvc_zero)
2135     real_from_integer (&num[n], VOIDmode, n, 0, 1);
2136
2137   return &num[n];
2138 }
2139
2140 /* Multiply R by 10**EXP.  */
2141
2142 static void
2143 times_pten (REAL_VALUE_TYPE *r, int exp)
2144 {
2145   REAL_VALUE_TYPE pten, *rr;
2146   bool negative = (exp < 0);
2147   int i;
2148
2149   if (negative)
2150     {
2151       exp = -exp;
2152       pten = *real_digit (1);
2153       rr = &pten;
2154     }
2155   else
2156     rr = r;
2157
2158   for (i = 0; exp > 0; ++i, exp >>= 1)
2159     if (exp & 1)
2160       do_multiply (rr, rr, ten_to_ptwo (i));
2161
2162   if (negative)
2163     do_divide (r, r, &pten);
2164 }
2165
2166 /* Fills R with +Inf.  */
2167
2168 void
2169 real_inf (REAL_VALUE_TYPE *r)
2170 {
2171   get_inf (r, 0);
2172 }
2173
2174 /* Fills R with a NaN whose significand is described by STR.  If QUIET,
2175    we force a QNaN, else we force an SNaN.  The string, if not empty,
2176    is parsed as a number and placed in the significand.  Return true
2177    if the string was successfully parsed.  */
2178
2179 bool
2180 real_nan (REAL_VALUE_TYPE *r, const char *str, int quiet,
2181           enum machine_mode mode)
2182 {
2183   const struct real_format *fmt;
2184
2185   fmt = REAL_MODE_FORMAT (mode);
2186   gcc_assert (fmt);
2187
2188   if (*str == 0)
2189     {
2190       if (quiet)
2191         get_canonical_qnan (r, 0);
2192       else
2193         get_canonical_snan (r, 0);
2194     }
2195   else
2196     {
2197       int base = 10, d;
2198
2199       memset (r, 0, sizeof (*r));
2200       r->cl = rvc_nan;
2201
2202       /* Parse akin to strtol into the significand of R.  */
2203
2204       while (ISSPACE (*str))
2205         str++;
2206       if (*str == '-')
2207         str++;
2208       else if (*str == '+')
2209         str++;
2210       if (*str == '0')
2211         {
2212           str++;
2213           if (*str == 'x' || *str == 'X')
2214             {
2215               base = 16;
2216               str++;
2217             }
2218           else
2219             base = 8;
2220         }
2221
2222       while ((d = hex_value (*str)) < base)
2223         {
2224           REAL_VALUE_TYPE u;
2225
2226           switch (base)
2227             {
2228             case 8:
2229               lshift_significand (r, r, 3);
2230               break;
2231             case 16:
2232               lshift_significand (r, r, 4);
2233               break;
2234             case 10:
2235               lshift_significand_1 (&u, r);
2236               lshift_significand (r, r, 3);
2237               add_significands (r, r, &u);
2238               break;
2239             default:
2240               gcc_unreachable ();
2241             }
2242
2243           get_zero (&u, 0);
2244           u.sig[0] = d;
2245           add_significands (r, r, &u);
2246
2247           str++;
2248         }
2249
2250       /* Must have consumed the entire string for success.  */
2251       if (*str != 0)
2252         return false;
2253
2254       /* Shift the significand into place such that the bits
2255          are in the most significant bits for the format.  */
2256       lshift_significand (r, r, SIGNIFICAND_BITS - fmt->pnan);
2257
2258       /* Our MSB is always unset for NaNs.  */
2259       r->sig[SIGSZ-1] &= ~SIG_MSB;
2260
2261       /* Force quiet or signalling NaN.  */
2262       r->signalling = !quiet;
2263     }
2264
2265   return true;
2266 }
2267
2268 /* Fills R with the largest finite value representable in mode MODE.
2269    If SIGN is nonzero, R is set to the most negative finite value.  */
2270
2271 void
2272 real_maxval (REAL_VALUE_TYPE *r, int sign, enum machine_mode mode)
2273 {
2274   const struct real_format *fmt;
2275   int np2;
2276
2277   fmt = REAL_MODE_FORMAT (mode);
2278   gcc_assert (fmt);
2279   memset (r, 0, sizeof (*r));
2280   
2281   if (fmt->b == 10)
2282     decimal_real_maxval (r, sign, mode);
2283   else
2284     {
2285       r->cl = rvc_normal;
2286       r->sign = sign;
2287       SET_REAL_EXP (r, fmt->emax);
2288
2289       np2 = SIGNIFICAND_BITS - fmt->p;
2290       memset (r->sig, -1, SIGSZ * sizeof (unsigned long));
2291       clear_significand_below (r, np2);
2292
2293       if (fmt->pnan < fmt->p)
2294         /* This is an IBM extended double format made up of two IEEE
2295            doubles.  The value of the long double is the sum of the
2296            values of the two parts.  The most significant part is
2297            required to be the value of the long double rounded to the
2298            nearest double.  Rounding means we need a slightly smaller
2299            value for LDBL_MAX.  */
2300         clear_significand_bit (r, SIGNIFICAND_BITS - fmt->pnan);
2301     }
2302 }
2303
2304 /* Fills R with 2**N.  */
2305
2306 void
2307 real_2expN (REAL_VALUE_TYPE *r, int n, enum machine_mode fmode)
2308 {
2309   memset (r, 0, sizeof (*r));
2310
2311   n++;
2312   if (n > MAX_EXP)
2313     r->cl = rvc_inf;
2314   else if (n < -MAX_EXP)
2315     ;
2316   else
2317     {
2318       r->cl = rvc_normal;
2319       SET_REAL_EXP (r, n);
2320       r->sig[SIGSZ-1] = SIG_MSB;
2321     }
2322   if (DECIMAL_FLOAT_MODE_P (fmode))
2323     decimal_real_convert (r, fmode, r);
2324 }
2325
2326 \f
2327 static void
2328 round_for_format (const struct real_format *fmt, REAL_VALUE_TYPE *r)
2329 {
2330   int p2, np2, i, w;
2331   unsigned long sticky;
2332   bool guard, lsb;
2333   int emin2m1, emax2;
2334
2335   if (r->decimal)
2336     {
2337       if (fmt->b == 10)
2338         {
2339           decimal_round_for_format (fmt, r);
2340           return;
2341         }
2342       /* FIXME. We can come here via fp_easy_constant
2343          (e.g. -O0 on '_Decimal32 x = 1.0 + 2.0dd'), but have not
2344          investigated whether this convert needs to be here, or
2345          something else is missing. */
2346       decimal_real_convert (r, DFmode, r);
2347     }
2348
2349   p2 = fmt->p;
2350   emin2m1 = fmt->emin - 1;
2351   emax2 = fmt->emax;
2352
2353   np2 = SIGNIFICAND_BITS - p2;
2354   switch (r->cl)
2355     {
2356     underflow:
2357       get_zero (r, r->sign);
2358     case rvc_zero:
2359       if (!fmt->has_signed_zero)
2360         r->sign = 0;
2361       return;
2362
2363     overflow:
2364       get_inf (r, r->sign);
2365     case rvc_inf:
2366       return;
2367
2368     case rvc_nan:
2369       clear_significand_below (r, np2);
2370       return;
2371
2372     case rvc_normal:
2373       break;
2374
2375     default:
2376       gcc_unreachable ();
2377     }
2378
2379   /* Check the range of the exponent.  If we're out of range,
2380      either underflow or overflow.  */
2381   if (REAL_EXP (r) > emax2)
2382     goto overflow;
2383   else if (REAL_EXP (r) <= emin2m1)
2384     {
2385       int diff;
2386
2387       if (!fmt->has_denorm)
2388         {
2389           /* Don't underflow completely until we've had a chance to round.  */
2390           if (REAL_EXP (r) < emin2m1)
2391             goto underflow;
2392         }
2393       else
2394         {
2395           diff = emin2m1 - REAL_EXP (r) + 1;
2396           if (diff > p2)
2397             goto underflow;
2398
2399           /* De-normalize the significand.  */
2400           r->sig[0] |= sticky_rshift_significand (r, r, diff);
2401           SET_REAL_EXP (r, REAL_EXP (r) + diff);
2402         }
2403     }
2404
2405   /* There are P2 true significand bits, followed by one guard bit,
2406      followed by one sticky bit, followed by stuff.  Fold nonzero
2407      stuff into the sticky bit.  */
2408
2409   sticky = 0;
2410   for (i = 0, w = (np2 - 1) / HOST_BITS_PER_LONG; i < w; ++i)
2411     sticky |= r->sig[i];
2412   sticky |=
2413     r->sig[w] & (((unsigned long)1 << ((np2 - 1) % HOST_BITS_PER_LONG)) - 1);
2414
2415   guard = test_significand_bit (r, np2 - 1);
2416   lsb = test_significand_bit (r, np2);
2417
2418   /* Round to even.  */
2419   if (guard && (sticky || lsb))
2420     {
2421       REAL_VALUE_TYPE u;
2422       get_zero (&u, 0);
2423       set_significand_bit (&u, np2);
2424
2425       if (add_significands (r, r, &u))
2426         {
2427           /* Overflow.  Means the significand had been all ones, and
2428              is now all zeros.  Need to increase the exponent, and
2429              possibly re-normalize it.  */
2430           SET_REAL_EXP (r, REAL_EXP (r) + 1);
2431           if (REAL_EXP (r) > emax2)
2432             goto overflow;
2433           r->sig[SIGSZ-1] = SIG_MSB;
2434         }
2435     }
2436
2437   /* Catch underflow that we deferred until after rounding.  */
2438   if (REAL_EXP (r) <= emin2m1)
2439     goto underflow;
2440
2441   /* Clear out trailing garbage.  */
2442   clear_significand_below (r, np2);
2443 }
2444
2445 /* Extend or truncate to a new mode.  */
2446
2447 void
2448 real_convert (REAL_VALUE_TYPE *r, enum machine_mode mode,
2449               const REAL_VALUE_TYPE *a)
2450 {
2451   const struct real_format *fmt;
2452
2453   fmt = REAL_MODE_FORMAT (mode);
2454   gcc_assert (fmt);
2455
2456   *r = *a;
2457
2458   if (a->decimal || fmt->b == 10)
2459     decimal_real_convert (r, mode, a);
2460
2461   round_for_format (fmt, r);
2462
2463   /* round_for_format de-normalizes denormals.  Undo just that part.  */
2464   if (r->cl == rvc_normal)
2465     normalize (r);
2466 }
2467
2468 /* Legacy.  Likewise, except return the struct directly.  */
2469
2470 REAL_VALUE_TYPE
2471 real_value_truncate (enum machine_mode mode, REAL_VALUE_TYPE a)
2472 {
2473   REAL_VALUE_TYPE r;
2474   real_convert (&r, mode, &a);
2475   return r;
2476 }
2477
2478 /* Return true if truncating to MODE is exact.  */
2479
2480 bool
2481 exact_real_truncate (enum machine_mode mode, const REAL_VALUE_TYPE *a)
2482 {
2483   const struct real_format *fmt;
2484   REAL_VALUE_TYPE t;
2485   int emin2m1;
2486
2487   fmt = REAL_MODE_FORMAT (mode);
2488   gcc_assert (fmt);
2489
2490   /* Don't allow conversion to denormals.  */
2491   emin2m1 = fmt->emin - 1;
2492   if (REAL_EXP (a) <= emin2m1)
2493     return false;
2494
2495   /* After conversion to the new mode, the value must be identical.  */
2496   real_convert (&t, mode, a);
2497   return real_identical (&t, a);
2498 }
2499
2500 /* Write R to the given target format.  Place the words of the result
2501    in target word order in BUF.  There are always 32 bits in each
2502    long, no matter the size of the host long.
2503
2504    Legacy: return word 0 for implementing REAL_VALUE_TO_TARGET_SINGLE.  */
2505
2506 long
2507 real_to_target_fmt (long *buf, const REAL_VALUE_TYPE *r_orig,
2508                     const struct real_format *fmt)
2509 {
2510   REAL_VALUE_TYPE r;
2511   long buf1;
2512
2513   r = *r_orig;
2514   round_for_format (fmt, &r);
2515
2516   if (!buf)
2517     buf = &buf1;
2518   (*fmt->encode) (fmt, buf, &r);
2519
2520   return *buf;
2521 }
2522
2523 /* Similar, but look up the format from MODE.  */
2524
2525 long
2526 real_to_target (long *buf, const REAL_VALUE_TYPE *r, enum machine_mode mode)
2527 {
2528   const struct real_format *fmt;
2529
2530   fmt = REAL_MODE_FORMAT (mode);
2531   gcc_assert (fmt);
2532
2533   return real_to_target_fmt (buf, r, fmt);
2534 }
2535
2536 /* Read R from the given target format.  Read the words of the result
2537    in target word order in BUF.  There are always 32 bits in each
2538    long, no matter the size of the host long.  */
2539
2540 void
2541 real_from_target_fmt (REAL_VALUE_TYPE *r, const long *buf,
2542                       const struct real_format *fmt)
2543 {
2544   (*fmt->decode) (fmt, r, buf);
2545 }
2546
2547 /* Similar, but look up the format from MODE.  */
2548
2549 void
2550 real_from_target (REAL_VALUE_TYPE *r, const long *buf, enum machine_mode mode)
2551 {
2552   const struct real_format *fmt;
2553
2554   fmt = REAL_MODE_FORMAT (mode);
2555   gcc_assert (fmt);
2556
2557   (*fmt->decode) (fmt, r, buf);
2558 }
2559
2560 /* Return the number of bits of the largest binary value that the
2561    significand of MODE will hold.  */
2562 /* ??? Legacy.  Should get access to real_format directly.  */
2563
2564 int
2565 significand_size (enum machine_mode mode)
2566 {
2567   const struct real_format *fmt;
2568
2569   fmt = REAL_MODE_FORMAT (mode);
2570   if (fmt == NULL)
2571     return 0;
2572
2573   if (fmt->b == 10)
2574     {
2575       /* Return the size in bits of the largest binary value that can be
2576          held by the decimal coefficient for this mode.  This is one more
2577          than the number of bits required to hold the largest coefficient
2578          of this mode.  */
2579       double log2_10 = 3.3219281;
2580       return fmt->p * log2_10;
2581     }
2582   return fmt->p;
2583 }
2584
2585 /* Return a hash value for the given real value.  */
2586 /* ??? The "unsigned int" return value is intended to be hashval_t,
2587    but I didn't want to pull hashtab.h into real.h.  */
2588
2589 unsigned int
2590 real_hash (const REAL_VALUE_TYPE *r)
2591 {
2592   unsigned int h;
2593   size_t i;
2594
2595   h = r->cl | (r->sign << 2);
2596   switch (r->cl)
2597     {
2598     case rvc_zero:
2599     case rvc_inf:
2600       return h;
2601
2602     case rvc_normal:
2603       h |= REAL_EXP (r) << 3;
2604       break;
2605
2606     case rvc_nan:
2607       if (r->signalling)
2608         h ^= (unsigned int)-1;
2609       if (r->canonical)
2610         return h;
2611       break;
2612
2613     default:
2614       gcc_unreachable ();
2615     }
2616
2617   if (sizeof(unsigned long) > sizeof(unsigned int))
2618     for (i = 0; i < SIGSZ; ++i)
2619       {
2620         unsigned long s = r->sig[i];
2621         h ^= s ^ (s >> (HOST_BITS_PER_LONG / 2));
2622       }
2623   else
2624     for (i = 0; i < SIGSZ; ++i)
2625       h ^= r->sig[i];
2626
2627   return h;
2628 }
2629 \f
2630 /* IEEE single-precision format.  */
2631
2632 static void encode_ieee_single (const struct real_format *fmt,
2633                                 long *, const REAL_VALUE_TYPE *);
2634 static void decode_ieee_single (const struct real_format *,
2635                                 REAL_VALUE_TYPE *, const long *);
2636
2637 static void
2638 encode_ieee_single (const struct real_format *fmt, long *buf,
2639                     const REAL_VALUE_TYPE *r)
2640 {
2641   unsigned long image, sig, exp;
2642   unsigned long sign = r->sign;
2643   bool denormal = (r->sig[SIGSZ-1] & SIG_MSB) == 0;
2644
2645   image = sign << 31;
2646   sig = (r->sig[SIGSZ-1] >> (HOST_BITS_PER_LONG - 24)) & 0x7fffff;
2647
2648   switch (r->cl)
2649     {
2650     case rvc_zero:
2651       break;
2652
2653     case rvc_inf:
2654       if (fmt->has_inf)
2655         image |= 255 << 23;
2656       else
2657         image |= 0x7fffffff;
2658       break;
2659
2660     case rvc_nan:
2661       if (fmt->has_nans)
2662         {
2663           if (r->canonical)
2664             sig = (fmt->canonical_nan_lsbs_set ? (1 << 22) - 1 : 0);
2665           if (r->signalling == fmt->qnan_msb_set)
2666             sig &= ~(1 << 22);
2667           else
2668             sig |= 1 << 22;
2669           if (sig == 0)
2670             sig = 1 << 21;
2671
2672           image |= 255 << 23;
2673           image |= sig;
2674         }
2675       else
2676         image |= 0x7fffffff;
2677       break;
2678
2679     case rvc_normal:
2680       /* Recall that IEEE numbers are interpreted as 1.F x 2**exp,
2681          whereas the intermediate representation is 0.F x 2**exp.
2682          Which means we're off by one.  */
2683       if (denormal)
2684         exp = 0;
2685       else
2686       exp = REAL_EXP (r) + 127 - 1;
2687       image |= exp << 23;
2688       image |= sig;
2689       break;
2690
2691     default:
2692       gcc_unreachable ();
2693     }
2694
2695   buf[0] = image;
2696 }
2697
2698 static void
2699 decode_ieee_single (const struct real_format *fmt, REAL_VALUE_TYPE *r,
2700                     const long *buf)
2701 {
2702   unsigned long image = buf[0] & 0xffffffff;
2703   bool sign = (image >> 31) & 1;
2704   int exp = (image >> 23) & 0xff;
2705
2706   memset (r, 0, sizeof (*r));
2707   image <<= HOST_BITS_PER_LONG - 24;
2708   image &= ~SIG_MSB;
2709
2710   if (exp == 0)
2711     {
2712       if (image && fmt->has_denorm)
2713         {
2714           r->cl = rvc_normal;
2715           r->sign = sign;
2716           SET_REAL_EXP (r, -126);
2717           r->sig[SIGSZ-1] = image << 1;
2718           normalize (r);
2719         }
2720       else if (fmt->has_signed_zero)
2721         r->sign = sign;
2722     }
2723   else if (exp == 255 && (fmt->has_nans || fmt->has_inf))
2724     {
2725       if (image)
2726         {
2727           r->cl = rvc_nan;
2728           r->sign = sign;
2729           r->signalling = (((image >> (HOST_BITS_PER_LONG - 2)) & 1)
2730                            ^ fmt->qnan_msb_set);
2731           r->sig[SIGSZ-1] = image;
2732         }
2733       else
2734         {
2735           r->cl = rvc_inf;
2736           r->sign = sign;
2737         }
2738     }
2739   else
2740     {
2741       r->cl = rvc_normal;
2742       r->sign = sign;
2743       SET_REAL_EXP (r, exp - 127 + 1);
2744       r->sig[SIGSZ-1] = image | SIG_MSB;
2745     }
2746 }
2747
2748 const struct real_format ieee_single_format =
2749   {
2750     encode_ieee_single,
2751     decode_ieee_single,
2752     2,
2753     24,
2754     24,
2755     -125,
2756     128,
2757     31,
2758     31,
2759     true,
2760     true,
2761     true,
2762     true,
2763     true,
2764     false
2765   };
2766
2767 const struct real_format mips_single_format =
2768   {
2769     encode_ieee_single,
2770     decode_ieee_single,
2771     2,
2772     24,
2773     24,
2774     -125,
2775     128,
2776     31,
2777     31,
2778     true,
2779     true,
2780     true,
2781     true,
2782     false,
2783     true
2784   };
2785
2786 const struct real_format motorola_single_format =
2787   {
2788     encode_ieee_single,
2789     decode_ieee_single,
2790     2,
2791     24,
2792     24,
2793     -125,
2794     128,
2795     31,
2796     31,
2797     true,
2798     true,
2799     true,
2800     true,
2801     true,
2802     true
2803   };
2804 \f
2805 /* IEEE double-precision format.  */
2806
2807 static void encode_ieee_double (const struct real_format *fmt,
2808                                 long *, const REAL_VALUE_TYPE *);
2809 static void decode_ieee_double (const struct real_format *,
2810                                 REAL_VALUE_TYPE *, const long *);
2811
2812 static void
2813 encode_ieee_double (const struct real_format *fmt, long *buf,
2814                     const REAL_VALUE_TYPE *r)
2815 {
2816   unsigned long image_lo, image_hi, sig_lo, sig_hi, exp;
2817   bool denormal = (r->sig[SIGSZ-1] & SIG_MSB) == 0;
2818
2819   image_hi = r->sign << 31;
2820   image_lo = 0;
2821
2822   if (HOST_BITS_PER_LONG == 64)
2823     {
2824       sig_hi = r->sig[SIGSZ-1];
2825       sig_lo = (sig_hi >> (64 - 53)) & 0xffffffff;
2826       sig_hi = (sig_hi >> (64 - 53 + 1) >> 31) & 0xfffff;
2827     }
2828   else
2829     {
2830       sig_hi = r->sig[SIGSZ-1];
2831       sig_lo = r->sig[SIGSZ-2];
2832       sig_lo = (sig_hi << 21) | (sig_lo >> 11);
2833       sig_hi = (sig_hi >> 11) & 0xfffff;
2834     }
2835
2836   switch (r->cl)
2837     {
2838     case rvc_zero:
2839       break;
2840
2841     case rvc_inf:
2842       if (fmt->has_inf)
2843         image_hi |= 2047 << 20;
2844       else
2845         {
2846           image_hi |= 0x7fffffff;
2847           image_lo = 0xffffffff;
2848         }
2849       break;
2850
2851     case rvc_nan:
2852       if (fmt->has_nans)
2853         {
2854           if (r->canonical)
2855             {
2856               if (fmt->canonical_nan_lsbs_set)
2857                 {
2858                   sig_hi = (1 << 19) - 1;
2859                   sig_lo = 0xffffffff;
2860                 }
2861               else
2862                 {
2863                   sig_hi = 0;
2864                   sig_lo = 0;
2865                 }
2866             }
2867           if (r->signalling == fmt->qnan_msb_set)
2868             sig_hi &= ~(1 << 19);
2869           else
2870             sig_hi |= 1 << 19;
2871           if (sig_hi == 0 && sig_lo == 0)
2872             sig_hi = 1 << 18;
2873
2874           image_hi |= 2047 << 20;
2875           image_hi |= sig_hi;
2876           image_lo = sig_lo;
2877         }
2878       else
2879         {
2880           image_hi |= 0x7fffffff;
2881           image_lo = 0xffffffff;
2882         }
2883       break;
2884
2885     case rvc_normal:
2886       /* Recall that IEEE numbers are interpreted as 1.F x 2**exp,
2887          whereas the intermediate representation is 0.F x 2**exp.
2888          Which means we're off by one.  */
2889       if (denormal)
2890         exp = 0;
2891       else
2892         exp = REAL_EXP (r) + 1023 - 1;
2893       image_hi |= exp << 20;
2894       image_hi |= sig_hi;
2895       image_lo = sig_lo;
2896       break;
2897
2898     default:
2899       gcc_unreachable ();
2900     }
2901
2902   if (FLOAT_WORDS_BIG_ENDIAN)
2903     buf[0] = image_hi, buf[1] = image_lo;
2904   else
2905     buf[0] = image_lo, buf[1] = image_hi;
2906 }
2907
2908 static void
2909 decode_ieee_double (const struct real_format *fmt, REAL_VALUE_TYPE *r,
2910                     const long *buf)
2911 {
2912   unsigned long image_hi, image_lo;
2913   bool sign;
2914   int exp;
2915
2916   if (FLOAT_WORDS_BIG_ENDIAN)
2917     image_hi = buf[0], image_lo = buf[1];
2918   else
2919     image_lo = buf[0], image_hi = buf[1];
2920   image_lo &= 0xffffffff;
2921   image_hi &= 0xffffffff;
2922
2923   sign = (image_hi >> 31) & 1;
2924   exp = (image_hi >> 20) & 0x7ff;
2925
2926   memset (r, 0, sizeof (*r));
2927
2928   image_hi <<= 32 - 21;
2929   image_hi |= image_lo >> 21;
2930   image_hi &= 0x7fffffff;
2931   image_lo <<= 32 - 21;
2932
2933   if (exp == 0)
2934     {
2935       if ((image_hi || image_lo) && fmt->has_denorm)
2936         {
2937           r->cl = rvc_normal;
2938           r->sign = sign;
2939           SET_REAL_EXP (r, -1022);
2940           if (HOST_BITS_PER_LONG == 32)
2941             {
2942               image_hi = (image_hi << 1) | (image_lo >> 31);
2943               image_lo <<= 1;
2944               r->sig[SIGSZ-1] = image_hi;
2945               r->sig[SIGSZ-2] = image_lo;
2946             }
2947           else
2948             {
2949               image_hi = (image_hi << 31 << 2) | (image_lo << 1);
2950               r->sig[SIGSZ-1] = image_hi;
2951             }
2952           normalize (r);
2953         }
2954       else if (fmt->has_signed_zero)
2955         r->sign = sign;
2956     }
2957   else if (exp == 2047 && (fmt->has_nans || fmt->has_inf))
2958     {
2959       if (image_hi || image_lo)
2960         {
2961           r->cl = rvc_nan;
2962           r->sign = sign;
2963           r->signalling = ((image_hi >> 30) & 1) ^ fmt->qnan_msb_set;
2964           if (HOST_BITS_PER_LONG == 32)
2965             {
2966               r->sig[SIGSZ-1] = image_hi;
2967               r->sig[SIGSZ-2] = image_lo;
2968             }
2969           else
2970             r->sig[SIGSZ-1] = (image_hi << 31 << 1) | image_lo;
2971         }
2972       else
2973         {
2974           r->cl = rvc_inf;
2975           r->sign = sign;
2976         }
2977     }
2978   else
2979     {
2980       r->cl = rvc_normal;
2981       r->sign = sign;
2982       SET_REAL_EXP (r, exp - 1023 + 1);
2983       if (HOST_BITS_PER_LONG == 32)
2984         {
2985           r->sig[SIGSZ-1] = image_hi | SIG_MSB;
2986           r->sig[SIGSZ-2] = image_lo;
2987         }
2988       else
2989         r->sig[SIGSZ-1] = (image_hi << 31 << 1) | image_lo | SIG_MSB;
2990     }
2991 }
2992
2993 const struct real_format ieee_double_format =
2994   {
2995     encode_ieee_double,
2996     decode_ieee_double,
2997     2,
2998     53,
2999     53,
3000     -1021,
3001     1024,
3002     63,
3003     63,
3004     true,
3005     true,
3006     true,
3007     true,
3008     true,
3009     false
3010   };
3011
3012 const struct real_format mips_double_format =
3013   {
3014     encode_ieee_double,
3015     decode_ieee_double,
3016     2,
3017     53,
3018     53,
3019     -1021,
3020     1024,
3021     63,
3022     63,
3023     true,
3024     true,
3025     true,
3026     true,
3027     false,
3028     true
3029   };
3030
3031 const struct real_format motorola_double_format =
3032   {
3033     encode_ieee_double,
3034     decode_ieee_double,
3035     2,
3036     53,
3037     53,
3038     -1021,
3039     1024,
3040     63,
3041     63,
3042     true,
3043     true,
3044     true,
3045     true,
3046     true,
3047     true
3048   };
3049 \f
3050 /* IEEE extended real format.  This comes in three flavors: Intel's as
3051    a 12 byte image, Intel's as a 16 byte image, and Motorola's.  Intel
3052    12- and 16-byte images may be big- or little endian; Motorola's is
3053    always big endian.  */
3054
3055 /* Helper subroutine which converts from the internal format to the
3056    12-byte little-endian Intel format.  Functions below adjust this
3057    for the other possible formats.  */
3058 static void
3059 encode_ieee_extended (const struct real_format *fmt, long *buf,
3060                       const REAL_VALUE_TYPE *r)
3061 {
3062   unsigned long image_hi, sig_hi, sig_lo;
3063   bool denormal = (r->sig[SIGSZ-1] & SIG_MSB) == 0;
3064
3065   image_hi = r->sign << 15;
3066   sig_hi = sig_lo = 0;
3067
3068   switch (r->cl)
3069     {
3070     case rvc_zero:
3071       break;
3072
3073     case rvc_inf:
3074       if (fmt->has_inf)
3075         {
3076           image_hi |= 32767;
3077
3078           /* Intel requires the explicit integer bit to be set, otherwise
3079              it considers the value a "pseudo-infinity".  Motorola docs
3080              say it doesn't care.  */
3081           sig_hi = 0x80000000;
3082         }
3083       else
3084         {
3085           image_hi |= 32767;
3086           sig_lo = sig_hi = 0xffffffff;
3087         }
3088       break;
3089
3090     case rvc_nan:
3091       if (fmt->has_nans)
3092         {
3093           image_hi |= 32767;
3094           if (r->canonical)
3095             {
3096               if (fmt->canonical_nan_lsbs_set)
3097                 {
3098                   sig_hi = (1 << 30) - 1;
3099                   sig_lo = 0xffffffff;
3100                 }
3101             }
3102           else if (HOST_BITS_PER_LONG == 32)
3103             {
3104               sig_hi = r->sig[SIGSZ-1];
3105               sig_lo = r->sig[SIGSZ-2];
3106             }
3107           else
3108             {
3109               sig_lo = r->sig[SIGSZ-1];
3110               sig_hi = sig_lo >> 31 >> 1;
3111               sig_lo &= 0xffffffff;
3112             }
3113           if (r->signalling == fmt->qnan_msb_set)
3114             sig_hi &= ~(1 << 30);
3115           else
3116             sig_hi |= 1 << 30;
3117           if ((sig_hi & 0x7fffffff) == 0 && sig_lo == 0)
3118             sig_hi = 1 << 29;
3119
3120           /* Intel requires the explicit integer bit to be set, otherwise
3121              it considers the value a "pseudo-nan".  Motorola docs say it
3122              doesn't care.  */
3123           sig_hi |= 0x80000000;
3124         }
3125       else
3126         {
3127           image_hi |= 32767;
3128           sig_lo = sig_hi = 0xffffffff;
3129         }
3130       break;
3131
3132     case rvc_normal:
3133       {
3134         int exp = REAL_EXP (r);
3135
3136         /* Recall that IEEE numbers are interpreted as 1.F x 2**exp,
3137            whereas the intermediate representation is 0.F x 2**exp.
3138            Which means we're off by one.
3139
3140            Except for Motorola, which consider exp=0 and explicit
3141            integer bit set to continue to be normalized.  In theory
3142            this discrepancy has been taken care of by the difference
3143            in fmt->emin in round_for_format.  */
3144
3145         if (denormal)
3146           exp = 0;
3147         else
3148           {
3149             exp += 16383 - 1;
3150             gcc_assert (exp >= 0);
3151           }
3152         image_hi |= exp;
3153
3154         if (HOST_BITS_PER_LONG == 32)
3155           {
3156             sig_hi = r->sig[SIGSZ-1];
3157             sig_lo = r->sig[SIGSZ-2];
3158           }
3159         else
3160           {
3161             sig_lo = r->sig[SIGSZ-1];
3162             sig_hi = sig_lo >> 31 >> 1;
3163             sig_lo &= 0xffffffff;
3164           }
3165       }
3166       break;
3167
3168     default:
3169       gcc_unreachable ();
3170     }
3171
3172   buf[0] = sig_lo, buf[1] = sig_hi, buf[2] = image_hi;
3173 }
3174
3175 /* Convert from the internal format to the 12-byte Motorola format
3176    for an IEEE extended real.  */
3177 static void
3178 encode_ieee_extended_motorola (const struct real_format *fmt, long *buf,
3179                                const REAL_VALUE_TYPE *r)
3180 {
3181   long intermed[3];
3182   encode_ieee_extended (fmt, intermed, r);
3183
3184   /* Motorola chips are assumed always to be big-endian.  Also, the
3185      padding in a Motorola extended real goes between the exponent and
3186      the mantissa.  At this point the mantissa is entirely within
3187      elements 0 and 1 of intermed, and the exponent entirely within
3188      element 2, so all we have to do is swap the order around, and
3189      shift element 2 left 16 bits.  */
3190   buf[0] = intermed[2] << 16;
3191   buf[1] = intermed[1];
3192   buf[2] = intermed[0];
3193 }
3194
3195 /* Convert from the internal format to the 12-byte Intel format for
3196    an IEEE extended real.  */
3197 static void
3198 encode_ieee_extended_intel_96 (const struct real_format *fmt, long *buf,
3199                                const REAL_VALUE_TYPE *r)
3200 {
3201   if (FLOAT_WORDS_BIG_ENDIAN)
3202     {
3203       /* All the padding in an Intel-format extended real goes at the high
3204          end, which in this case is after the mantissa, not the exponent.
3205          Therefore we must shift everything down 16 bits.  */
3206       long intermed[3];
3207       encode_ieee_extended (fmt, intermed, r);
3208       buf[0] = ((intermed[2] << 16) | ((unsigned long)(intermed[1] & 0xFFFF0000) >> 16));
3209       buf[1] = ((intermed[1] << 16) | ((unsigned long)(intermed[0] & 0xFFFF0000) >> 16));
3210       buf[2] =  (intermed[0] << 16);
3211     }
3212   else
3213     /* encode_ieee_extended produces what we want directly.  */
3214     encode_ieee_extended (fmt, buf, r);
3215 }
3216
3217 /* Convert from the internal format to the 16-byte Intel format for
3218    an IEEE extended real.  */
3219 static void
3220 encode_ieee_extended_intel_128 (const struct real_format *fmt, long *buf,
3221                                 const REAL_VALUE_TYPE *r)
3222 {
3223   /* All the padding in an Intel-format extended real goes at the high end.  */
3224   encode_ieee_extended_intel_96 (fmt, buf, r);
3225   buf[3] = 0;
3226 }
3227
3228 /* As above, we have a helper function which converts from 12-byte
3229    little-endian Intel format to internal format.  Functions below
3230    adjust for the other possible formats.  */
3231 static void
3232 decode_ieee_extended (const struct real_format *fmt, REAL_VALUE_TYPE *r,
3233                       const long *buf)
3234 {
3235   unsigned long image_hi, sig_hi, sig_lo;
3236   bool sign;
3237   int exp;
3238
3239   sig_lo = buf[0], sig_hi = buf[1], image_hi = buf[2];
3240   sig_lo &= 0xffffffff;
3241   sig_hi &= 0xffffffff;
3242   image_hi &= 0xffffffff;
3243
3244   sign = (image_hi >> 15) & 1;
3245   exp = image_hi & 0x7fff;
3246
3247   memset (r, 0, sizeof (*r));
3248
3249   if (exp == 0)
3250     {
3251       if ((sig_hi || sig_lo) && fmt->has_denorm)
3252         {
3253           r->cl = rvc_normal;
3254           r->sign = sign;
3255
3256           /* When the IEEE format contains a hidden bit, we know that
3257              it's zero at this point, and so shift up the significand
3258              and decrease the exponent to match.  In this case, Motorola
3259              defines the explicit integer bit to be valid, so we don't
3260              know whether the msb is set or not.  */
3261           SET_REAL_EXP (r, fmt->emin);
3262           if (HOST_BITS_PER_LONG == 32)
3263             {
3264               r->sig[SIGSZ-1] = sig_hi;
3265               r->sig[SIGSZ-2] = sig_lo;
3266             }
3267           else
3268             r->sig[SIGSZ-1] = (sig_hi << 31 << 1) | sig_lo;
3269
3270           normalize (r);
3271         }
3272       else if (fmt->has_signed_zero)
3273         r->sign = sign;
3274     }
3275   else if (exp == 32767 && (fmt->has_nans || fmt->has_inf))
3276     {
3277       /* See above re "pseudo-infinities" and "pseudo-nans".
3278          Short summary is that the MSB will likely always be
3279          set, and that we don't care about it.  */
3280       sig_hi &= 0x7fffffff;
3281
3282       if (sig_hi || sig_lo)
3283         {
3284           r->cl = rvc_nan;
3285           r->sign = sign;
3286           r->signalling = ((sig_hi >> 30) & 1) ^ fmt->qnan_msb_set;
3287           if (HOST_BITS_PER_LONG == 32)
3288             {
3289               r->sig[SIGSZ-1] = sig_hi;
3290               r->sig[SIGSZ-2] = sig_lo;
3291             }
3292           else
3293             r->sig[SIGSZ-1] = (sig_hi << 31 << 1) | sig_lo;
3294         }
3295       else
3296         {
3297           r->cl = rvc_inf;
3298           r->sign = sign;
3299         }
3300     }
3301   else
3302     {
3303       r->cl = rvc_normal;
3304       r->sign = sign;
3305       SET_REAL_EXP (r, exp - 16383 + 1);
3306       if (HOST_BITS_PER_LONG == 32)
3307         {
3308           r->sig[SIGSZ-1] = sig_hi;
3309           r->sig[SIGSZ-2] = sig_lo;
3310         }
3311       else
3312         r->sig[SIGSZ-1] = (sig_hi << 31 << 1) | sig_lo;
3313     }
3314 }
3315
3316 /* Convert from the internal format to the 12-byte Motorola format
3317    for an IEEE extended real.  */
3318 static void
3319 decode_ieee_extended_motorola (const struct real_format *fmt, REAL_VALUE_TYPE *r,
3320                                const long *buf)
3321 {
3322   long intermed[3];
3323
3324   /* Motorola chips are assumed always to be big-endian.  Also, the
3325      padding in a Motorola extended real goes between the exponent and
3326      the mantissa; remove it.  */
3327   intermed[0] = buf[2];
3328   intermed[1] = buf[1];
3329   intermed[2] = (unsigned long)buf[0] >> 16;
3330
3331   decode_ieee_extended (fmt, r, intermed);
3332 }
3333
3334 /* Convert from the internal format to the 12-byte Intel format for
3335    an IEEE extended real.  */
3336 static void
3337 decode_ieee_extended_intel_96 (const struct real_format *fmt, REAL_VALUE_TYPE *r,
3338                                const long *buf)
3339 {
3340   if (FLOAT_WORDS_BIG_ENDIAN)
3341     {
3342       /* All the padding in an Intel-format extended real goes at the high
3343          end, which in this case is after the mantissa, not the exponent.
3344          Therefore we must shift everything up 16 bits.  */
3345       long intermed[3];
3346
3347       intermed[0] = (((unsigned long)buf[2] >> 16) | (buf[1] << 16));
3348       intermed[1] = (((unsigned long)buf[1] >> 16) | (buf[0] << 16));
3349       intermed[2] =  ((unsigned long)buf[0] >> 16);
3350
3351       decode_ieee_extended (fmt, r, intermed);
3352     }
3353   else
3354     /* decode_ieee_extended produces what we want directly.  */
3355     decode_ieee_extended (fmt, r, buf);
3356 }
3357
3358 /* Convert from the internal format to the 16-byte Intel format for
3359    an IEEE extended real.  */
3360 static void
3361 decode_ieee_extended_intel_128 (const struct real_format *fmt, REAL_VALUE_TYPE *r,
3362                                 const long *buf)
3363 {
3364   /* All the padding in an Intel-format extended real goes at the high end.  */
3365   decode_ieee_extended_intel_96 (fmt, r, buf);
3366 }
3367
3368 const struct real_format ieee_extended_motorola_format =
3369   {
3370     encode_ieee_extended_motorola,
3371     decode_ieee_extended_motorola,
3372     2,
3373     64,
3374     64,
3375     -16382,
3376     16384,
3377     95,
3378     95,
3379     true,
3380     true,
3381     true,
3382     true,
3383     true,
3384     true
3385   };
3386
3387 const struct real_format ieee_extended_intel_96_format =
3388   {
3389     encode_ieee_extended_intel_96,
3390     decode_ieee_extended_intel_96,
3391     2,
3392     64,
3393     64,
3394     -16381,
3395     16384,
3396     79,
3397     79,
3398     true,
3399     true,
3400     true,
3401     true,
3402     true,
3403     false
3404   };
3405
3406 const struct real_format ieee_extended_intel_128_format =
3407   {
3408     encode_ieee_extended_intel_128,
3409     decode_ieee_extended_intel_128,
3410     2,
3411     64,
3412     64,
3413     -16381,
3414     16384,
3415     79,
3416     79,
3417     true,
3418     true,
3419     true,
3420     true,
3421     true,
3422     false
3423   };
3424
3425 /* The following caters to i386 systems that set the rounding precision
3426    to 53 bits instead of 64, e.g. FreeBSD.  */
3427 const struct real_format ieee_extended_intel_96_round_53_format =
3428   {
3429     encode_ieee_extended_intel_96,
3430     decode_ieee_extended_intel_96,
3431     2,
3432     53,
3433     53,
3434     -16381,
3435     16384,
3436     79,
3437     79,
3438     true,
3439     true,
3440     true,
3441     true,
3442     true,
3443     false
3444   };
3445 \f
3446 /* IBM 128-bit extended precision format: a pair of IEEE double precision
3447    numbers whose sum is equal to the extended precision value.  The number
3448    with greater magnitude is first.  This format has the same magnitude
3449    range as an IEEE double precision value, but effectively 106 bits of
3450    significand precision.  Infinity and NaN are represented by their IEEE
3451    double precision value stored in the first number, the second number is
3452    +0.0 or -0.0 for Infinity and don't-care for NaN.  */
3453
3454 static void encode_ibm_extended (const struct real_format *fmt,
3455                                  long *, const REAL_VALUE_TYPE *);
3456 static void decode_ibm_extended (const struct real_format *,
3457                                  REAL_VALUE_TYPE *, const long *);
3458
3459 static void
3460 encode_ibm_extended (const struct real_format *fmt, long *buf,
3461                      const REAL_VALUE_TYPE *r)
3462 {
3463   REAL_VALUE_TYPE u, normr, v;
3464   const struct real_format *base_fmt;
3465
3466   base_fmt = fmt->qnan_msb_set ? &ieee_double_format : &mips_double_format;
3467
3468   /* Renormlize R before doing any arithmetic on it.  */
3469   normr = *r;
3470   if (normr.cl == rvc_normal)
3471     normalize (&normr);
3472
3473   /* u = IEEE double precision portion of significand.  */
3474   u = normr;
3475   round_for_format (base_fmt, &u);
3476   encode_ieee_double (base_fmt, &buf[0], &u);
3477
3478   if (u.cl == rvc_normal)
3479     {
3480       do_add (&v, &normr, &u, 1);
3481       /* Call round_for_format since we might need to denormalize.  */
3482       round_for_format (base_fmt, &v);
3483       encode_ieee_double (base_fmt, &buf[2], &v);
3484     }
3485   else
3486     {
3487       /* Inf, NaN, 0 are all representable as doubles, so the
3488          least-significant part can be 0.0.  */
3489       buf[2] = 0;
3490       buf[3] = 0;
3491     }
3492 }
3493
3494 static void
3495 decode_ibm_extended (const struct real_format *fmt ATTRIBUTE_UNUSED, REAL_VALUE_TYPE *r,
3496                      const long *buf)
3497 {
3498   REAL_VALUE_TYPE u, v;
3499   const struct real_format *base_fmt;
3500
3501   base_fmt = fmt->qnan_msb_set ? &ieee_double_format : &mips_double_format;
3502   decode_ieee_double (base_fmt, &u, &buf[0]);
3503
3504   if (u.cl != rvc_zero && u.cl != rvc_inf && u.cl != rvc_nan)
3505     {
3506       decode_ieee_double (base_fmt, &v, &buf[2]);
3507       do_add (r, &u, &v, 0);
3508     }
3509   else
3510     *r = u;
3511 }
3512
3513 const struct real_format ibm_extended_format =
3514   {
3515     encode_ibm_extended,
3516     decode_ibm_extended,
3517     2,
3518     53 + 53,
3519     53,
3520     -1021 + 53,
3521     1024,
3522     127,
3523     -1,
3524     true,
3525     true,
3526     true,
3527     true,
3528     true,
3529     false
3530   };
3531
3532 const struct real_format mips_extended_format =
3533   {
3534     encode_ibm_extended,
3535     decode_ibm_extended,
3536     2,
3537     53 + 53,
3538     53,
3539     -1021 + 53,
3540     1024,
3541     127,
3542     -1,
3543     true,
3544     true,
3545     true,
3546     true,
3547     false,
3548     true
3549   };
3550
3551 \f
3552 /* IEEE quad precision format.  */
3553
3554 static void encode_ieee_quad (const struct real_format *fmt,
3555                               long *, const REAL_VALUE_TYPE *);
3556 static void decode_ieee_quad (const struct real_format *,
3557                               REAL_VALUE_TYPE *, const long *);
3558
3559 static void
3560 encode_ieee_quad (const struct real_format *fmt, long *buf,
3561                   const REAL_VALUE_TYPE *r)
3562 {
3563   unsigned long image3, image2, image1, image0, exp;
3564   bool denormal = (r->sig[SIGSZ-1] & SIG_MSB) == 0;
3565   REAL_VALUE_TYPE u;
3566
3567   image3 = r->sign << 31;
3568   image2 = 0;
3569   image1 = 0;
3570   image0 = 0;
3571
3572   rshift_significand (&u, r, SIGNIFICAND_BITS - 113);
3573
3574   switch (r->cl)
3575     {
3576     case rvc_zero:
3577       break;
3578
3579     case rvc_inf:
3580       if (fmt->has_inf)
3581         image3 |= 32767 << 16;
3582       else
3583         {
3584           image3 |= 0x7fffffff;
3585           image2 = 0xffffffff;
3586           image1 = 0xffffffff;
3587           image0 = 0xffffffff;
3588         }
3589       break;
3590
3591     case rvc_nan:
3592       if (fmt->has_nans)
3593         {
3594           image3 |= 32767 << 16;
3595
3596           if (r->canonical)
3597             {
3598               if (fmt->canonical_nan_lsbs_set)
3599                 {
3600                   image3 |= 0x7fff;
3601                   image2 = image1 = image0 = 0xffffffff;
3602                 }
3603             }
3604           else if (HOST_BITS_PER_LONG == 32)
3605             {
3606               image0 = u.sig[0];
3607               image1 = u.sig[1];
3608               image2 = u.sig[2];
3609               image3 |= u.sig[3] & 0xffff;
3610             }
3611           else
3612             {
3613               image0 = u.sig[0];
3614               image1 = image0 >> 31 >> 1;
3615               image2 = u.sig[1];
3616               image3 |= (image2 >> 31 >> 1) & 0xffff;
3617               image0 &= 0xffffffff;
3618               image2 &= 0xffffffff;
3619             }
3620           if (r->signalling == fmt->qnan_msb_set)
3621             image3 &= ~0x8000;
3622           else
3623             image3 |= 0x8000;
3624           if (((image3 & 0xffff) | image2 | image1 | image0) == 0)
3625             image3 |= 0x4000;
3626         }
3627       else
3628         {
3629           image3 |= 0x7fffffff;
3630           image2 = 0xffffffff;
3631           image1 = 0xffffffff;
3632           image0 = 0xffffffff;
3633         }
3634       break;
3635
3636     case rvc_normal:
3637       /* Recall that IEEE numbers are interpreted as 1.F x 2**exp,
3638          whereas the intermediate representation is 0.F x 2**exp.
3639          Which means we're off by one.  */
3640       if (denormal)
3641         exp = 0;
3642       else
3643         exp = REAL_EXP (r) + 16383 - 1;
3644       image3 |= exp << 16;
3645
3646       if (HOST_BITS_PER_LONG == 32)
3647         {
3648           image0 = u.sig[0];
3649           image1 = u.sig[1];
3650           image2 = u.sig[2];
3651           image3 |= u.sig[3] & 0xffff;
3652         }
3653       else
3654         {
3655           image0 = u.sig[0];
3656           image1 = image0 >> 31 >> 1;
3657           image2 = u.sig[1];
3658           image3 |= (image2 >> 31 >> 1) & 0xffff;
3659           image0 &= 0xffffffff;
3660           image2 &= 0xffffffff;
3661         }
3662       break;
3663
3664     default:
3665       gcc_unreachable ();
3666     }
3667
3668   if (FLOAT_WORDS_BIG_ENDIAN)
3669     {
3670       buf[0] = image3;
3671       buf[1] = image2;
3672       buf[2] = image1;
3673       buf[3] = image0;
3674     }
3675   else
3676     {
3677       buf[0] = image0;
3678       buf[1] = image1;
3679       buf[2] = image2;
3680       buf[3] = image3;
3681     }
3682 }
3683
3684 static void
3685 decode_ieee_quad (const struct real_format *fmt, REAL_VALUE_TYPE *r,
3686                   const long *buf)
3687 {
3688   unsigned long image3, image2, image1, image0;
3689   bool sign;
3690   int exp;
3691
3692   if (FLOAT_WORDS_BIG_ENDIAN)
3693     {
3694       image3 = buf[0];
3695       image2 = buf[1];
3696       image1 = buf[2];
3697       image0 = buf[3];
3698     }
3699   else
3700     {
3701       image0 = buf[0];
3702       image1 = buf[1];
3703       image2 = buf[2];
3704       image3 = buf[3];
3705     }
3706   image0 &= 0xffffffff;
3707   image1 &= 0xffffffff;
3708   image2 &= 0xffffffff;
3709
3710   sign = (image3 >> 31) & 1;
3711   exp = (image3 >> 16) & 0x7fff;
3712   image3 &= 0xffff;
3713
3714   memset (r, 0, sizeof (*r));
3715
3716   if (exp == 0)
3717     {
3718       if ((image3 | image2 | image1 | image0) && fmt->has_denorm)
3719         {
3720           r->cl = rvc_normal;
3721           r->sign = sign;
3722
3723           SET_REAL_EXP (r, -16382 + (SIGNIFICAND_BITS - 112));
3724           if (HOST_BITS_PER_LONG == 32)
3725             {
3726               r->sig[0] = image0;
3727               r->sig[1] = image1;
3728               r->sig[2] = image2;
3729               r->sig[3] = image3;
3730             }
3731           else
3732             {
3733               r->sig[0] = (image1 << 31 << 1) | image0;
3734               r->sig[1] = (image3 << 31 << 1) | image2;
3735             }
3736
3737           normalize (r);
3738         }
3739       else if (fmt->has_signed_zero)
3740         r->sign = sign;
3741     }
3742   else if (exp == 32767 && (fmt->has_nans || fmt->has_inf))
3743     {
3744       if (image3 | image2 | image1 | image0)
3745         {
3746           r->cl = rvc_nan;
3747           r->sign = sign;
3748           r->signalling = ((image3 >> 15) & 1) ^ fmt->qnan_msb_set;
3749
3750           if (HOST_BITS_PER_LONG == 32)
3751             {
3752               r->sig[0] = image0;
3753               r->sig[1] = image1;
3754               r->sig[2] = image2;
3755               r->sig[3] = image3;
3756             }
3757           else
3758             {
3759               r->sig[0] = (image1 << 31 << 1) | image0;
3760               r->sig[1] = (image3 << 31 << 1) | image2;
3761             }
3762           lshift_significand (r, r, SIGNIFICAND_BITS - 113);
3763         }
3764       else
3765         {
3766           r->cl = rvc_inf;
3767           r->sign = sign;
3768         }
3769     }
3770   else
3771     {
3772       r->cl = rvc_normal;
3773       r->sign = sign;
3774       SET_REAL_EXP (r, exp - 16383 + 1);
3775
3776       if (HOST_BITS_PER_LONG == 32)
3777         {
3778           r->sig[0] = image0;
3779           r->sig[1] = image1;
3780           r->sig[2] = image2;
3781           r->sig[3] = image3;
3782         }
3783       else
3784         {
3785           r->sig[0] = (image1 << 31 << 1) | image0;
3786