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