OSDN Git Service

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