packaging: Add python3-base dependency
[platform/upstream/gdb.git] / sim / common / sim-fpu.c
1 /* This is a software floating point library which can be used instead
2    of the floating point routines in libgcc1.c for targets without
3    hardware floating point.  */
4
5 /* Copyright 1994-2023 Free Software Foundation, Inc.
6
7 This program is free software; you can redistribute it and/or modify
8 it under the terms of the GNU General Public License as published by
9 the Free Software Foundation; either version 3 of the License, or
10 (at your option) any later version.
11
12 This program is distributed in the hope that it will be useful,
13 but WITHOUT ANY WARRANTY; without even the implied warranty of
14 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
15 GNU General Public License for more details.
16
17 You should have received a copy of the GNU General Public License
18 along with this program.  If not, see <http://www.gnu.org/licenses/>.  */
19
20 /* As a special exception, if you link this library with other files,
21    some of which are compiled with GCC, to produce an executable,
22    this library does not by itself cause the resulting executable
23    to be covered by the GNU General Public License.
24    This exception does not however invalidate any other reasons why
25    the executable file might be covered by the GNU General Public License.  */
26
27 /* This implements IEEE 754 format arithmetic, but does not provide a
28    mechanism for setting the rounding mode, or for generating or handling
29    exceptions.
30
31    The original code by Steve Chamberlain, hacked by Mark Eichin and Jim
32    Wilson, all of Cygnus Support.  */
33
34
35 #ifndef SIM_FPU_C
36 #define SIM_FPU_C
37
38 /* This must come before any other includes.  */
39 #include "defs.h"
40
41 #include <stdlib.h>
42
43 #include "sim-basics.h"
44 #include "sim-fpu.h"
45
46 #include "sim-io.h"
47 #include "sim-assert.h"
48
49 /* Debugging support.
50    If digits is -1, then print all digits.  */
51
52 static void
53 print_bits (uint64_t x,
54             int msbit,
55             int digits,
56             sim_fpu_print_func print,
57             void *arg)
58 {
59   uint64_t bit = LSBIT64 (msbit);
60   int i = 4;
61   while (bit && digits)
62     {
63       if (i == 0)
64         print (arg, ",");
65
66       if ((x & bit))
67         print (arg, "1");
68       else
69         print (arg, "0");
70       bit >>= 1;
71
72       if (digits > 0)
73         digits--;
74       i = (i + 1) % 4;
75     }
76 }
77
78
79
80 /* Quick and dirty conversion between a host double and host 64bit int.  */
81
82 typedef union
83 {
84   double d;
85   uint64_t i;
86 } sim_fpu_map;
87
88
89 /* A packed IEEE floating point number.
90
91    Form is <SIGN:1><BIASEDEXP:NR_EXPBITS><FRAC:NR_FRACBITS> for both
92    32 and 64 bit numbers.  This number is interpreted as:
93
94    Normalized (0 < BIASEDEXP && BIASEDEXP < EXPMAX):
95    (sign ? '-' : '+') 1.<FRAC> x 2 ^ (BIASEDEXP - EXPBIAS)
96
97    Denormalized (0 == BIASEDEXP && FRAC != 0):
98    (sign ? "-" : "+") 0.<FRAC> x 2 ^ (- EXPBIAS)
99
100    Zero (0 == BIASEDEXP && FRAC == 0):
101    (sign ? "-" : "+") 0.0
102
103    Infinity (BIASEDEXP == EXPMAX && FRAC == 0):
104    (sign ? "-" : "+") "infinity"
105
106    SignalingNaN (BIASEDEXP == EXPMAX && FRAC > 0 && FRAC < QUIET_NAN):
107    SNaN.FRAC
108
109    QuietNaN (BIASEDEXP == EXPMAX && FRAC > 0 && FRAC > QUIET_NAN):
110    QNaN.FRAC
111
112    */
113
114 #define NR_EXPBITS  (is_double ?   11 :   8)
115 #define NR_FRACBITS (is_double ?   52 : 23)
116 #define SIGNBIT     (is_double ? MSBIT64 (0) : MSBIT64 (32))
117
118 #define EXPMAX32    (255)
119 #define EXMPAX64    (2047)
120 #define EXPMAX      ((unsigned) (is_double ? EXMPAX64 : EXPMAX32))
121
122 #define EXPBIAS32   (127)
123 #define EXPBIAS64   (1023)
124 #define EXPBIAS     (is_double ? EXPBIAS64 : EXPBIAS32)
125
126 #define QUIET_NAN   LSBIT64 (NR_FRACBITS - 1)
127
128
129
130 /* An unpacked floating point number.
131
132    When unpacked, the fraction of both a 32 and 64 bit floating point
133    number is stored using the same format:
134
135    64 bit - <IMPLICIT_1:1><FRACBITS:52><GUARDS:8><PAD:00>
136    32 bit - <IMPLICIT_1:1><FRACBITS:23><GUARDS:7><PAD:30> */
137
138 #define NR_PAD32    (30)
139 #define NR_PAD64    (0)
140 #define NR_PAD      (is_double ? NR_PAD64 : NR_PAD32)
141 #define PADMASK     (is_double ? 0 : LSMASK64 (NR_PAD32 - 1, 0))
142
143 #define NR_GUARDS32 (7 + NR_PAD32)
144 #define NR_GUARDS64 (8 + NR_PAD64)
145 #define NR_GUARDS  (is_double ? NR_GUARDS64 : NR_GUARDS32)
146 #define GUARDMASK  LSMASK64 (NR_GUARDS - 1, 0)
147
148 #define GUARDMSB   LSBIT64  (NR_GUARDS - 1)
149 #define GUARDLSB   LSBIT64  (NR_PAD)
150 #define GUARDROUND LSMASK64 (NR_GUARDS - 2, 0)
151
152 #define NR_FRAC_GUARD   (60)
153 #define IMPLICIT_1 LSBIT64 (NR_FRAC_GUARD)
154 #define IMPLICIT_2 LSBIT64 (NR_FRAC_GUARD + 1)
155 #define IMPLICIT_4 LSBIT64 (NR_FRAC_GUARD + 2)
156 #define NR_SPARE 2
157
158 #define FRAC32MASK LSMASK64 (63, NR_FRAC_GUARD - 32 + 1)
159
160 #define NORMAL_EXPMIN (-(EXPBIAS)+1)
161
162 #define NORMAL_EXPMAX32 (EXPBIAS32)
163 #define NORMAL_EXPMAX64 (EXPBIAS64)
164 #define NORMAL_EXPMAX (EXPBIAS)
165
166
167 /* Integer constants */
168
169 #define MAX_INT32  ((int64_t) LSMASK64 (30, 0))
170 #define MAX_UINT32 LSMASK64 (31, 0)
171 #define MIN_INT32  ((int64_t) LSMASK64 (63, 31))
172
173 #define MAX_INT64  ((int64_t) LSMASK64 (62, 0))
174 #define MAX_UINT64 LSMASK64 (63, 0)
175 #define MIN_INT64  ((int64_t) LSMASK64 (63, 63))
176
177 #define MAX_INT   (is_64bit ? MAX_INT64  : MAX_INT32)
178 #define MIN_INT   (is_64bit ? MIN_INT64  : MIN_INT32)
179 #define MAX_UINT  (is_64bit ? MAX_UINT64 : MAX_UINT32)
180 #define NR_INTBITS (is_64bit ? 64 : 32)
181
182 /* Squeeze an unpacked sim_fpu struct into a 32/64 bit integer.  */
183 STATIC_INLINE_SIM_FPU (uint64_t)
184 pack_fpu (const sim_fpu *src,
185           int is_double)
186 {
187   int sign;
188   uint64_t exp;
189   uint64_t fraction;
190   uint64_t packed;
191
192   switch (src->class)
193     {
194       /* Create a NaN.  */
195     case sim_fpu_class_qnan:
196       sign = src->sign;
197       exp = EXPMAX;
198       /* Force fraction to correct class.  */
199       fraction = src->fraction;
200       fraction >>= NR_GUARDS;
201       if (sim_fpu_quiet_nan_inverted)
202         fraction |= QUIET_NAN - 1;
203       else
204         fraction |= QUIET_NAN;
205       break;
206     case sim_fpu_class_snan:
207       sign = src->sign;
208       exp = EXPMAX;
209       /* Force fraction to correct class.  */
210       fraction = src->fraction;
211       fraction >>= NR_GUARDS;
212       if (sim_fpu_quiet_nan_inverted)
213         fraction |= QUIET_NAN;
214       else
215         fraction &= ~QUIET_NAN;
216       break;
217     case sim_fpu_class_infinity:
218       sign = src->sign;
219       exp = EXPMAX;
220       fraction = 0;
221       break;
222     case sim_fpu_class_zero:
223       sign = src->sign;
224       exp = 0;
225       fraction = 0;
226       break;
227     case sim_fpu_class_number:
228     case sim_fpu_class_denorm:
229       ASSERT (src->fraction >= IMPLICIT_1);
230       ASSERT (src->fraction < IMPLICIT_2);
231       if (src->normal_exp < NORMAL_EXPMIN)
232         {
233           /* This number's exponent is too low to fit into the bits
234              available in the number We'll denormalize the number by
235              storing zero in the exponent and shift the fraction to
236              the right to make up for it. */
237           int nr_shift = NORMAL_EXPMIN - src->normal_exp;
238           if (nr_shift > NR_FRACBITS)
239             {
240               /* Underflow, just make the number zero.  */
241               sign = src->sign;
242               exp = 0;
243               fraction = 0;
244             }
245           else
246             {
247               sign = src->sign;
248               exp = 0;
249               /* Shift by the value.  */
250               fraction = src->fraction;
251               fraction >>= NR_GUARDS;
252               fraction >>= nr_shift;
253             }
254         }
255       else if (src->normal_exp > NORMAL_EXPMAX)
256         {
257           /* Infinity */
258           sign = src->sign;
259           exp = EXPMAX;
260           fraction = 0;
261         }
262       else
263         {
264           exp = (src->normal_exp + EXPBIAS);
265           sign = src->sign;
266           fraction = src->fraction;
267           /* FIXME: Need to round according to WITH_SIM_FPU_ROUNDING
268              or some such.  */
269           /* Round to nearest: If the guard bits are the all zero, but
270              the first, then we're half way between two numbers,
271              choose the one which makes the lsb of the answer 0.  */
272           if ((fraction & GUARDMASK) == GUARDMSB)
273             {
274               if ((fraction & (GUARDMSB << 1)))
275                 fraction += (GUARDMSB << 1);
276             }
277           else
278             {
279               /* Add a one to the guards to force round to nearest.  */
280               fraction += GUARDROUND;
281             }
282           if ((fraction & IMPLICIT_2)) /* Rounding resulted in carry.  */
283             {
284               exp += 1;
285               fraction >>= 1;
286             }
287           fraction >>= NR_GUARDS;
288           /* When exp == EXPMAX (overflow from carry) fraction must
289              have been made zero.  */
290           ASSERT ((exp == EXPMAX) <= ((fraction & ~IMPLICIT_1) == 0));
291         }
292       break;
293     default:
294       abort ();
295     }
296
297   packed = ((sign ? SIGNBIT : 0)
298              | (exp << NR_FRACBITS)
299              | LSMASKED64 (fraction, NR_FRACBITS - 1, 0));
300
301   /* Trace operation.  */
302 #if 0
303   if (is_double)
304     {
305     }
306   else
307     {
308       printf ("pack_fpu: ");
309       printf ("-> %c%0lX.%06lX\n",
310               LSMASKED32 (packed, 31, 31) ? '8' : '0',
311               (long) LSEXTRACTED32 (packed, 30, 23),
312               (long) LSEXTRACTED32 (packed, 23 - 1, 0));
313     }
314 #endif
315
316   return packed;
317 }
318
319
320 /* Unpack a 32/64 bit integer into a sim_fpu structure.  */
321 STATIC_INLINE_SIM_FPU (void)
322 unpack_fpu (sim_fpu *dst, uint64_t packed, int is_double)
323 {
324   uint64_t fraction = LSMASKED64 (packed, NR_FRACBITS - 1, 0);
325   unsigned exp = LSEXTRACTED64 (packed, NR_EXPBITS + NR_FRACBITS - 1, NR_FRACBITS);
326   int sign = (packed & SIGNBIT) != 0;
327
328   if (exp == 0)
329     {
330       /* Hmm.  Looks like 0 */
331       if (fraction == 0)
332         {
333           /* Tastes like zero.  */
334           dst->class = sim_fpu_class_zero;
335           dst->sign = sign;
336           dst->normal_exp = 0;
337         }
338       else
339         {
340           /* Zero exponent with non zero fraction - it's denormalized,
341              so there isn't a leading implicit one - we'll shift it so
342              it gets one.  */
343           dst->normal_exp = exp - EXPBIAS + 1;
344           dst->class = sim_fpu_class_denorm;
345           dst->sign = sign;
346           fraction <<= NR_GUARDS;
347           while (fraction < IMPLICIT_1)
348             {
349               fraction <<= 1;
350               dst->normal_exp--;
351             }
352           dst->fraction = fraction;
353         }
354     }
355   else if (exp == EXPMAX)
356     {
357       /* Huge exponent*/
358       if (fraction == 0)
359         {
360           /* Attached to a zero fraction - means infinity.  */
361           dst->class = sim_fpu_class_infinity;
362           dst->sign = sign;
363           /* dst->normal_exp = EXPBIAS; */
364           /* dst->fraction = 0; */
365         }
366       else
367         {
368           int qnan;
369
370           /* Non zero fraction, means NaN.  */
371           dst->sign = sign;
372           dst->fraction = (fraction << NR_GUARDS);
373           if (sim_fpu_quiet_nan_inverted)
374             qnan = (fraction & QUIET_NAN) == 0;
375           else
376             qnan = fraction >= QUIET_NAN;
377           if (qnan)
378             dst->class = sim_fpu_class_qnan;
379           else
380             dst->class = sim_fpu_class_snan;
381         }
382     }
383   else
384     {
385       /* Nothing strange about this number.  */
386       dst->class = sim_fpu_class_number;
387       dst->sign = sign;
388       dst->fraction = ((fraction << NR_GUARDS) | IMPLICIT_1);
389       dst->normal_exp = exp - EXPBIAS;
390     }
391
392   /* Trace operation.  */
393 #if 0
394   if (is_double)
395     {
396     }
397   else
398     {
399       printf ("unpack_fpu: %c%02lX.%06lX ->\n",
400               LSMASKED32 (packed, 31, 31) ? '8' : '0',
401               (long) LSEXTRACTED32 (packed, 30, 23),
402               (long) LSEXTRACTED32 (packed, 23 - 1, 0));
403     }
404 #endif
405
406   /* sanity checks */
407   {
408     sim_fpu_map val;
409     val.i = pack_fpu (dst, 1);
410     if (is_double)
411       {
412         ASSERT (val.i == packed);
413       }
414     else
415       {
416         uint32_t val = pack_fpu (dst, 0);
417         uint32_t org = packed;
418         ASSERT (val == org);
419       }
420   }
421 }
422
423
424 /* Convert a floating point into an integer.  */
425 STATIC_INLINE_SIM_FPU (int)
426 fpu2i (int64_t *i,
427        const sim_fpu *s,
428        int is_64bit,
429        sim_fpu_round round)
430 {
431   uint64_t tmp;
432   int shift;
433   int status = 0;
434   if (sim_fpu_is_zero (s))
435     {
436       *i = 0;
437       return 0;
438     }
439   if (sim_fpu_is_snan (s))
440     {
441       *i = MIN_INT; /* FIXME */
442       return sim_fpu_status_invalid_cvi;
443     }
444   if (sim_fpu_is_qnan (s))
445     {
446       *i = MIN_INT; /* FIXME */
447       return sim_fpu_status_invalid_cvi;
448     }
449   /* Map infinity onto MAX_INT...  */
450   if (sim_fpu_is_infinity (s))
451     {
452       *i = s->sign ? MIN_INT : MAX_INT;
453       return sim_fpu_status_invalid_cvi;
454     }
455   /* It is a number, but a small one.  */
456   if (s->normal_exp < 0)
457     {
458       *i = 0;
459       return sim_fpu_status_inexact;
460     }
461   /* Is the floating point MIN_INT or just close? */
462   if (s->sign && s->normal_exp == (NR_INTBITS - 1))
463     {
464       *i = MIN_INT;
465       ASSERT (s->fraction >= IMPLICIT_1);
466       if (s->fraction == IMPLICIT_1)
467         return 0; /* exact */
468       if (is_64bit) /* can't round */
469         return sim_fpu_status_invalid_cvi; /* must be overflow */
470       /* For a 32bit with MAX_INT, rounding is possible.  */
471       switch (round)
472         {
473         case sim_fpu_round_default:
474           abort ();
475         case sim_fpu_round_zero:
476           if ((s->fraction & FRAC32MASK) != IMPLICIT_1)
477             return sim_fpu_status_invalid_cvi;
478           else
479             return sim_fpu_status_inexact;
480           break;
481         case sim_fpu_round_near:
482           {
483             if ((s->fraction & FRAC32MASK) != IMPLICIT_1)
484               return sim_fpu_status_invalid_cvi;
485             else if ((s->fraction & !FRAC32MASK) >= (~FRAC32MASK >> 1))
486               return sim_fpu_status_invalid_cvi;
487             else
488               return sim_fpu_status_inexact;
489           }
490         case sim_fpu_round_up:
491           if ((s->fraction & FRAC32MASK) == IMPLICIT_1)
492             return sim_fpu_status_inexact;
493           else
494             return sim_fpu_status_invalid_cvi;
495         case sim_fpu_round_down:
496           return sim_fpu_status_invalid_cvi;
497         }
498     }
499   /* Would right shifting result in the FRAC being shifted into
500      (through) the integer's sign bit? */
501   if (s->normal_exp > (NR_INTBITS - 2))
502     {
503       *i = s->sign ? MIN_INT : MAX_INT;
504       return sim_fpu_status_invalid_cvi;
505     }
506   /* Normal number, shift it into place.  */
507   tmp = s->fraction;
508   shift = (s->normal_exp - (NR_FRAC_GUARD));
509   if (shift > 0)
510     {
511       tmp <<= shift;
512     }
513   else
514     {
515       shift = -shift;
516       if (tmp & ((SIGNED64 (1) << shift) - 1))
517         status |= sim_fpu_status_inexact;
518       tmp >>= shift;
519     }
520   *i = s->sign ? (-tmp) : (tmp);
521   return status;
522 }
523
524 /* Convert an integer into a floating point.  */
525 STATIC_INLINE_SIM_FPU (int)
526 i2fpu (sim_fpu *f, int64_t i, int is_64bit)
527 {
528   int status = 0;
529   if (i == 0)
530     {
531       f->class = sim_fpu_class_zero;
532       f->sign = 0;
533       f->normal_exp = 0;
534     }
535   else
536     {
537       f->class = sim_fpu_class_number;
538       f->sign = (i < 0);
539       f->normal_exp = NR_FRAC_GUARD;
540
541       if (f->sign)
542         {
543           /* Special case for minint, since there is no corresponding
544              +ve integer representation for it.  */
545           if (i == MIN_INT)
546             {
547               f->fraction = IMPLICIT_1;
548               f->normal_exp = NR_INTBITS - 1;
549             }
550           else
551             f->fraction = (-i);
552         }
553       else
554         f->fraction = i;
555
556       if (f->fraction >= IMPLICIT_2)
557         {
558           do
559             {
560               f->fraction = (f->fraction >> 1) | (f->fraction & 1);
561               f->normal_exp += 1;
562             }
563           while (f->fraction >= IMPLICIT_2);
564         }
565       else if (f->fraction < IMPLICIT_1)
566         {
567           do
568             {
569               f->fraction <<= 1;
570               f->normal_exp -= 1;
571             }
572           while (f->fraction < IMPLICIT_1);
573         }
574     }
575
576   /* trace operation */
577 #if 0
578   {
579     printf ("i2fpu: 0x%08lX ->\n", (long) i);
580   }
581 #endif
582
583   /* sanity check */
584   {
585     int64_t val;
586     fpu2i (&val, f, is_64bit, sim_fpu_round_zero);
587     if (i >= MIN_INT32 && i <= MAX_INT32)
588       {
589         ASSERT (val == i);
590       }
591   }
592
593   return status;
594 }
595
596
597 /* Convert a floating point into an integer.  */
598 STATIC_INLINE_SIM_FPU (int)
599 fpu2u (uint64_t *u, const sim_fpu *s, int is_64bit)
600 {
601   const int is_double = 1;
602   uint64_t tmp;
603   int shift;
604   if (sim_fpu_is_zero (s))
605     {
606       *u = 0;
607       return 0;
608     }
609   if (sim_fpu_is_nan (s))
610     {
611       *u = 0;
612       return 0;
613     }
614   /* It is a negative number.  */
615   if (s->sign)
616     {
617       *u = 0;
618       return 0;
619     }
620   /* Get reasonable MAX_USI_INT...  */
621   if (sim_fpu_is_infinity (s))
622     {
623       *u = MAX_UINT;
624       return 0;
625     }
626   /* It is a number, but a small one.  */
627   if (s->normal_exp < 0)
628     {
629       *u = 0;
630       return 0;
631     }
632   /* overflow */
633   if (s->normal_exp > (NR_INTBITS - 1))
634     {
635       *u = MAX_UINT;
636       return 0;
637     }
638   /* normal number */
639   tmp = (s->fraction & ~PADMASK);
640   shift = (s->normal_exp - (NR_FRACBITS + NR_GUARDS));
641   if (shift > 0)
642     {
643       tmp <<= shift;
644     }
645   else
646     {
647       shift = -shift;
648       tmp >>= shift;
649     }
650   *u = tmp;
651   return 0;
652 }
653
654 /* Convert an unsigned integer into a floating point.  */
655 STATIC_INLINE_SIM_FPU (int)
656 u2fpu (sim_fpu *f, uint64_t u, int is_64bit)
657 {
658   if (u == 0)
659     {
660       f->class = sim_fpu_class_zero;
661       f->sign = 0;
662       f->normal_exp = 0;
663     }
664   else
665     {
666       f->class = sim_fpu_class_number;
667       f->sign = 0;
668       f->normal_exp = NR_FRAC_GUARD;
669       f->fraction = u;
670
671       while (f->fraction < IMPLICIT_1)
672         {
673           f->fraction <<= 1;
674           f->normal_exp -= 1;
675         }
676     }
677   return 0;
678 }
679
680
681 /* register <-> sim_fpu */
682
683 INLINE_SIM_FPU (void)
684 sim_fpu_32to (sim_fpu *f, uint32_t s)
685 {
686   unpack_fpu (f, s, 0);
687 }
688
689
690 INLINE_SIM_FPU (void)
691 sim_fpu_232to (sim_fpu *f, uint32_t h, uint32_t l)
692 {
693   uint64_t s = h;
694   s = (s << 32) | l;
695   unpack_fpu (f, s, 1);
696 }
697
698
699 INLINE_SIM_FPU (void)
700 sim_fpu_64to (sim_fpu *f, uint64_t s)
701 {
702   unpack_fpu (f, s, 1);
703 }
704
705
706 INLINE_SIM_FPU (void)
707 sim_fpu_to32 (uint32_t *s,
708               const sim_fpu *f)
709 {
710   *s = pack_fpu (f, 0);
711 }
712
713
714 INLINE_SIM_FPU (void)
715 sim_fpu_to232 (uint32_t *h, uint32_t *l,
716                const sim_fpu *f)
717 {
718   uint64_t s = pack_fpu (f, 1);
719   *l = s;
720   *h = (s >> 32);
721 }
722
723
724 INLINE_SIM_FPU (void)
725 sim_fpu_to64 (uint64_t *u,
726               const sim_fpu *f)
727 {
728   *u = pack_fpu (f, 1);
729 }
730
731
732 INLINE_SIM_FPU (void)
733 sim_fpu_fractionto (sim_fpu *f,
734                     int sign,
735                     int normal_exp,
736                     uint64_t fraction,
737                     int precision)
738 {
739   int shift = (NR_FRAC_GUARD - precision);
740   f->class = sim_fpu_class_number;
741   f->sign = sign;
742   f->normal_exp = normal_exp;
743   /* Shift the fraction to where sim-fpu expects it.  */
744   if (shift >= 0)
745     f->fraction = (fraction << shift);
746   else
747     f->fraction = (fraction >> -shift);
748   f->fraction |= IMPLICIT_1;
749 }
750
751
752 INLINE_SIM_FPU (uint64_t)
753 sim_fpu_tofraction (const sim_fpu *d,
754                     int precision)
755 {
756   /* We have NR_FRAC_GUARD bits, we want only PRECISION bits.  */
757   int shift = (NR_FRAC_GUARD - precision);
758   uint64_t fraction = (d->fraction & ~IMPLICIT_1);
759   if (shift >= 0)
760     return fraction >> shift;
761   else
762     return fraction << -shift;
763 }
764
765
766 /* Rounding */
767
768 STATIC_INLINE_SIM_FPU (int)
769 do_normal_overflow (sim_fpu *f,
770                     int is_double,
771                     sim_fpu_round round)
772 {
773   switch (round)
774     {
775     case sim_fpu_round_default:
776       return 0;
777     case sim_fpu_round_near:
778       f->class = sim_fpu_class_infinity;
779       break;
780     case sim_fpu_round_up:
781       if (!f->sign)
782         f->class = sim_fpu_class_infinity;
783       break;
784     case sim_fpu_round_down:
785       if (f->sign)
786         f->class = sim_fpu_class_infinity;
787       break;
788     case sim_fpu_round_zero:
789       break;
790     }
791   f->normal_exp = NORMAL_EXPMAX;
792   f->fraction = LSMASK64 (NR_FRAC_GUARD, NR_GUARDS);
793   return (sim_fpu_status_overflow | sim_fpu_status_inexact);
794 }
795
796 STATIC_INLINE_SIM_FPU (int)
797 do_normal_underflow (sim_fpu *f,
798                      int is_double,
799                      sim_fpu_round round)
800 {
801   switch (round)
802     {
803     case sim_fpu_round_default:
804       return 0;
805     case sim_fpu_round_near:
806       f->class = sim_fpu_class_zero;
807       break;
808     case sim_fpu_round_up:
809       if (f->sign)
810         f->class = sim_fpu_class_zero;
811       break;
812     case sim_fpu_round_down:
813       if (!f->sign)
814         f->class = sim_fpu_class_zero;
815       break;
816     case sim_fpu_round_zero:
817       f->class = sim_fpu_class_zero;
818       break;
819     }
820   f->normal_exp = NORMAL_EXPMIN - NR_FRACBITS;
821   f->fraction = IMPLICIT_1;
822   return (sim_fpu_status_inexact | sim_fpu_status_underflow);
823 }
824
825
826
827 /* Round a number using NR_GUARDS.
828    Will return the rounded number or F->FRACTION == 0 when underflow.  */
829
830 STATIC_INLINE_SIM_FPU (int)
831 do_normal_round (sim_fpu *f,
832                  int nr_guards,
833                  sim_fpu_round round)
834 {
835   uint64_t guardmask = LSMASK64 (nr_guards - 1, 0);
836   uint64_t guardmsb = LSBIT64 (nr_guards - 1);
837   uint64_t fraclsb = guardmsb << 1;
838   if ((f->fraction & guardmask))
839     {
840       int status = sim_fpu_status_inexact;
841       switch (round)
842         {
843         case sim_fpu_round_default:
844           return 0;
845         case sim_fpu_round_near:
846           if ((f->fraction & guardmsb))
847             {
848               if ((f->fraction & fraclsb))
849                 {
850                   status |= sim_fpu_status_rounded;
851                 }
852               else if ((f->fraction & (guardmask >> 1)))
853                 {
854                   status |= sim_fpu_status_rounded;
855                 }
856             }
857           break;
858         case sim_fpu_round_up:
859           if (!f->sign)
860             status |= sim_fpu_status_rounded;
861           break;
862         case sim_fpu_round_down:
863           if (f->sign)
864             status |= sim_fpu_status_rounded;
865           break;
866         case sim_fpu_round_zero:
867           break;
868         }
869       f->fraction &= ~guardmask;
870       /* Round if needed, handle resulting overflow.  */
871       if ((status & sim_fpu_status_rounded))
872         {
873           f->fraction += fraclsb;
874           if ((f->fraction & IMPLICIT_2))
875             {
876               f->fraction >>= 1;
877               f->normal_exp += 1;
878             }
879         }
880       return status;
881     }
882   else
883     return 0;
884 }
885
886
887 STATIC_INLINE_SIM_FPU (int)
888 do_round (sim_fpu *f,
889           int is_double,
890           sim_fpu_round round,
891           sim_fpu_denorm denorm)
892 {
893   switch (f->class)
894     {
895     case sim_fpu_class_qnan:
896     case sim_fpu_class_zero:
897     case sim_fpu_class_infinity:
898       return 0;
899       break;
900     case sim_fpu_class_snan:
901       /* Quieten a SignalingNaN.  */
902       f->class = sim_fpu_class_qnan;
903       return sim_fpu_status_invalid_snan;
904       break;
905     case sim_fpu_class_number:
906     case sim_fpu_class_denorm:
907       {
908         int status;
909         ASSERT (f->fraction < IMPLICIT_2);
910         ASSERT (f->fraction >= IMPLICIT_1);
911         if (f->normal_exp < NORMAL_EXPMIN)
912           {
913             /* This number's exponent is too low to fit into the bits
914                available in the number.  Round off any bits that will be
915                discarded as a result of denormalization.  Edge case is
916                the implicit bit shifted to GUARD0 and then rounded
917                up. */
918             int shift = NORMAL_EXPMIN - f->normal_exp;
919             if (shift + NR_GUARDS <= NR_FRAC_GUARD + 1
920                 && !(denorm & sim_fpu_denorm_zero))
921               {
922                 status = do_normal_round (f, shift + NR_GUARDS, round);
923                 if (f->fraction == 0) /* Rounding underflowed.  */
924                   {
925                     status |= do_normal_underflow (f, is_double, round);
926                   }
927                 else if (f->normal_exp < NORMAL_EXPMIN) /* still underflow? */
928                   {
929                     status |= sim_fpu_status_denorm;
930                     /* Any loss of precision when denormalizing is
931                        underflow. Some processors check for underflow
932                        before rounding, some after! */
933                     if (status & sim_fpu_status_inexact)
934                       status |= sim_fpu_status_underflow;
935                     /* Flag that resultant value has been denormalized.  */
936                     f->class = sim_fpu_class_denorm;
937                   }
938                 else if ((denorm & sim_fpu_denorm_underflow_inexact))
939                   {
940                     if ((status & sim_fpu_status_inexact))
941                       status |= sim_fpu_status_underflow;
942                   }
943               }
944             else
945               {
946                 status = do_normal_underflow (f, is_double, round);
947               }
948           }
949         else if (f->normal_exp > NORMAL_EXPMAX)
950           {
951             /* Infinity */
952             status = do_normal_overflow (f, is_double, round);
953           }
954         else
955           {
956             status = do_normal_round (f, NR_GUARDS, round);
957             if (f->fraction == 0)
958               /* f->class = sim_fpu_class_zero; */
959               status |= do_normal_underflow (f, is_double, round);
960             else if (f->normal_exp > NORMAL_EXPMAX)
961               /* Oops! rounding caused overflow.  */
962               status |= do_normal_overflow (f, is_double, round);
963           }
964         ASSERT ((f->class == sim_fpu_class_number
965                  || f->class == sim_fpu_class_denorm)
966                 <= (f->fraction < IMPLICIT_2 && f->fraction >= IMPLICIT_1));
967         return status;
968       }
969     }
970   return 0;
971 }
972
973 INLINE_SIM_FPU (int)
974 sim_fpu_round_32 (sim_fpu *f,
975                   sim_fpu_round round,
976                   sim_fpu_denorm denorm)
977 {
978   return do_round (f, 0, round, denorm);
979 }
980
981 INLINE_SIM_FPU (int)
982 sim_fpu_round_64 (sim_fpu *f,
983                   sim_fpu_round round,
984                   sim_fpu_denorm denorm)
985 {
986   return do_round (f, 1, round, denorm);
987 }
988
989 /* NaN handling for binary operations.  */
990
991 INLINE_SIM_FPU (int)
992 sim_fpu_op_nan (sim_fpu *f, const sim_fpu *l, const sim_fpu *r)
993 {
994   if (sim_fpu_is_snan (l) || sim_fpu_is_snan (r))
995     {
996       *f = sim_fpu_is_snan (l) ? *l : *r;
997       f->class = sim_fpu_class_qnan;
998       return sim_fpu_status_invalid_snan;
999     }
1000    ASSERT (sim_fpu_is_nan (l) || sim_fpu_is_nan (r));
1001    if (sim_fpu_is_qnan (l))
1002      *f = *l;
1003    else /* if (sim_fpu_is_qnan (r)) */
1004      *f = *r;
1005   return 0;
1006 }
1007
1008 /* NaN handling specific to min/max operations.  */
1009
1010 INLINE_SIM_FPU (int)
1011 sim_fpu_minmax_nan (sim_fpu *f, const sim_fpu *l, const sim_fpu *r)
1012 {
1013   if (sim_fpu_is_snan (l)
1014       || sim_fpu_is_snan (r)
1015       || sim_fpu_is_ieee754_1985 ())
1016     return sim_fpu_op_nan (f, l, r);
1017   else
1018     /* if sim_fpu_is_ieee754_2008()
1019        && ((sim_fpu_is_qnan (l) || sim_fpu_is_qnan (r)))  */
1020     {
1021       /* In IEEE754-2008:
1022          "minNum/maxNum is ... the canonicalized number if one
1023          operand is a number and the other a quiet NaN."  */
1024       if (sim_fpu_is_qnan (l))
1025         *f = *r;
1026       else /* if (sim_fpu_is_qnan (r))  */
1027         *f = *l;
1028       return 0;
1029     }
1030 }
1031
1032 /* Arithmetic ops */
1033
1034 INLINE_SIM_FPU (int)
1035 sim_fpu_add (sim_fpu *f,
1036              const sim_fpu *l,
1037              const sim_fpu *r)
1038 {
1039   if (sim_fpu_is_nan (l) || sim_fpu_is_nan (r))
1040     return sim_fpu_op_nan (f, l, r);
1041   if (sim_fpu_is_infinity (l))
1042     {
1043       if (sim_fpu_is_infinity (r)
1044           && l->sign != r->sign)
1045         {
1046           *f = sim_fpu_qnan;
1047           return sim_fpu_status_invalid_isi;
1048         }
1049       *f = *l;
1050       return 0;
1051     }
1052   if (sim_fpu_is_infinity (r))
1053     {
1054       *f = *r;
1055       return 0;
1056     }
1057   if (sim_fpu_is_zero (l))
1058     {
1059       if (sim_fpu_is_zero (r))
1060         {
1061           *f = sim_fpu_zero;
1062           f->sign = l->sign & r->sign;
1063         }
1064       else
1065         *f = *r;
1066       return 0;
1067     }
1068   if (sim_fpu_is_zero (r))
1069     {
1070       *f = *l;
1071       return 0;
1072     }
1073   {
1074     int status = 0;
1075     int shift = l->normal_exp - r->normal_exp;
1076     uint64_t lfraction;
1077     uint64_t rfraction;
1078     /* use exp of larger */
1079     if (shift >= NR_FRAC_GUARD)
1080       {
1081         /* left has much bigger magnitude */
1082         *f = *l;
1083         return sim_fpu_status_inexact;
1084       }
1085     if (shift <= - NR_FRAC_GUARD)
1086       {
1087         /* right has much bigger magnitude */
1088         *f = *r;
1089         return sim_fpu_status_inexact;
1090       }
1091     lfraction = l->fraction;
1092     rfraction = r->fraction;
1093     if (shift > 0)
1094       {
1095         f->normal_exp = l->normal_exp;
1096         if (rfraction & LSMASK64 (shift - 1, 0))
1097           {
1098             status |= sim_fpu_status_inexact;
1099             rfraction |= LSBIT64 (shift); /* Stick LSBit.  */
1100           }
1101         rfraction >>= shift;
1102       }
1103     else if (shift < 0)
1104       {
1105         f->normal_exp = r->normal_exp;
1106         if (lfraction & LSMASK64 (- shift - 1, 0))
1107           {
1108             status |= sim_fpu_status_inexact;
1109             lfraction |= LSBIT64 (- shift); /* Stick LSBit.  */
1110           }
1111         lfraction >>= -shift;
1112       }
1113     else
1114       {
1115         f->normal_exp = r->normal_exp;
1116       }
1117
1118     /* Perform the addition.  */
1119     if (l->sign)
1120       lfraction = - lfraction;
1121     if (r->sign)
1122       rfraction = - rfraction;
1123     f->fraction = lfraction + rfraction;
1124
1125     /* zero? */
1126     if (f->fraction == 0)
1127       {
1128         *f = sim_fpu_zero;
1129         return 0;
1130       }
1131
1132     /* sign? */
1133     f->class = sim_fpu_class_number;
1134     if (((int64_t) f->fraction) >= 0)
1135       f->sign = 0;
1136     else
1137       {
1138         f->sign = 1;
1139         f->fraction = - f->fraction;
1140       }
1141
1142     /* Normalize it.  */
1143     if ((f->fraction & IMPLICIT_2))
1144       {
1145         f->fraction = (f->fraction >> 1) | (f->fraction & 1);
1146         f->normal_exp ++;
1147       }
1148     else if (f->fraction < IMPLICIT_1)
1149       {
1150         do
1151           {
1152             f->fraction <<= 1;
1153             f->normal_exp --;
1154           }
1155         while (f->fraction < IMPLICIT_1);
1156       }
1157     ASSERT (f->fraction >= IMPLICIT_1 && f->fraction < IMPLICIT_2);
1158     return status;
1159   }
1160 }
1161
1162
1163 INLINE_SIM_FPU (int)
1164 sim_fpu_sub (sim_fpu *f,
1165              const sim_fpu *l,
1166              const sim_fpu *r)
1167 {
1168   if (sim_fpu_is_nan (l) || sim_fpu_is_nan (r))
1169     return sim_fpu_op_nan (f, l, r);
1170   if (sim_fpu_is_infinity (l))
1171     {
1172       if (sim_fpu_is_infinity (r)
1173           && l->sign == r->sign)
1174         {
1175           *f = sim_fpu_qnan;
1176           return sim_fpu_status_invalid_isi;
1177         }
1178       *f = *l;
1179       return 0;
1180     }
1181   if (sim_fpu_is_infinity (r))
1182     {
1183       *f = *r;
1184       f->sign = !r->sign;
1185       return 0;
1186     }
1187   if (sim_fpu_is_zero (l))
1188     {
1189       if (sim_fpu_is_zero (r))
1190         {
1191           *f = sim_fpu_zero;
1192           f->sign = l->sign & !r->sign;
1193         }
1194       else
1195         {
1196           *f = *r;
1197           f->sign = !r->sign;
1198         }
1199       return 0;
1200     }
1201   if (sim_fpu_is_zero (r))
1202     {
1203       *f = *l;
1204       return 0;
1205     }
1206   {
1207     int status = 0;
1208     int shift = l->normal_exp - r->normal_exp;
1209     uint64_t lfraction;
1210     uint64_t rfraction;
1211     /* use exp of larger */
1212     if (shift >= NR_FRAC_GUARD)
1213       {
1214         /* left has much bigger magnitude */
1215         *f = *l;
1216         return sim_fpu_status_inexact;
1217       }
1218     if (shift <= - NR_FRAC_GUARD)
1219       {
1220         /* right has much bigger magnitude */
1221         *f = *r;
1222         f->sign = !r->sign;
1223         return sim_fpu_status_inexact;
1224       }
1225     lfraction = l->fraction;
1226     rfraction = r->fraction;
1227     if (shift > 0)
1228       {
1229         f->normal_exp = l->normal_exp;
1230         if (rfraction & LSMASK64 (shift - 1, 0))
1231           {
1232             status |= sim_fpu_status_inexact;
1233             rfraction |= LSBIT64 (shift); /* Stick LSBit.  */
1234           }
1235         rfraction >>= shift;
1236       }
1237     else if (shift < 0)
1238       {
1239         f->normal_exp = r->normal_exp;
1240         if (lfraction & LSMASK64 (- shift - 1, 0))
1241           {
1242             status |= sim_fpu_status_inexact;
1243             lfraction |= LSBIT64 (- shift); /* Stick LSBit.  */
1244           }
1245         lfraction >>= -shift;
1246       }
1247     else
1248       {
1249         f->normal_exp = r->normal_exp;
1250       }
1251
1252     /* Perform the subtraction.  */
1253     if (l->sign)
1254       lfraction = - lfraction;
1255     if (!r->sign)
1256       rfraction = - rfraction;
1257     f->fraction = lfraction + rfraction;
1258
1259     /* zero? */
1260     if (f->fraction == 0)
1261       {
1262         *f = sim_fpu_zero;
1263         return 0;
1264       }
1265
1266     /* sign? */
1267     f->class = sim_fpu_class_number;
1268     if (((int64_t) f->fraction) >= 0)
1269       f->sign = 0;
1270     else
1271       {
1272         f->sign = 1;
1273         f->fraction = - f->fraction;
1274       }
1275
1276     /* Normalize it.  */
1277     if ((f->fraction & IMPLICIT_2))
1278       {
1279         f->fraction = (f->fraction >> 1) | (f->fraction & 1);
1280         f->normal_exp ++;
1281       }
1282     else if (f->fraction < IMPLICIT_1)
1283       {
1284         do
1285           {
1286             f->fraction <<= 1;
1287             f->normal_exp --;
1288           }
1289         while (f->fraction < IMPLICIT_1);
1290       }
1291     ASSERT (f->fraction >= IMPLICIT_1 && f->fraction < IMPLICIT_2);
1292     return status;
1293   }
1294 }
1295
1296
1297 INLINE_SIM_FPU (int)
1298 sim_fpu_mul (sim_fpu *f,
1299              const sim_fpu *l,
1300              const sim_fpu *r)
1301 {
1302   if (sim_fpu_is_nan (l) || sim_fpu_is_nan (r))
1303     return sim_fpu_op_nan (f, l, r);
1304   if (sim_fpu_is_infinity (l))
1305     {
1306       if (sim_fpu_is_zero (r))
1307         {
1308           *f = sim_fpu_qnan;
1309           return sim_fpu_status_invalid_imz;
1310         }
1311       *f = *l;
1312       f->sign = l->sign ^ r->sign;
1313       return 0;
1314     }
1315   if (sim_fpu_is_infinity (r))
1316     {
1317       if (sim_fpu_is_zero (l))
1318         {
1319           *f = sim_fpu_qnan;
1320           return sim_fpu_status_invalid_imz;
1321         }
1322       *f = *r;
1323       f->sign = l->sign ^ r->sign;
1324       return 0;
1325     }
1326   if (sim_fpu_is_zero (l) || sim_fpu_is_zero (r))
1327     {
1328       *f = sim_fpu_zero;
1329       f->sign = l->sign ^ r->sign;
1330       return 0;
1331     }
1332   /* Calculate the mantissa by multiplying both 64bit numbers to get a
1333      128 bit number.  */
1334   {
1335     uint64_t low;
1336     uint64_t high;
1337     uint64_t nl = l->fraction & 0xffffffff;
1338     uint64_t nh = l->fraction >> 32;
1339     uint64_t ml = r->fraction & 0xffffffff;
1340     uint64_t mh = r->fraction >>32;
1341     uint64_t pp_ll = ml * nl;
1342     uint64_t pp_hl = mh * nl;
1343     uint64_t pp_lh = ml * nh;
1344     uint64_t pp_hh = mh * nh;
1345     uint64_t res2 = 0;
1346     uint64_t res0 = 0;
1347     uint64_t ps_hh__ = pp_hl + pp_lh;
1348     if (ps_hh__ < pp_hl)
1349       res2 += UNSIGNED64 (0x100000000);
1350     pp_hl = (ps_hh__ << 32) & UNSIGNED64 (0xffffffff00000000);
1351     res0 = pp_ll + pp_hl;
1352     if (res0 < pp_ll)
1353       res2++;
1354     res2 += ((ps_hh__ >> 32) & 0xffffffff) + pp_hh;
1355     high = res2;
1356     low = res0;
1357
1358     f->normal_exp = l->normal_exp + r->normal_exp;
1359     f->sign = l->sign ^ r->sign;
1360     f->class = sim_fpu_class_number;
1361
1362     /* Input is bounded by [1,2)   ;   [2^60,2^61)
1363        Output is bounded by [1,4)  ;   [2^120,2^122) */
1364
1365     /* Adjust the exponent according to where the decimal point ended
1366        up in the high 64 bit word.  In the source the decimal point
1367        was at NR_FRAC_GUARD. */
1368     f->normal_exp += NR_FRAC_GUARD + 64 - (NR_FRAC_GUARD * 2);
1369
1370     /* The high word is bounded according to the above.  Consequently
1371        it has never overflowed into IMPLICIT_2. */
1372     ASSERT (high < LSBIT64 (((NR_FRAC_GUARD + 1) * 2) - 64));
1373     ASSERT (high >= LSBIT64 ((NR_FRAC_GUARD * 2) - 64));
1374     ASSERT (LSBIT64 (((NR_FRAC_GUARD + 1) * 2) - 64) < IMPLICIT_1);
1375
1376     /* Normalize.  */
1377     do
1378       {
1379         f->normal_exp--;
1380         high <<= 1;
1381         if (low & LSBIT64 (63))
1382           high |= 1;
1383         low <<= 1;
1384       }
1385     while (high < IMPLICIT_1);
1386
1387     ASSERT (high >= IMPLICIT_1 && high < IMPLICIT_2);
1388     if (low != 0)
1389       {
1390         f->fraction = (high | 1); /* sticky */
1391         return sim_fpu_status_inexact;
1392       }
1393     else
1394       {
1395         f->fraction = high;
1396         return 0;
1397       }
1398     return 0;
1399   }
1400 }
1401
1402 INLINE_SIM_FPU (int)
1403 sim_fpu_div (sim_fpu *f,
1404              const sim_fpu *l,
1405              const sim_fpu *r)
1406 {
1407   if (sim_fpu_is_nan (l) || sim_fpu_is_nan (r))
1408     return sim_fpu_op_nan (f, l, r);
1409   if (sim_fpu_is_infinity (l))
1410     {
1411       if (sim_fpu_is_infinity (r))
1412         {
1413           *f = sim_fpu_qnan;
1414           return sim_fpu_status_invalid_idi;
1415         }
1416       else
1417         {
1418           *f = *l;
1419           f->sign = l->sign ^ r->sign;
1420           return 0;
1421         }
1422     }
1423   if (sim_fpu_is_zero (l))
1424     {
1425       if (sim_fpu_is_zero (r))
1426         {
1427           *f = sim_fpu_qnan;
1428           return sim_fpu_status_invalid_zdz;
1429         }
1430       else
1431         {
1432           *f = *l;
1433           f->sign = l->sign ^ r->sign;
1434           return 0;
1435         }
1436     }
1437   if (sim_fpu_is_infinity (r))
1438     {
1439       *f = sim_fpu_zero;
1440       f->sign = l->sign ^ r->sign;
1441       return 0;
1442     }
1443   if (sim_fpu_is_zero (r))
1444     {
1445       f->class = sim_fpu_class_infinity;
1446       f->sign = l->sign ^ r->sign;
1447       return sim_fpu_status_invalid_div0;
1448     }
1449
1450   /* Calculate the mantissa by multiplying both 64bit numbers to get a
1451      128 bit number.  */
1452   {
1453     /* quotient =  ( ( numerator / denominator)
1454                       x 2^(numerator exponent -  denominator exponent)
1455      */
1456     uint64_t numerator;
1457     uint64_t denominator;
1458     uint64_t quotient;
1459     uint64_t bit;
1460
1461     f->class = sim_fpu_class_number;
1462     f->sign = l->sign ^ r->sign;
1463     f->normal_exp = l->normal_exp - r->normal_exp;
1464
1465     numerator = l->fraction;
1466     denominator = r->fraction;
1467
1468     /* Fraction will be less than 1.0 */
1469     if (numerator < denominator)
1470       {
1471         numerator <<= 1;
1472         f->normal_exp--;
1473       }
1474     ASSERT (numerator >= denominator);
1475
1476     /* Gain extra precision, already used one spare bit.  */
1477     numerator <<=    NR_SPARE;
1478     denominator <<=  NR_SPARE;
1479
1480     /* Does divide one bit at a time.  Optimize???  */
1481     quotient = 0;
1482     bit = (IMPLICIT_1 << NR_SPARE);
1483     while (bit)
1484       {
1485         if (numerator >= denominator)
1486           {
1487             quotient |= bit;
1488             numerator -= denominator;
1489           }
1490         bit >>= 1;
1491         numerator <<= 1;
1492       }
1493
1494     /* Discard (but save) the extra bits.  */
1495     if ((quotient & LSMASK64 (NR_SPARE -1, 0)))
1496       quotient = (quotient >> NR_SPARE) | 1;
1497     else
1498       quotient = (quotient >> NR_SPARE);
1499
1500     f->fraction = quotient;
1501     ASSERT (f->fraction >= IMPLICIT_1 && f->fraction < IMPLICIT_2);
1502     if (numerator != 0)
1503       {
1504         f->fraction |= 1; /* Stick remaining bits.  */
1505         return sim_fpu_status_inexact;
1506       }
1507     else
1508       return 0;
1509   }
1510 }
1511
1512
1513 INLINE_SIM_FPU (int)
1514 sim_fpu_rem (sim_fpu *f,
1515              const sim_fpu *l,
1516              const sim_fpu *r)
1517 {
1518   if (sim_fpu_is_nan (l) || sim_fpu_is_nan (r))
1519     return sim_fpu_op_nan (f, l, r);
1520   if (sim_fpu_is_infinity (l))
1521     {
1522       *f = sim_fpu_qnan;
1523       return sim_fpu_status_invalid_irx;
1524     }
1525   if (sim_fpu_is_zero (r))
1526     {
1527       *f = sim_fpu_qnan;
1528       return sim_fpu_status_invalid_div0;
1529     }
1530   if (sim_fpu_is_zero (l))
1531     {
1532       *f = *l;
1533       return 0;
1534     }
1535   if (sim_fpu_is_infinity (r))
1536     {
1537       *f = *l;
1538       return 0;
1539     }
1540   {
1541     sim_fpu n, tmp;
1542
1543     /* Remainder is calculated as l-n*r, where n is l/r rounded to the
1544        nearest integer.  The variable n is rounded half even.  */
1545
1546     sim_fpu_div (&n, l, r);
1547     sim_fpu_round_64 (&n, 0, 0);
1548
1549     if (n.normal_exp < -1) /* If n looks like zero just return l.  */
1550       {
1551         *f = *l;
1552         return 0;
1553       }
1554     else if (n.class == sim_fpu_class_number
1555              && n.normal_exp <= (NR_FRAC_GUARD)) /* If not too large round.  */
1556       do_normal_round (&n, (NR_FRAC_GUARD) - n.normal_exp, sim_fpu_round_near);
1557
1558     /* Mark 0's as zero so multiply can detect zero.  */
1559     if (n.fraction == 0)
1560       n.class = sim_fpu_class_zero;
1561
1562     /* Calculate n*r.  */
1563     sim_fpu_mul (&tmp, &n, r);
1564     sim_fpu_round_64 (&tmp, 0, 0);
1565
1566     /* Finally calculate l-n*r.  */
1567     sim_fpu_sub (f, l, &tmp);
1568
1569     return 0;
1570   }
1571 }
1572
1573
1574 INLINE_SIM_FPU (int)
1575 sim_fpu_max (sim_fpu *f,
1576              const sim_fpu *l,
1577              const sim_fpu *r)
1578 {
1579   if (sim_fpu_is_nan (l) || sim_fpu_is_nan (r))
1580     return sim_fpu_minmax_nan (f, l, r);
1581   if (sim_fpu_is_infinity (l))
1582     {
1583       if (sim_fpu_is_infinity (r)
1584           && l->sign == r->sign)
1585         {
1586           *f = sim_fpu_qnan;
1587           return sim_fpu_status_invalid_isi;
1588         }
1589       if (l->sign)
1590         *f = *r; /* -inf < anything */
1591       else
1592         *f = *l; /* +inf > anything */
1593       return 0;
1594     }
1595   if (sim_fpu_is_infinity (r))
1596     {
1597       if (r->sign)
1598         *f = *l; /* anything > -inf */
1599       else
1600         *f = *r; /* anything < +inf */
1601       return 0;
1602     }
1603   if (l->sign > r->sign)
1604     {
1605       *f = *r; /* -ve < +ve */
1606       return 0;
1607     }
1608   if (l->sign < r->sign)
1609     {
1610       *f = *l; /* +ve > -ve */
1611       return 0;
1612     }
1613   ASSERT (l->sign == r->sign);
1614   if (l->normal_exp > r->normal_exp
1615       || (l->normal_exp == r->normal_exp
1616           && l->fraction > r->fraction))
1617     {
1618       /* |l| > |r| */
1619       if (l->sign)
1620         *f = *r; /* -ve < -ve */
1621       else
1622         *f = *l; /* +ve > +ve */
1623       return 0;
1624     }
1625   else
1626     {
1627       /* |l| <= |r| */
1628       if (l->sign)
1629         *f = *l; /* -ve > -ve */
1630       else
1631         *f = *r; /* +ve < +ve */
1632       return 0;
1633     }
1634 }
1635
1636
1637 INLINE_SIM_FPU (int)
1638 sim_fpu_min (sim_fpu *f,
1639              const sim_fpu *l,
1640              const sim_fpu *r)
1641 {
1642   if (sim_fpu_is_nan (l) || sim_fpu_is_nan (r))
1643     return sim_fpu_minmax_nan (f, l, r);
1644   if (sim_fpu_is_infinity (l))
1645     {
1646       if (sim_fpu_is_infinity (r)
1647           && l->sign == r->sign)
1648         {
1649           *f = sim_fpu_qnan;
1650           return sim_fpu_status_invalid_isi;
1651         }
1652       if (l->sign)
1653         *f = *l; /* -inf < anything */
1654       else
1655         *f = *r; /* +inf > anthing */
1656       return 0;
1657     }
1658   if (sim_fpu_is_infinity (r))
1659     {
1660       if (r->sign)
1661         *f = *r; /* anything > -inf */
1662       else
1663         *f = *l; /* anything < +inf */
1664       return 0;
1665     }
1666   if (l->sign > r->sign)
1667     {
1668       *f = *l; /* -ve < +ve */
1669       return 0;
1670     }
1671   if (l->sign < r->sign)
1672     {
1673       *f = *r; /* +ve > -ve */
1674       return 0;
1675     }
1676   ASSERT (l->sign == r->sign);
1677   if (l->normal_exp > r->normal_exp
1678       || (l->normal_exp == r->normal_exp
1679           && l->fraction > r->fraction))
1680     {
1681       /* |l| > |r| */
1682       if (l->sign)
1683         *f = *l; /* -ve < -ve */
1684       else
1685         *f = *r; /* +ve > +ve */
1686       return 0;
1687     }
1688   else
1689     {
1690       /* |l| <= |r| */
1691       if (l->sign)
1692         *f = *r; /* -ve > -ve */
1693       else
1694         *f = *l; /* +ve < +ve */
1695       return 0;
1696     }
1697 }
1698
1699
1700 INLINE_SIM_FPU (int)
1701 sim_fpu_neg (sim_fpu *f,
1702              const sim_fpu *r)
1703 {
1704   if (sim_fpu_is_ieee754_1985 () && sim_fpu_is_snan (r))
1705     {
1706       *f = *r;
1707       f->class = sim_fpu_class_qnan;
1708       return sim_fpu_status_invalid_snan;
1709     }
1710   if (sim_fpu_is_qnan (r))
1711     {
1712       *f = *r;
1713       return 0;
1714     }
1715   *f = *r;
1716   f->sign = !r->sign;
1717   return 0;
1718 }
1719
1720
1721 INLINE_SIM_FPU (int)
1722 sim_fpu_abs (sim_fpu *f,
1723              const sim_fpu *r)
1724 {
1725   *f = *r;
1726   f->sign = 0;
1727   if (sim_fpu_is_ieee754_1985 () && sim_fpu_is_snan (r))
1728     {
1729       f->class = sim_fpu_class_qnan;
1730       return sim_fpu_status_invalid_snan;
1731     }
1732   return 0;
1733 }
1734
1735
1736 INLINE_SIM_FPU (int)
1737 sim_fpu_inv (sim_fpu *f,
1738              const sim_fpu *r)
1739 {
1740   return sim_fpu_div (f, &sim_fpu_one, r);
1741 }
1742
1743
1744 INLINE_SIM_FPU (int)
1745 sim_fpu_sqrt (sim_fpu *f,
1746               const sim_fpu *r)
1747 {
1748   if (sim_fpu_is_snan (r))
1749     {
1750       *f = sim_fpu_qnan;
1751       return sim_fpu_status_invalid_snan;
1752     }
1753   if (sim_fpu_is_qnan (r))
1754     {
1755       *f = sim_fpu_qnan;
1756       return 0;
1757     }
1758   if (sim_fpu_is_zero (r))
1759     {
1760       f->class = sim_fpu_class_zero;
1761       f->sign = r->sign;
1762       f->normal_exp = 0;
1763       return 0;
1764     }
1765   if (sim_fpu_is_infinity (r))
1766     {
1767       if (r->sign)
1768         {
1769           *f = sim_fpu_qnan;
1770           return sim_fpu_status_invalid_sqrt;
1771         }
1772       else
1773         {
1774           f->class = sim_fpu_class_infinity;
1775           f->sign = 0;
1776           f->sign = 0;
1777           return 0;
1778         }
1779     }
1780   if (r->sign)
1781     {
1782       *f = sim_fpu_qnan;
1783       return sim_fpu_status_invalid_sqrt;
1784     }
1785
1786   /* @(#)e_sqrt.c 5.1 93/09/24 */
1787   /*
1788    * ====================================================
1789    * Copyright (C) 1993 by Sun Microsystems, Inc. All rights reserved.
1790    *
1791    * Developed at SunPro, a Sun Microsystems, Inc. business.
1792    * Permission to use, copy, modify, and distribute this
1793    * software is freely granted, provided that this notice
1794    * is preserved.
1795    * ====================================================
1796    */
1797
1798   /* __ieee754_sqrt(x)
1799    * Return correctly rounded sqrt.
1800    *           ------------------------------------------
1801    *           |  Use the hardware sqrt if you have one |
1802    *           ------------------------------------------
1803    * Method:
1804    *   Bit by bit method using integer arithmetic. (Slow, but portable)
1805    *   1. Normalization
1806    *    Scale x to y in [1,4) with even powers of 2:
1807    *    find an integer k such that  1 <= (y=x*2^(2k)) < 4, then
1808    *            sqrt(x) = 2^k * sqrt(y)
1809    -
1810    - Since:
1811    -   sqrt ( x*2^(2m) )     = sqrt(x).2^m    ; m even
1812    -   sqrt ( x*2^(2m + 1) ) = sqrt(2.x).2^m  ; m odd
1813    - Define:
1814    -   y = ((m even) ? x : 2.x)
1815    - Then:
1816    -   y in [1, 4)                            ; [IMPLICIT_1,IMPLICIT_4)
1817    - And:
1818    -   sqrt (y) in [1, 2)                     ; [IMPLICIT_1,IMPLICIT_2)
1819    -
1820    *   2. Bit by bit computation
1821    *    Let q  = sqrt(y) truncated to i bit after binary point (q = 1),
1822    *         i                                                   0
1823    *                                     i+1         2
1824    *        s  = 2*q , and      y  =  2   * ( y - q  ).         (1)
1825    *         i      i            i                 i
1826    *
1827    *    To compute q    from q , one checks whether
1828    *                i+1       i
1829    *
1830    *                          -(i+1) 2
1831    *                    (q + 2      ) <= y.                     (2)
1832    *                              i
1833    *                                                          -(i+1)
1834    *    If (2) is false, then q   = q ; otherwise q   = q  + 2      .
1835    *                           i+1   i             i+1   i
1836    *
1837    *    With some algebraic manipulation, it is not difficult to see
1838    *    that (2) is equivalent to
1839    *                             -(i+1)
1840    *                    s  +  2       <= y                      (3)
1841    *                     i                i
1842    *
1843    *    The advantage of (3) is that s  and y  can be computed by
1844    *                                  i      i
1845    *    the following recurrence formula:
1846    *        if (3) is false
1847    *
1848    *        s     =  s  ,       y    = y   ;                    (4)
1849    *         i+1      i          i+1    i
1850    *
1851    -
1852    -                      NOTE: y    = 2*y
1853    -                             i+1      i
1854    -
1855    *        otherwise,
1856    *                       -i                      -(i+1)
1857    *        s     =  s  + 2  ,  y    = y  -  s  - 2             (5)
1858    *         i+1      i          i+1    i     i
1859    *
1860    -
1861    -                                                   -(i+1)
1862    -                      NOTE: y    = 2 (y  -  s  -  2      )
1863    -                             i+1       i     i
1864    -
1865    *    One may easily use induction to prove (4) and (5).
1866    *    Note. Since the left hand side of (3) contain only i+2 bits,
1867    *          it does not necessary to do a full (53-bit) comparison
1868    *          in (3).
1869    *   3. Final rounding
1870    *    After generating the 53 bits result, we compute one more bit.
1871    *    Together with the remainder, we can decide whether the
1872    *    result is exact, bigger than 1/2ulp, or less than 1/2ulp
1873    *    (it will never equal to 1/2ulp).
1874    *    The rounding mode can be detected by checking whether
1875    *    huge + tiny is equal to huge, and whether huge - tiny is
1876    *    equal to huge for some floating point number "huge" and "tiny".
1877    *
1878    * Special cases:
1879    *    sqrt(+-0) = +-0         ... exact
1880    *    sqrt(inf) = inf
1881    *    sqrt(-ve) = NaN         ... with invalid signal
1882    *    sqrt(NaN) = NaN         ... with invalid signal for signalling NaN
1883    *
1884    * Other methods : see the appended file at the end of the program below.
1885    *---------------
1886    */
1887
1888   {
1889     /* Generate sqrt(x) bit by bit.  */
1890     uint64_t y;
1891     uint64_t q;
1892     uint64_t s;
1893     uint64_t b;
1894
1895     f->class = sim_fpu_class_number;
1896     f->sign = 0;
1897     y = r->fraction;
1898     f->normal_exp = (r->normal_exp >> 1);       /* exp = [exp/2] */
1899
1900     /* Odd exp, double x to make it even.  */
1901     ASSERT (y >= IMPLICIT_1 && y < IMPLICIT_4);
1902     if ((r->normal_exp & 1))
1903       {
1904         y += y;
1905       }
1906     ASSERT (y >= IMPLICIT_1 && y < (IMPLICIT_2 << 1));
1907
1908     /* Let loop determine first value of s (either 1 or 2) */
1909     b = IMPLICIT_1;
1910     q = 0;
1911     s = 0;
1912
1913     while (b)
1914       {
1915         uint64_t t = s + b;
1916         if (t <= y)
1917           {
1918             s |= (b << 1);
1919             y -= t;
1920             q |= b;
1921           }
1922         y <<= 1;
1923         b >>= 1;
1924       }
1925
1926     ASSERT (q >= IMPLICIT_1 && q < IMPLICIT_2);
1927     f->fraction = q;
1928     if (y != 0)
1929       {
1930         f->fraction |= 1; /* Stick remaining bits.  */
1931         return sim_fpu_status_inexact;
1932       }
1933     else
1934       return 0;
1935   }
1936 }
1937
1938
1939 /* int/long <-> sim_fpu */
1940
1941 INLINE_SIM_FPU (int)
1942 sim_fpu_i32to (sim_fpu *f,
1943                int32_t i,
1944                sim_fpu_round round)
1945 {
1946   i2fpu (f, i, 0);
1947   return 0;
1948 }
1949
1950 INLINE_SIM_FPU (int)
1951 sim_fpu_u32to (sim_fpu *f,
1952                uint32_t u,
1953                sim_fpu_round round)
1954 {
1955   u2fpu (f, u, 0);
1956   return 0;
1957 }
1958
1959 INLINE_SIM_FPU (int)
1960 sim_fpu_i64to (sim_fpu *f,
1961                int64_t i,
1962                sim_fpu_round round)
1963 {
1964   i2fpu (f, i, 1);
1965   return 0;
1966 }
1967
1968 INLINE_SIM_FPU (int)
1969 sim_fpu_u64to (sim_fpu *f,
1970                uint64_t u,
1971                sim_fpu_round round)
1972 {
1973   u2fpu (f, u, 1);
1974   return 0;
1975 }
1976
1977
1978 INLINE_SIM_FPU (int)
1979 sim_fpu_to32i (int32_t *i,
1980                const sim_fpu *f,
1981                sim_fpu_round round)
1982 {
1983   int64_t i64;
1984   int status = fpu2i (&i64, f, 0, round);
1985   *i = i64;
1986   return status;
1987 }
1988
1989 INLINE_SIM_FPU (int)
1990 sim_fpu_to32u (uint32_t *u,
1991                const sim_fpu *f,
1992                sim_fpu_round round)
1993 {
1994   uint64_t u64;
1995   int status = fpu2u (&u64, f, 0);
1996   *u = u64;
1997   return status;
1998 }
1999
2000 INLINE_SIM_FPU (int)
2001 sim_fpu_to64i (int64_t *i,
2002                const sim_fpu *f,
2003                sim_fpu_round round)
2004 {
2005   return fpu2i (i, f, 1, round);
2006 }
2007
2008
2009 INLINE_SIM_FPU (int)
2010 sim_fpu_to64u (uint64_t *u,
2011                const sim_fpu *f,
2012                sim_fpu_round round)
2013 {
2014   return fpu2u (u, f, 1);
2015 }
2016
2017
2018
2019 /* sim_fpu -> host format */
2020
2021 #if 0
2022 INLINE_SIM_FPU (float)
2023 sim_fpu_2f (const sim_fpu *f)
2024 {
2025   return fval.d;
2026 }
2027 #endif
2028
2029
2030 INLINE_SIM_FPU (double)
2031 sim_fpu_2d (const sim_fpu *s)
2032 {
2033   sim_fpu_map val;
2034   if (sim_fpu_is_snan (s))
2035     {
2036       /* gag SNaN's */
2037       sim_fpu n = *s;
2038       n.class = sim_fpu_class_qnan;
2039       val.i = pack_fpu (&n, 1);
2040     }
2041   else
2042     {
2043       val.i = pack_fpu (s, 1);
2044     }
2045   return val.d;
2046 }
2047
2048
2049 #if 0
2050 INLINE_SIM_FPU (void)
2051 sim_fpu_f2 (sim_fpu *f,
2052             float s)
2053 {
2054   sim_fpu_map val;
2055   val.d = s;
2056   unpack_fpu (f, val.i, 1);
2057 }
2058 #endif
2059
2060
2061 INLINE_SIM_FPU (void)
2062 sim_fpu_d2 (sim_fpu *f,
2063             double d)
2064 {
2065   sim_fpu_map val;
2066   val.d = d;
2067   unpack_fpu (f, val.i, 1);
2068 }
2069
2070
2071 /* General */
2072
2073 INLINE_SIM_FPU (int)
2074 sim_fpu_is_nan (const sim_fpu *d)
2075 {
2076   switch (d->class)
2077     {
2078     case sim_fpu_class_qnan:
2079     case sim_fpu_class_snan:
2080       return 1;
2081     default:
2082       return 0;
2083     }
2084 }
2085
2086 INLINE_SIM_FPU (int)
2087 sim_fpu_is_qnan (const sim_fpu *d)
2088 {
2089   switch (d->class)
2090     {
2091     case sim_fpu_class_qnan:
2092       return 1;
2093     default:
2094       return 0;
2095     }
2096 }
2097
2098 INLINE_SIM_FPU (int)
2099 sim_fpu_is_snan (const sim_fpu *d)
2100 {
2101   switch (d->class)
2102     {
2103     case sim_fpu_class_snan:
2104       return 1;
2105     default:
2106       return 0;
2107     }
2108 }
2109
2110 INLINE_SIM_FPU (int)
2111 sim_fpu_is_zero (const sim_fpu *d)
2112 {
2113   switch (d->class)
2114     {
2115     case sim_fpu_class_zero:
2116       return 1;
2117     default:
2118       return 0;
2119     }
2120 }
2121
2122 INLINE_SIM_FPU (int)
2123 sim_fpu_is_infinity (const sim_fpu *d)
2124 {
2125   switch (d->class)
2126     {
2127     case sim_fpu_class_infinity:
2128       return 1;
2129     default:
2130       return 0;
2131     }
2132 }
2133
2134 INLINE_SIM_FPU (int)
2135 sim_fpu_is_number (const sim_fpu *d)
2136 {
2137   switch (d->class)
2138     {
2139     case sim_fpu_class_denorm:
2140     case sim_fpu_class_number:
2141       return 1;
2142     default:
2143       return 0;
2144     }
2145 }
2146
2147 INLINE_SIM_FPU (int)
2148 sim_fpu_is_denorm (const sim_fpu *d)
2149 {
2150   switch (d->class)
2151     {
2152     case sim_fpu_class_denorm:
2153       return 1;
2154     default:
2155       return 0;
2156     }
2157 }
2158
2159
2160 INLINE_SIM_FPU (int)
2161 sim_fpu_sign (const sim_fpu *d)
2162 {
2163   return d->sign;
2164 }
2165
2166
2167 INLINE_SIM_FPU (int)
2168 sim_fpu_exp (const sim_fpu *d)
2169 {
2170   return d->normal_exp;
2171 }
2172
2173
2174 INLINE_SIM_FPU (uint64_t)
2175 sim_fpu_fraction (const sim_fpu *d)
2176 {
2177   return d->fraction;
2178 }
2179
2180
2181 INLINE_SIM_FPU (uint64_t)
2182 sim_fpu_guard (const sim_fpu *d, int is_double)
2183 {
2184   uint64_t rv;
2185   uint64_t guardmask = LSMASK64 (NR_GUARDS - 1, 0);
2186   rv = (d->fraction & guardmask) >> NR_PAD;
2187   return rv;
2188 }
2189
2190
2191 INLINE_SIM_FPU (int)
2192 sim_fpu_is (const sim_fpu *d)
2193 {
2194   switch (d->class)
2195     {
2196     case sim_fpu_class_qnan:
2197       return SIM_FPU_IS_QNAN;
2198     case sim_fpu_class_snan:
2199       return SIM_FPU_IS_SNAN;
2200     case sim_fpu_class_infinity:
2201       if (d->sign)
2202         return SIM_FPU_IS_NINF;
2203       else
2204         return SIM_FPU_IS_PINF;
2205     case sim_fpu_class_number:
2206       if (d->sign)
2207         return SIM_FPU_IS_NNUMBER;
2208       else
2209         return SIM_FPU_IS_PNUMBER;
2210     case sim_fpu_class_denorm:
2211       if (d->sign)
2212         return SIM_FPU_IS_NDENORM;
2213       else
2214         return SIM_FPU_IS_PDENORM;
2215     case sim_fpu_class_zero:
2216       if (d->sign)
2217         return SIM_FPU_IS_NZERO;
2218       else
2219         return SIM_FPU_IS_PZERO;
2220     default:
2221       return -1;
2222       abort ();
2223     }
2224 }
2225
2226 INLINE_SIM_FPU (int)
2227 sim_fpu_cmp (const sim_fpu *l, const sim_fpu *r)
2228 {
2229   sim_fpu res;
2230   sim_fpu_sub (&res, l, r);
2231   return sim_fpu_is (&res);
2232 }
2233
2234 INLINE_SIM_FPU (int)
2235 sim_fpu_is_lt (const sim_fpu *l, const sim_fpu *r)
2236 {
2237   int status;
2238   sim_fpu_lt (&status, l, r);
2239   return status;
2240 }
2241
2242 INLINE_SIM_FPU (int)
2243 sim_fpu_is_le (const sim_fpu *l, const sim_fpu *r)
2244 {
2245   int is;
2246   sim_fpu_le (&is, l, r);
2247   return is;
2248 }
2249
2250 INLINE_SIM_FPU (int)
2251 sim_fpu_is_eq (const sim_fpu *l, const sim_fpu *r)
2252 {
2253   int is;
2254   sim_fpu_eq (&is, l, r);
2255   return is;
2256 }
2257
2258 INLINE_SIM_FPU (int)
2259 sim_fpu_is_ne (const sim_fpu *l, const sim_fpu *r)
2260 {
2261   int is;
2262   sim_fpu_ne (&is, l, r);
2263   return is;
2264 }
2265
2266 INLINE_SIM_FPU (int)
2267 sim_fpu_is_ge (const sim_fpu *l, const sim_fpu *r)
2268 {
2269   int is;
2270   sim_fpu_ge (&is, l, r);
2271   return is;
2272 }
2273
2274 INLINE_SIM_FPU (int)
2275 sim_fpu_is_gt (const sim_fpu *l, const sim_fpu *r)
2276 {
2277   int is;
2278   sim_fpu_gt (&is, l, r);
2279   return is;
2280 }
2281
2282 INLINE_SIM_FPU (int)
2283 sim_fpu_is_un (const sim_fpu *l, const sim_fpu *r)
2284 {
2285   int is;
2286   sim_fpu_un (&is, l, r);
2287   return is;
2288 }
2289
2290 INLINE_SIM_FPU (int)
2291 sim_fpu_is_or (const sim_fpu *l, const sim_fpu *r)
2292 {
2293   int is;
2294   sim_fpu_or (&is, l, r);
2295   return is;
2296 }
2297
2298 /* Compare operators */
2299
2300 INLINE_SIM_FPU (int)
2301 sim_fpu_lt (int *is,
2302             const sim_fpu *l,
2303             const sim_fpu *r)
2304 {
2305   if (!sim_fpu_is_nan (l) && !sim_fpu_is_nan (r))
2306     {
2307       sim_fpu_map lval;
2308       sim_fpu_map rval;
2309       lval.i = pack_fpu (l, 1);
2310       rval.i = pack_fpu (r, 1);
2311       (*is) = (lval.d < rval.d);
2312       return 0;
2313     }
2314   else if (sim_fpu_is_snan (l) || sim_fpu_is_snan (r))
2315     {
2316       *is = 0;
2317       return sim_fpu_status_invalid_snan;
2318     }
2319   else
2320     {
2321       *is = 0;
2322       return sim_fpu_status_invalid_qnan;
2323     }
2324 }
2325
2326 INLINE_SIM_FPU (int)
2327 sim_fpu_le (int *is,
2328             const sim_fpu *l,
2329             const sim_fpu *r)
2330 {
2331   if (!sim_fpu_is_nan (l) && !sim_fpu_is_nan (r))
2332     {
2333       sim_fpu_map lval;
2334       sim_fpu_map rval;
2335       lval.i = pack_fpu (l, 1);
2336       rval.i = pack_fpu (r, 1);
2337       *is = (lval.d <= rval.d);
2338       return 0;
2339     }
2340   else if (sim_fpu_is_snan (l) || sim_fpu_is_snan (r))
2341     {
2342       *is = 0;
2343       return sim_fpu_status_invalid_snan;
2344     }
2345   else
2346     {
2347       *is = 0;
2348       return sim_fpu_status_invalid_qnan;
2349     }
2350 }
2351
2352 INLINE_SIM_FPU (int)
2353 sim_fpu_eq (int *is,
2354             const sim_fpu *l,
2355             const sim_fpu *r)
2356 {
2357   if (!sim_fpu_is_nan (l) && !sim_fpu_is_nan (r))
2358     {
2359       sim_fpu_map lval;
2360       sim_fpu_map rval;
2361       lval.i = pack_fpu (l, 1);
2362       rval.i = pack_fpu (r, 1);
2363       (*is) = (lval.d == rval.d);
2364       return 0;
2365     }
2366   else if (sim_fpu_is_snan (l) || sim_fpu_is_snan (r))
2367     {
2368       *is = 0;
2369       return sim_fpu_status_invalid_snan;
2370     }
2371   else
2372     {
2373       *is = 0;
2374       return sim_fpu_status_invalid_qnan;
2375     }
2376 }
2377
2378 INLINE_SIM_FPU (int)
2379 sim_fpu_ne (int *is,
2380             const sim_fpu *l,
2381             const sim_fpu *r)
2382 {
2383   if (!sim_fpu_is_nan (l) && !sim_fpu_is_nan (r))
2384     {
2385       sim_fpu_map lval;
2386       sim_fpu_map rval;
2387       lval.i = pack_fpu (l, 1);
2388       rval.i = pack_fpu (r, 1);
2389       (*is) = (lval.d != rval.d);
2390       return 0;
2391     }
2392   else if (sim_fpu_is_snan (l) || sim_fpu_is_snan (r))
2393     {
2394       *is = 0;
2395       return sim_fpu_status_invalid_snan;
2396     }
2397   else
2398     {
2399       *is = 0;
2400       return sim_fpu_status_invalid_qnan;
2401     }
2402 }
2403
2404 INLINE_SIM_FPU (int)
2405 sim_fpu_ge (int *is,
2406             const sim_fpu *l,
2407             const sim_fpu *r)
2408 {
2409   return sim_fpu_le (is, r, l);
2410 }
2411
2412 INLINE_SIM_FPU (int)
2413 sim_fpu_gt (int *is,
2414             const sim_fpu *l,
2415             const sim_fpu *r)
2416 {
2417   return sim_fpu_lt (is, r, l);
2418 }
2419
2420 INLINE_SIM_FPU (int)
2421 sim_fpu_un (int *is, const sim_fpu *l, const sim_fpu *r)
2422 {
2423   if (sim_fpu_is_nan (l) || sim_fpu_is_nan (r))
2424    {
2425     *is = 1;
2426     return 0;
2427    }
2428
2429   *is = 0;
2430   return 0;
2431 }
2432
2433 INLINE_SIM_FPU (int)
2434 sim_fpu_or (int *is, const sim_fpu *l, const sim_fpu *r)
2435 {
2436   sim_fpu_un (is, l, r);
2437
2438   /* Invert result.  */
2439   *is = !*is;
2440   return 0;
2441 }
2442
2443 INLINE_SIM_FPU(int)
2444 sim_fpu_classify (const sim_fpu *f)
2445 {
2446   switch (f->class)
2447     {
2448     case sim_fpu_class_snan: return SIM_FPU_IS_SNAN;
2449     case sim_fpu_class_qnan: return SIM_FPU_IS_QNAN;
2450     case sim_fpu_class_infinity:
2451       return f->sign ? SIM_FPU_IS_NINF : SIM_FPU_IS_PINF;
2452     case sim_fpu_class_zero:
2453       return f->sign ? SIM_FPU_IS_NZERO : SIM_FPU_IS_PZERO;
2454     case sim_fpu_class_number:
2455       return f->sign ? SIM_FPU_IS_NNUMBER : SIM_FPU_IS_PNUMBER;
2456     case sim_fpu_class_denorm:
2457       return f->sign ? SIM_FPU_IS_NDENORM : SIM_FPU_IS_PDENORM;
2458     default:
2459       fprintf (stderr, "Bad switch\n");
2460       abort ();
2461     }
2462   return 0;
2463 }
2464
2465 /* A number of useful constants */
2466
2467 #if EXTERN_SIM_FPU_P
2468 sim_fpu_state _sim_fpu = {
2469   .quiet_nan_inverted = false,
2470   .current_mode = sim_fpu_ieee754_1985,
2471 };
2472
2473 const sim_fpu sim_fpu_zero = {
2474   sim_fpu_class_zero, 0, 0, 0
2475 };
2476 const sim_fpu sim_fpu_qnan = {
2477   sim_fpu_class_qnan, 0, 0, 0
2478 };
2479 const sim_fpu sim_fpu_one = {
2480   sim_fpu_class_number, 0, IMPLICIT_1, 0
2481 };
2482 const sim_fpu sim_fpu_two = {
2483   sim_fpu_class_number, 0, IMPLICIT_1, 1
2484 };
2485 const sim_fpu sim_fpu_max32 = {
2486   sim_fpu_class_number, 0, LSMASK64 (NR_FRAC_GUARD, NR_GUARDS32), NORMAL_EXPMAX32
2487 };
2488 const sim_fpu sim_fpu_max64 = {
2489   sim_fpu_class_number, 0, LSMASK64 (NR_FRAC_GUARD, NR_GUARDS64), NORMAL_EXPMAX64
2490 };
2491 #endif
2492
2493 /* Specification swapping behaviour */
2494 INLINE_SIM_FPU (bool)
2495 sim_fpu_is_ieee754_1985 (void)
2496 {
2497   return (sim_fpu_current_mode == sim_fpu_ieee754_1985);
2498 }
2499
2500 INLINE_SIM_FPU (bool)
2501 sim_fpu_is_ieee754_2008 (void)
2502 {
2503   return (sim_fpu_current_mode == sim_fpu_ieee754_2008);
2504 }
2505
2506 INLINE_SIM_FPU (void)
2507 sim_fpu_set_mode (const sim_fpu_mode m)
2508 {
2509   sim_fpu_current_mode = m;
2510 }
2511
2512 /* For debugging */
2513
2514 INLINE_SIM_FPU (void)
2515 sim_fpu_print_fpu (const sim_fpu *f,
2516                    sim_fpu_print_func *print,
2517                    void *arg)
2518 {
2519   sim_fpu_printn_fpu (f, print, -1, arg);
2520 }
2521
2522 INLINE_SIM_FPU (void)
2523 sim_fpu_printn_fpu (const sim_fpu *f,
2524                    sim_fpu_print_func *print,
2525                    int digits,
2526                    void *arg)
2527 {
2528   print (arg, "%s", f->sign ? "-" : "+");
2529   switch (f->class)
2530     {
2531     case sim_fpu_class_qnan:
2532       print (arg, "0.");
2533       print_bits (f->fraction, NR_FRAC_GUARD - 1, digits, print, arg);
2534       print (arg, "*QuietNaN");
2535       break;
2536     case sim_fpu_class_snan:
2537       print (arg, "0.");
2538       print_bits (f->fraction, NR_FRAC_GUARD - 1, digits, print, arg);
2539       print (arg, "*SignalNaN");
2540       break;
2541     case sim_fpu_class_zero:
2542       print (arg, "0.0");
2543       break;
2544     case sim_fpu_class_infinity:
2545       print (arg, "INF");
2546       break;
2547     case sim_fpu_class_number:
2548     case sim_fpu_class_denorm:
2549       print (arg, "1.");
2550       print_bits (f->fraction, NR_FRAC_GUARD - 1, digits, print, arg);
2551       print (arg, "*2^%+d", f->normal_exp);
2552       ASSERT (f->fraction >= IMPLICIT_1);
2553       ASSERT (f->fraction < IMPLICIT_2);
2554     }
2555 }
2556
2557
2558 INLINE_SIM_FPU (void)
2559 sim_fpu_print_status (int status,
2560                       sim_fpu_print_func *print,
2561                       void *arg)
2562 {
2563   int i = 1;
2564   const char *prefix = "";
2565   while (status >= i)
2566     {
2567       switch ((sim_fpu_status) (status & i))
2568         {
2569         case sim_fpu_status_denorm:
2570           print (arg, "%sD", prefix);
2571           break;
2572         case sim_fpu_status_invalid_snan:
2573           print (arg, "%sSNaN", prefix);
2574           break;
2575         case sim_fpu_status_invalid_qnan:
2576           print (arg, "%sQNaN", prefix);
2577           break;
2578         case sim_fpu_status_invalid_isi:
2579           print (arg, "%sISI", prefix);
2580           break;
2581         case sim_fpu_status_invalid_idi:
2582           print (arg, "%sIDI", prefix);
2583           break;
2584         case sim_fpu_status_invalid_zdz:
2585           print (arg, "%sZDZ", prefix);
2586           break;
2587         case sim_fpu_status_invalid_imz:
2588           print (arg, "%sIMZ", prefix);
2589           break;
2590         case sim_fpu_status_invalid_cvi:
2591           print (arg, "%sCVI", prefix);
2592           break;
2593         case sim_fpu_status_invalid_cmp:
2594           print (arg, "%sCMP", prefix);
2595           break;
2596         case sim_fpu_status_invalid_sqrt:
2597           print (arg, "%sSQRT", prefix);
2598           break;
2599         case sim_fpu_status_invalid_irx:
2600           print (arg, "%sIRX", prefix);
2601           break;
2602         case sim_fpu_status_inexact:
2603           print (arg, "%sX", prefix);
2604           break;
2605         case sim_fpu_status_overflow:
2606           print (arg, "%sO", prefix);
2607           break;
2608         case sim_fpu_status_underflow:
2609           print (arg, "%sU", prefix);
2610           break;
2611         case sim_fpu_status_invalid_div0:
2612           print (arg, "%s/", prefix);
2613           break;
2614         case sim_fpu_status_rounded:
2615           print (arg, "%sR", prefix);
2616           break;
2617         }
2618       i <<= 1;
2619       prefix = ",";
2620     }
2621 }
2622
2623 #endif