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