PowerPC64 .branch_lt address
[binutils-gdb.git] / sim / common / sim-fpu.c
index 99e922c00b65c879dfcd0eb961cb0a18a2009d2a..ccaff9c766187aa4fd1b4c3176e4bd13d26e3b77 100644 (file)
@@ -2,29 +2,20 @@
    of the floating point routines in libgcc1.c for targets without
    hardware floating point.  */
 
-/* Copyright (C) 1994,1997 Free Software Foundation, Inc.
-
-This file is free software; you can redistribute it and/or modify it
-under the terms of the GNU General Public License as published by the
-Free Software Foundation; either version 2, or (at your option) any
-later version.
-
-In addition to the permissions in the GNU General Public License, the
-Free Software Foundation gives you unlimited permission to link the
-compiled version of this file with other programs, and to distribute
-those programs without any restriction coming from the use of this
-file.  (The General Public License restrictions do apply in other
-respects; for example, they cover modification of the file, and
-distribution when not linked into another program.)
-
-This file is distributed in the hope that it will be useful, but
-WITHOUT ANY WARRANTY; without even the implied warranty of
-MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
-General Public License for more details.
+/* Copyright 1994-2022 Free Software Foundation, Inc.
+
+This program is free software; you can redistribute it and/or modify
+it under the terms of the GNU General Public License as published by
+the Free Software Foundation; either version 3 of the License, or
+(at your option) any later version.
+
+This program is distributed in the hope that it will be useful,
+but WITHOUT ANY WARRANTY; without even the implied warranty of
+MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+GNU General Public License for more details.
 
 You should have received a copy of the GNU General Public License
-along with this program; see the file COPYING.  If not, write to
-the Free Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA.  */
+along with this program.  If not, see <http://www.gnu.org/licenses/>.  */
 
 /* As a special exception, if you link this library with other files,
    some of which are compiled with GCC, to produce an executable,
@@ -44,44 +35,55 @@ the Free Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA.  */
 #ifndef SIM_FPU_C
 #define SIM_FPU_C
 
+/* This must come before any other includes.  */
+#include "defs.h"
+
+#include <stdlib.h>
+
 #include "sim-basics.h"
 #include "sim-fpu.h"
 
 #include "sim-io.h"
 #include "sim-assert.h"
 
-
-/* Debugging support. */
+/* Debugging support.
+   If digits is -1, then print all digits.  */
 
 static void
-print_bits (unsigned64 x,
+print_bits (uint64_t x,
            int msbit,
+           int digits,
            sim_fpu_print_func print,
            void *arg)
 {
-  unsigned64 bit = LSBIT64 (msbit);
+  uint64_t bit = LSBIT64 (msbit);
   int i = 4;
-  while (bit)
+  while (bit && digits)
     {
       if (i == 0)
        print (arg, ",");
+
       if ((x & bit))
        print (arg, "1");
       else
        print (arg, "0");
       bit >>= 1;
+
+      if (digits > 0)
+       digits--;
       i = (i + 1) % 4;
     }
 }
 
 
 
-/* Quick and dirty conversion between a host double and host 64bit int */
+/* Quick and dirty conversion between a host double and host 64bit int */
 
-typedef union {
+typedef union
+{
   double d;
-  unsigned64 i;
-} sim_fpu_map;  
+  uint64_t i;
+} sim_fpu_map;
 
 
 /* A packed IEEE floating point number.
@@ -97,7 +99,7 @@ typedef union {
 
    Zero (0 == BIASEDEXP && FRAC == 0):
    (sign ? "-" : "+") 0.0
-   
+
    Infinity (BIASEDEXP == EXPMAX && FRAC == 0):
    (sign ? "-" : "+") "infinity"
 
@@ -164,47 +166,53 @@ typedef union {
 
 /* Integer constants */
 
-#define MAX_INT32  ((signed64) LSMASK64 (30, 0))
+#define MAX_INT32  ((int64_t) LSMASK64 (30, 0))
 #define MAX_UINT32 LSMASK64 (31, 0)
-#define MIN_INT32  ((signed64) LSMASK64 (63, 31))
+#define MIN_INT32  ((int64_t) LSMASK64 (63, 31))
 
-#define MAX_INT64  ((signed64) LSMASK64 (62, 0))
+#define MAX_INT64  ((int64_t) LSMASK64 (62, 0))
 #define MAX_UINT64 LSMASK64 (63, 0)
-#define MIN_INT64  ((signed64) LSMASK64 (63, 63))
+#define MIN_INT64  ((int64_t) LSMASK64 (63, 63))
 
 #define MAX_INT   (is_64bit ? MAX_INT64  : MAX_INT32)
 #define MIN_INT   (is_64bit ? MIN_INT64  : MIN_INT32)
 #define MAX_UINT  (is_64bit ? MAX_UINT64 : MAX_UINT32)
 #define NR_INTBITS (is_64bit ? 64 : 32)
 
-/* Squeese an unpacked sim_fpu struct into a 32/64 bit integer */
-STATIC_INLINE_SIM_FPU (unsigned64)
+/* Squeeze an unpacked sim_fpu struct into a 32/64 bit integer.  */
+STATIC_INLINE_SIM_FPU (uint64_t)
 pack_fpu (const sim_fpu *src,
          int is_double)
 {
   int sign;
-  unsigned64 exp;
-  unsigned64 fraction;
-  unsigned64 packed;
+  uint64_t exp;
+  uint64_t fraction;
+  uint64_t packed;
 
   switch (src->class)
     {
-      /* create a NaN */
+      /* Create a NaN.  */
     case sim_fpu_class_qnan:
       sign = src->sign;
       exp = EXPMAX;
-      /* force fraction to correct class */
+      /* Force fraction to correct class.  */
       fraction = src->fraction;
       fraction >>= NR_GUARDS;
-      fraction |= QUIET_NAN;
+      if (sim_fpu_quiet_nan_inverted)
+       fraction |= QUIET_NAN - 1;
+      else
+       fraction |= QUIET_NAN;
       break;
     case sim_fpu_class_snan:
       sign = src->sign;
       exp = EXPMAX;
-      /* force fraction to correct class */
+      /* Force fraction to correct class.  */
       fraction = src->fraction;
       fraction >>= NR_GUARDS;
-      fraction &= ~QUIET_NAN;
+      if (sim_fpu_quiet_nan_inverted)
+        fraction |= QUIET_NAN;
+      else
+       fraction &= ~QUIET_NAN;
       break;
     case sim_fpu_class_infinity:
       sign = src->sign;
@@ -229,7 +237,7 @@ pack_fpu (const sim_fpu *src,
          int nr_shift = NORMAL_EXPMIN - src->normal_exp;
          if (nr_shift > NR_FRACBITS)
            {
-             /* underflow, just make the number zero */
+             /* Underflow, just make the number zero.  */
              sign = src->sign;
              exp = 0;
              fraction = 0;
@@ -238,7 +246,7 @@ pack_fpu (const sim_fpu *src,
            {
              sign = src->sign;
              exp = 0;
-             /* Shift by the value */
+             /* Shift by the value */
              fraction = src->fraction;
              fraction >>= NR_GUARDS;
              fraction >>= nr_shift;
@@ -249,7 +257,7 @@ pack_fpu (const sim_fpu *src,
          /* Infinity */
          sign = src->sign;
          exp = EXPMAX;
-         fraction = 0; 
+         fraction = 0;
        }
       else
        {
@@ -257,7 +265,7 @@ pack_fpu (const sim_fpu *src,
          sign = src->sign;
          fraction = src->fraction;
          /* FIXME: Need to round according to WITH_SIM_FPU_ROUNDING
-             or some such */
+             or some such */
          /* Round to nearest: If the guard bits are the all zero, but
             the first, then we're half way between two numbers,
             choose the one which makes the lsb of the answer 0.  */
@@ -268,17 +276,17 @@ pack_fpu (const sim_fpu *src,
            }
          else
            {
-             /* Add a one to the guards to force round to nearest */
+             /* Add a one to the guards to force round to nearest */
              fraction += GUARDROUND;
            }
-         if ((fraction & IMPLICIT_2)) /* rounding resulted in carry */
+         if ((fraction & IMPLICIT_2)) /* Rounding resulted in carry.  */
            {
              exp += 1;
              fraction >>= 1;
            }
          fraction >>= NR_GUARDS;
          /* When exp == EXPMAX (overflow from carry) fraction must
-            have been made zero */
+            have been made zero */
          ASSERT ((exp == EXPMAX) <= ((fraction & ~IMPLICIT_1) == 0));
        }
       break;
@@ -290,7 +298,7 @@ pack_fpu (const sim_fpu *src,
             | (exp << NR_FRACBITS)
             | LSMASKED64 (fraction, NR_FRACBITS - 1, 0));
 
-  /* trace operation */
+  /* Trace operation.  */
 #if 0
   if (is_double)
     {
@@ -304,16 +312,16 @@ pack_fpu (const sim_fpu *src,
              (long) LSEXTRACTED32 (packed, 23 - 1, 0));
     }
 #endif
-  
+
   return packed;
 }
 
 
-/* Unpack a 32/64 bit integer into a sim_fpu structure */
+/* Unpack a 32/64 bit integer into a sim_fpu structure */
 STATIC_INLINE_SIM_FPU (void)
-unpack_fpu (sim_fpu *dst, unsigned64 packed, int is_double)
+unpack_fpu (sim_fpu *dst, uint64_t packed, int is_double)
 {
-  unsigned64 fraction = LSMASKED64 (packed, NR_FRACBITS - 1, 0);
+  uint64_t fraction = LSMASKED64 (packed, NR_FRACBITS - 1, 0);
   unsigned exp = LSEXTRACTED64 (packed, NR_EXPBITS + NR_FRACBITS - 1, NR_FRACBITS);
   int sign = (packed & SIGNBIT) != 0;
 
@@ -322,9 +330,10 @@ unpack_fpu (sim_fpu *dst, unsigned64 packed, int is_double)
       /* Hmm.  Looks like 0 */
       if (fraction == 0)
        {
-         /* tastes like zero */
+         /* Tastes like zero.  */
          dst->class = sim_fpu_class_zero;
          dst->sign = sign;
+         dst->normal_exp = 0;
        }
       else
        {
@@ -348,7 +357,7 @@ unpack_fpu (sim_fpu *dst, unsigned64 packed, int is_double)
       /* Huge exponent*/
       if (fraction == 0)
        {
-         /* Attached to a zero fraction - means infinity */
+         /* Attached to a zero fraction - means infinity */
          dst->class = sim_fpu_class_infinity;
          dst->sign = sign;
          /* dst->normal_exp = EXPBIAS; */
@@ -356,10 +365,16 @@ unpack_fpu (sim_fpu *dst, unsigned64 packed, int is_double)
        }
       else
        {
-         /* Non zero fraction, means NaN */
+         int qnan;
+
+         /* Non zero fraction, means NaN.  */
          dst->sign = sign;
          dst->fraction = (fraction << NR_GUARDS);
-         if (fraction >= QUIET_NAN)
+         if (sim_fpu_quiet_nan_inverted)
+           qnan = (fraction & QUIET_NAN) == 0;
+         else
+           qnan = fraction >= QUIET_NAN;
+         if (qnan)
            dst->class = sim_fpu_class_qnan;
          else
            dst->class = sim_fpu_class_snan;
@@ -367,14 +382,14 @@ unpack_fpu (sim_fpu *dst, unsigned64 packed, int is_double)
     }
   else
     {
-      /* Nothing strange about this number */
+      /* Nothing strange about this number */
       dst->class = sim_fpu_class_number;
       dst->sign = sign;
       dst->fraction = ((fraction << NR_GUARDS) | IMPLICIT_1);
       dst->normal_exp = exp - EXPBIAS;
     }
 
-  /* trace operation */
+  /* Trace operation.  */
 #if 0
   if (is_double)
     {
@@ -398,22 +413,22 @@ unpack_fpu (sim_fpu *dst, unsigned64 packed, int is_double)
       }
     else
       {
-       unsigned32 val = pack_fpu (dst, 0);
-       unsigned32 org = packed;
+       uint32_t val = pack_fpu (dst, 0);
+       uint32_t org = packed;
        ASSERT (val == org);
       }
   }
 }
 
 
-/* Convert a floating point into an integer */
+/* Convert a floating point into an integer */
 STATIC_INLINE_SIM_FPU (int)
-fpu2i (signed64 *i,
+fpu2i (int64_t *i,
        const sim_fpu *s,
        int is_64bit,
        sim_fpu_round round)
 {
-  unsigned64 tmp;
+  uint64_t tmp;
   int shift;
   int status = 0;
   if (sim_fpu_is_zero (s))
@@ -431,13 +446,13 @@ fpu2i (signed64 *i,
       *i = MIN_INT; /* FIXME */
       return sim_fpu_status_invalid_cvi;
     }
-  /* map infinity onto MAX_INT... */
+  /* Map infinity onto MAX_INT...  */
   if (sim_fpu_is_infinity (s))
     {
       *i = s->sign ? MIN_INT : MAX_INT;
       return sim_fpu_status_invalid_cvi;
     }
-  /* it is a number, but a small one */
+  /* It is a number, but a small one.  */
   if (s->normal_exp < 0)
     {
       *i = 0;
@@ -452,7 +467,7 @@ fpu2i (signed64 *i,
        return 0; /* exact */
       if (is_64bit) /* can't round */
        return sim_fpu_status_invalid_cvi; /* must be overflow */
-      /* For a 32bit with MAX_INT, rounding is possible */
+      /* For a 32bit with MAX_INT, rounding is possible */
       switch (round)
        {
        case sim_fpu_round_default:
@@ -488,7 +503,7 @@ fpu2i (signed64 *i,
       *i = s->sign ? MIN_INT : MAX_INT;
       return sim_fpu_status_invalid_cvi;
     }
-  /* normal number shift it into place */
+  /* Normal number, shift it into place.  */
   tmp = s->fraction;
   shift = (s->normal_exp - (NR_FRAC_GUARD));
   if (shift > 0)
@@ -506,15 +521,16 @@ fpu2i (signed64 *i,
   return status;
 }
 
-/* convert an integer into a floating point */
+/* Convert an integer into a floating point.  */
 STATIC_INLINE_SIM_FPU (int)
-i2fpu (sim_fpu *f, signed64 i, int is_64bit)
+i2fpu (sim_fpu *f, int64_t i, int is_64bit)
 {
   int status = 0;
   if (i == 0)
     {
       f->class = sim_fpu_class_zero;
       f->sign = 0;
+      f->normal_exp = 0;
     }
   else
     {
@@ -522,10 +538,10 @@ i2fpu (sim_fpu *f, signed64 i, int is_64bit)
       f->sign = (i < 0);
       f->normal_exp = NR_FRAC_GUARD;
 
-      if (f->sign) 
+      if (f->sign)
        {
          /* Special case for minint, since there is no corresponding
-            +ve integer representation for it */
+            +ve integer representation for it */
          if (i == MIN_INT)
            {
              f->fraction = IMPLICIT_1;
@@ -539,9 +555,9 @@ i2fpu (sim_fpu *f, signed64 i, int is_64bit)
 
       if (f->fraction >= IMPLICIT_2)
        {
-         do 
+         do
            {
-             f->fraction >>= 1;
+             f->fraction = (f->fraction >> 1) | (f->fraction & 1);
              f->normal_exp += 1;
            }
          while (f->fraction >= IMPLICIT_2);
@@ -566,7 +582,7 @@ i2fpu (sim_fpu *f, signed64 i, int is_64bit)
 
   /* sanity check */
   {
-    signed64 val;
+    int64_t val;
     fpu2i (&val, f, is_64bit, sim_fpu_round_zero);
     if (i >= MIN_INT32 && i <= MAX_INT32)
       {
@@ -578,12 +594,12 @@ i2fpu (sim_fpu *f, signed64 i, int is_64bit)
 }
 
 
-/* Convert a floating point into an integer */
+/* Convert a floating point into an integer */
 STATIC_INLINE_SIM_FPU (int)
-fpu2u (unsigned64 *u, const sim_fpu *s, int is_64bit)
+fpu2u (uint64_t *u, const sim_fpu *s, int is_64bit)
 {
   const int is_double = 1;
-  unsigned64 tmp;
+  uint64_t tmp;
   int shift;
   if (sim_fpu_is_zero (s))
     {
@@ -595,19 +611,19 @@ fpu2u (unsigned64 *u, const sim_fpu *s, int is_64bit)
       *u = 0;
       return 0;
     }
-  /* it is a negative number */
+  /* It is a negative number.  */
   if (s->sign)
     {
       *u = 0;
       return 0;
     }
-  /* get reasonable MAX_USI_INT... */
+  /* Get reasonable MAX_USI_INT...  */
   if (sim_fpu_is_infinity (s))
     {
       *u = MAX_UINT;
       return 0;
     }
-  /* it is a number, but a small one */
+  /* It is a number, but a small one.  */
   if (s->normal_exp < 0)
     {
       *u = 0;
@@ -635,14 +651,15 @@ fpu2u (unsigned64 *u, const sim_fpu *s, int is_64bit)
   return 0;
 }
 
-/* Convert an unsigned integer into a floating point */
+/* Convert an unsigned integer into a floating point */
 STATIC_INLINE_SIM_FPU (int)
-u2fpu (sim_fpu *f, unsigned64 u, int is_64bit)
+u2fpu (sim_fpu *f, uint64_t u, int is_64bit)
 {
   if (u == 0)
     {
       f->class = sim_fpu_class_zero;
       f->sign = 0;
+      f->normal_exp = 0;
     }
   else
     {
@@ -664,30 +681,30 @@ u2fpu (sim_fpu *f, unsigned64 u, int is_64bit)
 /* register <-> sim_fpu */
 
 INLINE_SIM_FPU (void)
-sim_fpu_32to (sim_fpu *f, unsigned32 s)
+sim_fpu_32to (sim_fpu *f, uint32_t s)
 {
   unpack_fpu (f, s, 0);
 }
 
 
 INLINE_SIM_FPU (void)
-sim_fpu_232to (sim_fpu *f, unsigned32 h, unsigned32 l)
+sim_fpu_232to (sim_fpu *f, uint32_t h, uint32_t l)
 {
-  unsigned64 s = h;
+  uint64_t s = h;
   s = (s << 32) | l;
   unpack_fpu (f, s, 1);
 }
 
 
 INLINE_SIM_FPU (void)
-sim_fpu_64to (sim_fpu *f, unsigned64 s)
+sim_fpu_64to (sim_fpu *f, uint64_t s)
 {
   unpack_fpu (f, s, 1);
 }
 
 
 INLINE_SIM_FPU (void)
-sim_fpu_to32 (unsigned32 *s,
+sim_fpu_to32 (uint32_t *s,
              const sim_fpu *f)
 {
   *s = pack_fpu (f, 0);
@@ -695,23 +712,57 @@ sim_fpu_to32 (unsigned32 *s,
 
 
 INLINE_SIM_FPU (void)
-sim_fpu_to232 (unsigned32 *h, unsigned32 *l,
+sim_fpu_to232 (uint32_t *h, uint32_t *l,
               const sim_fpu *f)
 {
-  unsigned64 s = pack_fpu (f, 1);
+  uint64_t s = pack_fpu (f, 1);
   *l = s;
   *h = (s >> 32);
 }
 
 
 INLINE_SIM_FPU (void)
-sim_fpu_to64 (unsigned64 *u,
+sim_fpu_to64 (uint64_t *u,
              const sim_fpu *f)
 {
   *u = pack_fpu (f, 1);
 }
 
 
+INLINE_SIM_FPU (void)
+sim_fpu_fractionto (sim_fpu *f,
+                   int sign,
+                   int normal_exp,
+                   uint64_t fraction,
+                   int precision)
+{
+  int shift = (NR_FRAC_GUARD - precision);
+  f->class = sim_fpu_class_number;
+  f->sign = sign;
+  f->normal_exp = normal_exp;
+  /* Shift the fraction to where sim-fpu expects it.  */
+  if (shift >= 0)
+    f->fraction = (fraction << shift);
+  else
+    f->fraction = (fraction >> -shift);
+  f->fraction |= IMPLICIT_1;
+}
+
+
+INLINE_SIM_FPU (uint64_t)
+sim_fpu_tofraction (const sim_fpu *d,
+                   int precision)
+{
+  /* We have NR_FRAC_GUARD bits, we want only PRECISION bits.  */
+  int shift = (NR_FRAC_GUARD - precision);
+  uint64_t fraction = (d->fraction & ~IMPLICIT_1);
+  if (shift >= 0)
+    return fraction >> shift;
+  else
+    return fraction << -shift;
+}
+
+
 /* Rounding */
 
 STATIC_INLINE_SIM_FPU (int)
@@ -774,16 +825,16 @@ do_normal_underflow (sim_fpu *f,
 
 
 /* Round a number using NR_GUARDS.
-   Will return the rounded number or F->FRACTION == 0 when underflow */
+   Will return the rounded number or F->FRACTION == 0 when underflow */
 
 STATIC_INLINE_SIM_FPU (int)
 do_normal_round (sim_fpu *f,
                 int nr_guards,
                 sim_fpu_round round)
 {
-  unsigned64 guardmask = LSMASK64 (nr_guards - 1, 0);
-  unsigned64 guardmsb = LSBIT64 (nr_guards - 1);
-  unsigned64 fraclsb = guardmsb << 1;
+  uint64_t guardmask = LSMASK64 (nr_guards - 1, 0);
+  uint64_t guardmsb = LSBIT64 (nr_guards - 1);
+  uint64_t fraclsb = guardmsb << 1;
   if ((f->fraction & guardmask))
     {
       int status = sim_fpu_status_inexact;
@@ -816,7 +867,7 @@ do_normal_round (sim_fpu *f,
          break;
        }
       f->fraction &= ~guardmask;
-      /* round if needed, handle resulting overflow */
+      /* Round if needed, handle resulting overflow.  */
       if ((status & sim_fpu_status_rounded))
        {
          f->fraction += fraclsb;
@@ -847,7 +898,7 @@ do_round (sim_fpu *f,
       return 0;
       break;
     case sim_fpu_class_snan:
-      /* Quieten a SignalingNaN */ 
+      /* Quieten a SignalingNaN.  */
       f->class = sim_fpu_class_qnan;
       return sim_fpu_status_invalid_snan;
       break;
@@ -869,7 +920,7 @@ do_round (sim_fpu *f,
                && !(denorm & sim_fpu_denorm_zero))
              {
                status = do_normal_round (f, shift + NR_GUARDS, round);
-               if (f->fraction == 0) /* rounding underflowed */
+               if (f->fraction == 0) /* Rounding underflowed.  */
                  {
                    status |= do_normal_underflow (f, is_double, round);
                  }
@@ -881,7 +932,7 @@ do_round (sim_fpu *f,
                       before rounding, some after! */
                    if (status & sim_fpu_status_inexact)
                      status |= sim_fpu_status_underflow;
-                   /* Flag that resultant value has been denormalized */
+                   /* Flag that resultant value has been denormalized */
                    f->class = sim_fpu_class_denorm;
                  }
                else if ((denorm & sim_fpu_denorm_underflow_inexact))
@@ -907,7 +958,7 @@ do_round (sim_fpu *f,
              /* f->class = sim_fpu_class_zero; */
              status |= do_normal_underflow (f, is_double, round);
            else if (f->normal_exp > NORMAL_EXPMAX)
-             /* oops! rounding caused overflow */
+             /* Oops! rounding caused overflow.  */
              status |= do_normal_overflow (f, is_double, round);
          }
        ASSERT ((f->class == sim_fpu_class_number
@@ -935,37 +986,58 @@ sim_fpu_round_64 (sim_fpu *f,
   return do_round (f, 1, round, denorm);
 }
 
-
-
-/* Arithmetic ops */
+/* NaN handling for binary operations.  */
 
 INLINE_SIM_FPU (int)
-sim_fpu_add (sim_fpu *f,
-            const sim_fpu *l,
-            const sim_fpu *r)
+sim_fpu_op_nan (sim_fpu *f, const sim_fpu *l, const sim_fpu *r)
 {
-  if (sim_fpu_is_snan (l))
-    {
-      *f = *l;
-      f->class = sim_fpu_class_qnan;
-      return sim_fpu_status_invalid_snan;
-    }
-  if (sim_fpu_is_snan (r))
+  if (sim_fpu_is_snan (l) || sim_fpu_is_snan (r))
     {
-      *f = *r;
+      *f = sim_fpu_is_snan (l) ? *l : *r;
       f->class = sim_fpu_class_qnan;
       return sim_fpu_status_invalid_snan;
     }
-  if (sim_fpu_is_qnan (l))
-    {
-      *f = *l;
-      return 0;
-    }
-  if (sim_fpu_is_qnan (r))
+   ASSERT (sim_fpu_is_nan (l) || sim_fpu_is_nan (r));
+   if (sim_fpu_is_qnan (l))
+     *f = *l;
+   else /* if (sim_fpu_is_qnan (r)) */
+     *f = *r;
+  return 0;
+}
+
+/* NaN handling specific to min/max operations.  */
+
+INLINE_SIM_FPU (int)
+sim_fpu_minmax_nan (sim_fpu *f, const sim_fpu *l, const sim_fpu *r)
+{
+  if (sim_fpu_is_snan (l)
+      || sim_fpu_is_snan (r)
+      || sim_fpu_is_ieee754_1985 ())
+    return sim_fpu_op_nan (f, l, r);
+  else
+    /* if sim_fpu_is_ieee754_2008()
+       && ((sim_fpu_is_qnan (l) || sim_fpu_is_qnan (r)))  */
     {
-      *f = *r;
+      /* In IEEE754-2008:
+        "minNum/maxNum is ... the canonicalized number if one
+        operand is a number and the other a quiet NaN."  */
+      if (sim_fpu_is_qnan (l))
+       *f = *r;
+      else /* if (sim_fpu_is_qnan (r))  */
+       *f = *l;
       return 0;
     }
+}
+
+/* Arithmetic ops */
+
+INLINE_SIM_FPU (int)
+sim_fpu_add (sim_fpu *f,
+            const sim_fpu *l,
+            const sim_fpu *r)
+{
+  if (sim_fpu_is_nan (l) || sim_fpu_is_nan (r))
+    return sim_fpu_op_nan (f, l, r);
   if (sim_fpu_is_infinity (l))
     {
       if (sim_fpu_is_infinity (r)
@@ -1001,18 +1073,18 @@ sim_fpu_add (sim_fpu *f,
   {
     int status = 0;
     int shift = l->normal_exp - r->normal_exp;
-    unsigned64 lfraction;
-    unsigned64 rfraction;
+    uint64_t lfraction;
+    uint64_t rfraction;
     /* use exp of larger */
     if (shift >= NR_FRAC_GUARD)
       {
-       /* left has much bigger magnitute */
+       /* left has much bigger magnitude */
        *f = *l;
        return sim_fpu_status_inexact;
       }
     if (shift <= - NR_FRAC_GUARD)
       {
-       /* right has much bigger magnitute */
+       /* right has much bigger magnitude */
        *f = *r;
        return sim_fpu_status_inexact;
       }
@@ -1024,7 +1096,7 @@ sim_fpu_add (sim_fpu *f,
        if (rfraction & LSMASK64 (shift - 1, 0))
          {
            status |= sim_fpu_status_inexact;
-           rfraction |= LSBIT64 (shift); /* stick LSBit */
+           rfraction |= LSBIT64 (shift); /* Stick LSBit.  */
          }
        rfraction >>= shift;
       }
@@ -1034,7 +1106,7 @@ sim_fpu_add (sim_fpu *f,
        if (lfraction & LSMASK64 (- shift - 1, 0))
          {
            status |= sim_fpu_status_inexact;
-           lfraction |= LSBIT64 (- shift); /* stick LSBit */
+           lfraction |= LSBIT64 (- shift); /* Stick LSBit.  */
          }
        lfraction >>= -shift;
       }
@@ -1043,7 +1115,7 @@ sim_fpu_add (sim_fpu *f,
        f->normal_exp = r->normal_exp;
       }
 
-    /* perform the addition */
+    /* Perform the addition.  */
     if (l->sign)
       lfraction = - lfraction;
     if (r->sign)
@@ -1059,7 +1131,7 @@ sim_fpu_add (sim_fpu *f,
 
     /* sign? */
     f->class = sim_fpu_class_number;
-    if ((signed64) f->fraction >= 0)
+    if (((int64_t) f->fraction) >= 0)
       f->sign = 0;
     else
       {
@@ -1067,7 +1139,7 @@ sim_fpu_add (sim_fpu *f,
        f->fraction = - f->fraction;
       }
 
-    /* normalize it */
+    /* Normalize it.  */
     if ((f->fraction & IMPLICIT_2))
       {
        f->fraction = (f->fraction >> 1) | (f->fraction & 1);
@@ -1093,28 +1165,8 @@ sim_fpu_sub (sim_fpu *f,
             const sim_fpu *l,
             const sim_fpu *r)
 {
-  if (sim_fpu_is_snan (l))
-    {
-      *f = *l;
-      f->class = sim_fpu_class_qnan;
-      return sim_fpu_status_invalid_snan;
-    }
-  if (sim_fpu_is_snan (r))
-    {
-      *f = *r;
-      f->class = sim_fpu_class_qnan;
-      return sim_fpu_status_invalid_snan;
-    }
-  if (sim_fpu_is_qnan (l))
-    {
-      *f = *l;
-      return 0;
-    }
-  if (sim_fpu_is_qnan (r))
-    {
-      *f = *r;
-      return 0;
-    }
+  if (sim_fpu_is_nan (l) || sim_fpu_is_nan (r))
+    return sim_fpu_op_nan (f, l, r);
   if (sim_fpu_is_infinity (l))
     {
       if (sim_fpu_is_infinity (r)
@@ -1154,18 +1206,18 @@ sim_fpu_sub (sim_fpu *f,
   {
     int status = 0;
     int shift = l->normal_exp - r->normal_exp;
-    unsigned64 lfraction;
-    unsigned64 rfraction;
+    uint64_t lfraction;
+    uint64_t rfraction;
     /* use exp of larger */
     if (shift >= NR_FRAC_GUARD)
       {
-       /* left has much bigger magnitute */
+       /* left has much bigger magnitude */
        *f = *l;
        return sim_fpu_status_inexact;
       }
     if (shift <= - NR_FRAC_GUARD)
       {
-       /* right has much bigger magnitute */
+       /* right has much bigger magnitude */
        *f = *r;
        f->sign = !r->sign;
        return sim_fpu_status_inexact;
@@ -1178,7 +1230,7 @@ sim_fpu_sub (sim_fpu *f,
        if (rfraction & LSMASK64 (shift - 1, 0))
          {
            status |= sim_fpu_status_inexact;
-           rfraction |= LSBIT64 (shift); /* stick LSBit */
+           rfraction |= LSBIT64 (shift); /* Stick LSBit.  */
          }
        rfraction >>= shift;
       }
@@ -1188,7 +1240,7 @@ sim_fpu_sub (sim_fpu *f,
        if (lfraction & LSMASK64 (- shift - 1, 0))
          {
            status |= sim_fpu_status_inexact;
-           lfraction |= LSBIT64 (- shift); /* stick LSBit */
+           lfraction |= LSBIT64 (- shift); /* Stick LSBit.  */
          }
        lfraction >>= -shift;
       }
@@ -1197,7 +1249,7 @@ sim_fpu_sub (sim_fpu *f,
        f->normal_exp = r->normal_exp;
       }
 
-    /* perform the subtraction */
+    /* Perform the subtraction.  */
     if (l->sign)
       lfraction = - lfraction;
     if (!r->sign)
@@ -1213,7 +1265,7 @@ sim_fpu_sub (sim_fpu *f,
 
     /* sign? */
     f->class = sim_fpu_class_number;
-    if ((signed64) f->fraction >= 0)
+    if (((int64_t) f->fraction) >= 0)
       f->sign = 0;
     else
       {
@@ -1221,7 +1273,7 @@ sim_fpu_sub (sim_fpu *f,
        f->fraction = - f->fraction;
       }
 
-    /* normalize it */
+    /* Normalize it.  */
     if ((f->fraction & IMPLICIT_2))
       {
        f->fraction = (f->fraction >> 1) | (f->fraction & 1);
@@ -1247,28 +1299,8 @@ sim_fpu_mul (sim_fpu *f,
             const sim_fpu *l,
             const sim_fpu *r)
 {
-  if (sim_fpu_is_snan (l))
-    {
-      *f = *l;
-      f->class = sim_fpu_class_qnan;
-      return sim_fpu_status_invalid_snan;
-    }
-  if (sim_fpu_is_snan (r))
-    {
-      *f = *r;
-      f->class = sim_fpu_class_qnan;
-      return sim_fpu_status_invalid_snan;
-    }
-  if (sim_fpu_is_qnan (l))
-    {
-      *f = *l;
-      return 0;
-    }
-  if (sim_fpu_is_qnan (r))
-    {
-      *f = *r;
-      return 0;
-    }
+  if (sim_fpu_is_nan (l) || sim_fpu_is_nan (r))
+    return sim_fpu_op_nan (f, l, r);
   if (sim_fpu_is_infinity (l))
     {
       if (sim_fpu_is_zero (r))
@@ -1298,21 +1330,21 @@ sim_fpu_mul (sim_fpu *f,
       return 0;
     }
   /* Calculate the mantissa by multiplying both 64bit numbers to get a
-     128 bit number */
+     128 bit number */
   {
-    unsigned64 low;
-    unsigned64 high;
-    unsigned64 nl = l->fraction & 0xffffffff;
-    unsigned64 nh = l->fraction >> 32;
-    unsigned64 ml = r->fraction & 0xffffffff;
-    unsigned64 mh = r->fraction >>32;
-    unsigned64 pp_ll = ml * nl;
-    unsigned64 pp_hl = mh * nl;
-    unsigned64 pp_lh = ml * nh;
-    unsigned64 pp_hh = mh * nh;
-    unsigned64 res2 = 0;
-    unsigned64 res0 = 0;
-    unsigned64 ps_hh__ = pp_hl + pp_lh;
+    uint64_t low;
+    uint64_t high;
+    uint64_t nl = l->fraction & 0xffffffff;
+    uint64_t nh = l->fraction >> 32;
+    uint64_t ml = r->fraction & 0xffffffff;
+    uint64_t mh = r->fraction >>32;
+    uint64_t pp_ll = ml * nl;
+    uint64_t pp_hl = mh * nl;
+    uint64_t pp_lh = ml * nh;
+    uint64_t pp_hh = mh * nh;
+    uint64_t res2 = 0;
+    uint64_t res0 = 0;
+    uint64_t ps_hh__ = pp_hl + pp_lh;
     if (ps_hh__ < pp_hl)
       res2 += UNSIGNED64 (0x100000000);
     pp_hl = (ps_hh__ << 32) & UNSIGNED64 (0xffffffff00000000);
@@ -1322,14 +1354,14 @@ sim_fpu_mul (sim_fpu *f,
     res2 += ((ps_hh__ >> 32) & 0xffffffff) + pp_hh;
     high = res2;
     low = res0;
-    
+
     f->normal_exp = l->normal_exp + r->normal_exp;
     f->sign = l->sign ^ r->sign;
     f->class = sim_fpu_class_number;
 
     /* Input is bounded by [1,2)   ;   [2^60,2^61)
        Output is bounded by [1,4)  ;   [2^120,2^122) */
+
     /* Adjust the exponent according to where the decimal point ended
        up in the high 64 bit word.  In the source the decimal point
        was at NR_FRAC_GUARD. */
@@ -1341,15 +1373,7 @@ sim_fpu_mul (sim_fpu *f,
     ASSERT (high >= LSBIT64 ((NR_FRAC_GUARD * 2) - 64));
     ASSERT (LSBIT64 (((NR_FRAC_GUARD + 1) * 2) - 64) < IMPLICIT_1);
 
-#if 0
-    printf ("\n");
-    print_bits (high, 63, (sim_fpu_print_func*)fprintf, stdout);
-    printf (";");
-    print_bits (low, 63, (sim_fpu_print_func*)fprintf, stdout);
-    printf ("\n");
-#endif
-
-    /* normalize */
+    /* Normalize.  */
     do
       {
        f->normal_exp--;
@@ -1360,13 +1384,6 @@ sim_fpu_mul (sim_fpu *f,
       }
     while (high < IMPLICIT_1);
 
-#if 0
-    print_bits (high, 63, (sim_fpu_print_func*)fprintf, stdout);
-    printf (";");
-    print_bits (low, 63, (sim_fpu_print_func*)fprintf, stdout);
-    printf ("\n");
-#endif
-
     ASSERT (high >= IMPLICIT_1 && high < IMPLICIT_2);
     if (low != 0)
       {
@@ -1387,30 +1404,8 @@ sim_fpu_div (sim_fpu *f,
             const sim_fpu *l,
             const sim_fpu *r)
 {
-  if (sim_fpu_is_snan (l))
-    {
-      *f = *l;
-      f->class = sim_fpu_class_qnan;
-      return sim_fpu_status_invalid_snan;
-    }
-  if (sim_fpu_is_snan (r))
-    {
-      *f = *r;
-      f->class = sim_fpu_class_qnan;
-      return sim_fpu_status_invalid_snan;
-    }
-  if (sim_fpu_is_qnan (l))
-    {
-      *f = *l;
-      f->class = sim_fpu_class_qnan;
-      return 0;
-    }
-  if (sim_fpu_is_qnan (r))
-    {
-      *f = *r;
-      f->class = sim_fpu_class_qnan;
-      return 0;
-    }
+  if (sim_fpu_is_nan (l) || sim_fpu_is_nan (r))
+    return sim_fpu_op_nan (f, l, r);
   if (sim_fpu_is_infinity (l))
     {
       if (sim_fpu_is_infinity (r))
@@ -1453,15 +1448,15 @@ sim_fpu_div (sim_fpu *f,
     }
 
   /* Calculate the mantissa by multiplying both 64bit numbers to get a
-     128 bit number */
+     128 bit number */
   {
     /* quotient =  ( ( numerator / denominator)
                       x 2^(numerator exponent -  denominator exponent)
      */
-    unsigned64 numerator;
-    unsigned64 denominator;
-    unsigned64 quotient;
-    unsigned64 bit;
+    uint64_t numerator;
+    uint64_t denominator;
+    uint64_t quotient;
+    uint64_t bit;
 
     f->class = sim_fpu_class_number;
     f->sign = l->sign ^ r->sign;
@@ -1477,8 +1472,8 @@ sim_fpu_div (sim_fpu *f,
        f->normal_exp--;
       }
     ASSERT (numerator >= denominator);
-    
-    /* Gain extra precision, already used one spare bit */
+
+    /* Gain extra precision, already used one spare bit */
     numerator <<=    NR_SPARE;
     denominator <<=  NR_SPARE;
 
@@ -1496,17 +1491,7 @@ sim_fpu_div (sim_fpu *f,
        numerator <<= 1;
       }
 
-#if 0
-    printf ("\n");
-    print_bits (quotient, 63, (sim_fpu_print_func*)fprintf, stdout);
-    printf ("\n");
-    print_bits (numerator, 63, (sim_fpu_print_func*)fprintf, stdout);
-    printf ("\n");
-    print_bits (denominator, 63, (sim_fpu_print_func*)fprintf, stdout);
-    printf ("\n");
-#endif
-
-    /* discard (but save) the extra bits */
+    /* Discard (but save) the extra bits.  */
     if ((quotient & LSMASK64 (NR_SPARE -1, 0)))
       quotient = (quotient >> NR_SPARE) | 1;
     else
@@ -1516,7 +1501,7 @@ sim_fpu_div (sim_fpu *f,
     ASSERT (f->fraction >= IMPLICIT_1 && f->fraction < IMPLICIT_2);
     if (numerator != 0)
       {
-       f->fraction |= 1; /* stick remaining bits */
+       f->fraction |= 1; /* Stick remaining bits.  */
        return sim_fpu_status_inexact;
       }
     else
@@ -1526,32 +1511,73 @@ sim_fpu_div (sim_fpu *f,
 
 
 INLINE_SIM_FPU (int)
-sim_fpu_max (sim_fpu *f,
+sim_fpu_rem (sim_fpu *f,
             const sim_fpu *l,
             const sim_fpu *r)
 {
-  if (sim_fpu_is_snan (l))
+  if (sim_fpu_is_nan (l) || sim_fpu_is_nan (r))
+    return sim_fpu_op_nan (f, l, r);
+  if (sim_fpu_is_infinity (l))
     {
-      *f = *l;
-      f->class = sim_fpu_class_qnan;
-      return sim_fpu_status_invalid_snan;
+      *f = sim_fpu_qnan;
+      return sim_fpu_status_invalid_irx;
     }
-  if (sim_fpu_is_snan (r))
+  if (sim_fpu_is_zero (r))
     {
-      *f = *r;
-      f->class = sim_fpu_class_qnan;
-      return sim_fpu_status_invalid_snan;
+      *f = sim_fpu_qnan;
+      return sim_fpu_status_invalid_div0;
     }
-  if (sim_fpu_is_qnan (l))
+  if (sim_fpu_is_zero (l))
     {
       *f = *l;
       return 0;
     }
-  if (sim_fpu_is_qnan (r))
+  if (sim_fpu_is_infinity (r))
     {
-      *f = *r;
+      *f = *l;
       return 0;
     }
+  {
+    sim_fpu n, tmp;
+
+    /* Remainder is calculated as l-n*r, where n is l/r rounded to the
+       nearest integer.  The variable n is rounded half even.  */
+
+    sim_fpu_div (&n, l, r);
+    sim_fpu_round_64 (&n, 0, 0);
+
+    if (n.normal_exp < -1) /* If n looks like zero just return l.  */
+      {
+       *f = *l;
+       return 0;
+      }
+    else if (n.class == sim_fpu_class_number
+            && n.normal_exp <= (NR_FRAC_GUARD)) /* If not too large round.  */
+      do_normal_round (&n, (NR_FRAC_GUARD) - n.normal_exp, sim_fpu_round_near);
+
+    /* Mark 0's as zero so multiply can detect zero.  */
+    if (n.fraction == 0)
+      n.class = sim_fpu_class_zero;
+
+    /* Calculate n*r.  */
+    sim_fpu_mul (&tmp, &n, r);
+    sim_fpu_round_64 (&tmp, 0, 0);
+
+    /* Finally calculate l-n*r.  */
+    sim_fpu_sub (f, l, &tmp);
+
+    return 0;
+  }
+}
+
+
+INLINE_SIM_FPU (int)
+sim_fpu_max (sim_fpu *f,
+            const sim_fpu *l,
+            const sim_fpu *r)
+{
+  if (sim_fpu_is_nan (l) || sim_fpu_is_nan (r))
+    return sim_fpu_minmax_nan (f, l, r);
   if (sim_fpu_is_infinity (l))
     {
       if (sim_fpu_is_infinity (r)
@@ -1563,7 +1589,7 @@ sim_fpu_max (sim_fpu *f,
       if (l->sign)
        *f = *r; /* -inf < anything */
       else
-       *f = *l; /* +inf > anthing */
+       *f = *l; /* +inf > anything */
       return 0;
     }
   if (sim_fpu_is_infinity (r))
@@ -1571,7 +1597,7 @@ sim_fpu_max (sim_fpu *f,
       if (r->sign)
        *f = *l; /* anything > -inf */
       else
-       *f = *r; /* anthing < +inf */
+       *f = *r; /* anything < +inf */
       return 0;
     }
   if (l->sign > r->sign)
@@ -1586,8 +1612,8 @@ sim_fpu_max (sim_fpu *f,
     }
   ASSERT (l->sign == r->sign);
   if (l->normal_exp > r->normal_exp
-      || (l->normal_exp == r->normal_exp && 
-         l->fraction > r->fraction))
+      || (l->normal_exp == r->normal_exp
+         && l->fraction > r->fraction))
     {
       /* |l| > |r| */
       if (l->sign)
@@ -1613,28 +1639,8 @@ sim_fpu_min (sim_fpu *f,
             const sim_fpu *l,
             const sim_fpu *r)
 {
-  if (sim_fpu_is_snan (l))
-    {
-      *f = *l;
-      f->class = sim_fpu_class_qnan;
-      return sim_fpu_status_invalid_snan;
-    }
-  if (sim_fpu_is_snan (r))
-    {
-      *f = *r;
-      f->class = sim_fpu_class_qnan;
-      return sim_fpu_status_invalid_snan;
-    }
-  if (sim_fpu_is_qnan (l))
-    {
-      *f = *l;
-      return 0;
-    }
-  if (sim_fpu_is_qnan (r))
-    {
-      *f = *r;
-      return 0;
-    }
+  if (sim_fpu_is_nan (l) || sim_fpu_is_nan (r))
+    return sim_fpu_minmax_nan (f, l, r);
   if (sim_fpu_is_infinity (l))
     {
       if (sim_fpu_is_infinity (r)
@@ -1669,8 +1675,8 @@ sim_fpu_min (sim_fpu *f,
     }
   ASSERT (l->sign == r->sign);
   if (l->normal_exp > r->normal_exp
-      || (l->normal_exp == r->normal_exp && 
-         l->fraction > r->fraction))
+      || (l->normal_exp == r->normal_exp
+         && l->fraction > r->fraction))
     {
       /* |l| > |r| */
       if (l->sign)
@@ -1695,7 +1701,7 @@ INLINE_SIM_FPU (int)
 sim_fpu_neg (sim_fpu *f,
             const sim_fpu *r)
 {
-  if (sim_fpu_is_snan (r))
+  if (sim_fpu_is_ieee754_1985 () && sim_fpu_is_snan (r))
     {
       *f = *r;
       f->class = sim_fpu_class_qnan;
@@ -1716,19 +1722,13 @@ INLINE_SIM_FPU (int)
 sim_fpu_abs (sim_fpu *f,
             const sim_fpu *r)
 {
-  if (sim_fpu_is_snan (r))
+  *f = *r;
+  f->sign = 0;
+  if (sim_fpu_is_ieee754_1985 () && sim_fpu_is_snan (r))
     {
-      *f = *r;
       f->class = sim_fpu_class_qnan;
       return sim_fpu_status_invalid_snan;
     }
-  if (sim_fpu_is_qnan (r))
-    {
-      *f = *r;
-      return 0;
-    }
-  *f = *r;
-  f->sign = 0;
   return 0;
 }
 
@@ -1737,33 +1737,7 @@ INLINE_SIM_FPU (int)
 sim_fpu_inv (sim_fpu *f,
             const sim_fpu *r)
 {
-  if (sim_fpu_is_snan (r))
-    {
-      *f = *r;
-      f->class = sim_fpu_class_qnan;
-      return sim_fpu_status_invalid_snan;
-    }
-  if (sim_fpu_is_qnan (r))
-    {
-      *f = *r;
-      f->class = sim_fpu_class_qnan;
-      return 0;
-    }
-  if (sim_fpu_is_infinity (r))
-    {
-      *f = sim_fpu_zero;
-      f->sign = r->sign;
-      return 0;
-    }
-  if (sim_fpu_is_zero (r))
-    {
-      f->class = sim_fpu_class_infinity;
-      f->sign = r->sign;
-      return sim_fpu_status_invalid_div0;
-    }
-  *f = *r;
-  f->normal_exp = - r->normal_exp;
-  return 0;
+  return sim_fpu_div (f, &sim_fpu_one, r);
 }
 
 
@@ -1785,6 +1759,7 @@ sim_fpu_sqrt (sim_fpu *f,
     {
       f->class = sim_fpu_class_zero;
       f->sign = r->sign;
+      f->normal_exp = 0;
       return 0;
     }
   if (sim_fpu_is_infinity (r))
@@ -1815,20 +1790,20 @@ sim_fpu_sqrt (sim_fpu *f,
    *
    * Developed at SunPro, a Sun Microsystems, Inc. business.
    * Permission to use, copy, modify, and distribute this
-   * software is freely granted, provided that this notice 
+   * software is freely granted, provided that this notice
    * is preserved.
    * ====================================================
    */
-  
+
   /* __ieee754_sqrt(x)
    * Return correctly rounded sqrt.
    *           ------------------------------------------
    *           |  Use the hardware sqrt if you have one |
    *           ------------------------------------------
-   * Method: 
-   *   Bit by bit method using integer arithmetic. (Slow, but portable) 
+   * Method:
+   *   Bit by bit method using integer arithmetic. (Slow, but portable)
    *   1. Normalization
-   *   Scale x to y in [1,4) with even powers of 2: 
+   *   Scale x to y in [1,4) with even powers of 2:
    *   find an integer k such that  1 <= (y=x*2^(2k)) < 4, then
    *           sqrt(x) = 2^k * sqrt(y)
    -
@@ -1848,9 +1823,9 @@ sim_fpu_sqrt (sim_fpu *f,
    *                                     i+1         2
    *       s  = 2*q , and      y  =  2   * ( y - q  ).         (1)
    *        i      i            i                 i
-   *                                                        
-   *   To compute q    from q , one checks whether 
-   *               i+1       i                       
+   *
+   *   To compute q    from q , one checks whether
+   *               i+1       i
    *
    *                         -(i+1) 2
    *                   (q + 2      ) <= y.                     (2)
@@ -1859,13 +1834,13 @@ sim_fpu_sqrt (sim_fpu *f,
    *   If (2) is false, then q   = q ; otherwise q   = q  + 2      .
    *                          i+1   i             i+1   i
    *
-   *   With some algebric manipulation, it is not difficult to see
-   *   that (2) is equivalent to 
+   *   With some algebraic manipulation, it is not difficult to see
+   *   that (2) is equivalent to
    *                             -(i+1)
    *                   s  +  2       <= y                      (3)
    *                    i                i
    *
-   *   The advantage of (3) is that s  and y  can be computed by 
+   *   The advantage of (3) is that s  and y  can be computed by
    *                                 i      i
    *   the following recurrence formula:
    *       if (3) is false
@@ -1881,15 +1856,15 @@ sim_fpu_sqrt (sim_fpu *f,
    *                       -i                      -(i+1)
    *       s     =  s  + 2  ,  y    = y  -  s  - 2             (5)
    *         i+1      i          i+1    i     i
-   *                           
+   *
    -
    -                                                   -(i+1)
-   -                      NOTE: y    = 2 (y  -  s  -  2      )                 
+   -                      NOTE: y    = 2 (y  -  s  -  2      )
    -                             i+1       i     i
    -
-   *   One may easily use induction to prove (4) and (5). 
+   *   One may easily use induction to prove (4) and (5).
    *   Note. Since the left hand side of (3) contain only i+2 bits,
-   *         it does not necessary to do a full (53-bit) comparison 
+   *         it does not necessary to do a full (53-bit) comparison
    *         in (3).
    *   3. Final rounding
    *   After generating the 53 bits result, we compute one more bit.
@@ -1899,30 +1874,30 @@ sim_fpu_sqrt (sim_fpu *f,
    *   The rounding mode can be detected by checking whether
    *   huge + tiny is equal to huge, and whether huge - tiny is
    *   equal to huge for some floating point number "huge" and "tiny".
-   *           
+   *
    * Special cases:
    *   sqrt(+-0) = +-0         ... exact
    *   sqrt(inf) = inf
    *   sqrt(-ve) = NaN         ... with invalid signal
-   *   sqrt(NaN) = NaN         ... with invalid signal for signaling NaN
+   *   sqrt(NaN) = NaN         ... with invalid signal for signalling NaN
    *
    * Other methods : see the appended file at the end of the program below.
    *---------------
    */
 
   {
-    /* generate sqrt(x) bit by bit */
-    unsigned64 y;
-    unsigned64 q;
-    unsigned64 s;
-    unsigned64 b;
+    /* Generate sqrt(x) bit by bit.  */
+    uint64_t y;
+    uint64_t q;
+    uint64_t s;
+    uint64_t b;
 
     f->class = sim_fpu_class_number;
     f->sign = 0;
     y = r->fraction;
     f->normal_exp = (r->normal_exp >> 1);      /* exp = [exp/2] */
 
-    /* odd exp, double x to make it even */
+    /* Odd exp, double x to make it even.  */
     ASSERT (y >= IMPLICIT_1 && y < IMPLICIT_4);
     if ((r->normal_exp & 1))
       {
@@ -1934,10 +1909,10 @@ sim_fpu_sqrt (sim_fpu *f,
     b = IMPLICIT_1;
     q = 0;
     s = 0;
-    
+
     while (b)
       {
-       unsigned64 t = s + b;
+       uint64_t t = s + b;
        if (t <= y)
          {
            s |= (b << 1);
@@ -1952,7 +1927,7 @@ sim_fpu_sqrt (sim_fpu *f,
     f->fraction = q;
     if (y != 0)
       {
-       f->fraction |= 1; /* stick remaining bits */
+       f->fraction |= 1; /* Stick remaining bits.  */
        return sim_fpu_status_inexact;
       }
     else
@@ -1965,7 +1940,7 @@ sim_fpu_sqrt (sim_fpu *f,
 
 INLINE_SIM_FPU (int)
 sim_fpu_i32to (sim_fpu *f,
-              signed32 i,
+              int32_t i,
               sim_fpu_round round)
 {
   i2fpu (f, i, 0);
@@ -1974,7 +1949,7 @@ sim_fpu_i32to (sim_fpu *f,
 
 INLINE_SIM_FPU (int)
 sim_fpu_u32to (sim_fpu *f,
-              unsigned32 u,
+              uint32_t u,
               sim_fpu_round round)
 {
   u2fpu (f, u, 0);
@@ -1983,7 +1958,7 @@ sim_fpu_u32to (sim_fpu *f,
 
 INLINE_SIM_FPU (int)
 sim_fpu_i64to (sim_fpu *f,
-              signed64 i,
+              int64_t i,
               sim_fpu_round round)
 {
   i2fpu (f, i, 1);
@@ -1992,7 +1967,7 @@ sim_fpu_i64to (sim_fpu *f,
 
 INLINE_SIM_FPU (int)
 sim_fpu_u64to (sim_fpu *f,
-              unsigned64 u,
+              uint64_t u,
               sim_fpu_round round)
 {
   u2fpu (f, u, 1);
@@ -2001,29 +1976,29 @@ sim_fpu_u64to (sim_fpu *f,
 
 
 INLINE_SIM_FPU (int)
-sim_fpu_to32i (signed32 *i,
+sim_fpu_to32i (int32_t *i,
               const sim_fpu *f,
               sim_fpu_round round)
 {
-  signed64 i64;
+  int64_t i64;
   int status = fpu2i (&i64, f, 0, round);
   *i = i64;
   return status;
 }
 
 INLINE_SIM_FPU (int)
-sim_fpu_to32u (unsigned32 *u,
+sim_fpu_to32u (uint32_t *u,
               const sim_fpu *f,
               sim_fpu_round round)
 {
-  unsigned64 u64;
+  uint64_t u64;
   int status = fpu2u (&u64, f, 0);
   *u = u64;
   return status;
 }
 
 INLINE_SIM_FPU (int)
-sim_fpu_to64i (signed64 *i,
+sim_fpu_to64i (int64_t *i,
               const sim_fpu *f,
               sim_fpu_round round)
 {
@@ -2032,7 +2007,7 @@ sim_fpu_to64i (signed64 *i,
 
 
 INLINE_SIM_FPU (int)
-sim_fpu_to64u (unsigned64 *u,
+sim_fpu_to64u (uint64_t *u,
               const sim_fpu *f,
               sim_fpu_round round)
 {
@@ -2056,7 +2031,17 @@ INLINE_SIM_FPU (double)
 sim_fpu_2d (const sim_fpu *s)
 {
   sim_fpu_map val;
-  val.i = pack_fpu (s, 1);
+  if (sim_fpu_is_snan (s))
+    {
+      /* gag SNaN's */
+      sim_fpu n = *s;
+      n.class = sim_fpu_class_qnan;
+      val.i = pack_fpu (&n, 1);
+    }
+  else
+    {
+      val.i = pack_fpu (s, 1);
+    }
   return val.d;
 }
 
@@ -2186,6 +2171,23 @@ sim_fpu_exp (const sim_fpu *d)
 }
 
 
+INLINE_SIM_FPU (uint64_t)
+sim_fpu_fraction (const sim_fpu *d)
+{
+  return d->fraction;
+}
+
+
+INLINE_SIM_FPU (uint64_t)
+sim_fpu_guard (const sim_fpu *d, int is_double)
+{
+  uint64_t rv;
+  uint64_t guardmask = LSMASK64 (NR_GUARDS - 1, 0);
+  rv = (d->fraction & guardmask) >> NR_PAD;
+  return rv;
+}
+
+
 INLINE_SIM_FPU (int)
 sim_fpu_is (const sim_fpu *d)
 {
@@ -2196,8 +2198,10 @@ sim_fpu_is (const sim_fpu *d)
     case sim_fpu_class_snan:
       return SIM_FPU_IS_SNAN;
     case sim_fpu_class_infinity:
-      return SIM_FPU_IS_NINF;
-      return SIM_FPU_IS_PINF;
+      if (d->sign)
+       return SIM_FPU_IS_NINF;
+      else
+       return SIM_FPU_IS_PINF;
     case sim_fpu_class_number:
       if (d->sign)
        return SIM_FPU_IS_NNUMBER;
@@ -2275,6 +2279,21 @@ sim_fpu_is_gt (const sim_fpu *l, const sim_fpu *r)
   return is;
 }
 
+INLINE_SIM_FPU (int)
+sim_fpu_is_un (const sim_fpu *l, const sim_fpu *r)
+{
+  int is;
+  sim_fpu_un (&is, l, r);
+  return is;
+}
+
+INLINE_SIM_FPU (int)
+sim_fpu_is_or (const sim_fpu *l, const sim_fpu *r)
+{
+  int is;
+  sim_fpu_or (&is, l, r);
+  return is;
+}
 
 /* Compare operators */
 
@@ -2398,28 +2417,97 @@ sim_fpu_gt (int *is,
   return sim_fpu_lt (is, r, l);
 }
 
+INLINE_SIM_FPU (int)
+sim_fpu_un (int *is, const sim_fpu *l, const sim_fpu *r)
+{
+  if (sim_fpu_is_nan (l) || sim_fpu_is_nan (r))
+   {
+    *is = 1;
+    return 0;
+   }
+
+  *is = 0;
+  return 0;
+}
+
+INLINE_SIM_FPU (int)
+sim_fpu_or (int *is, const sim_fpu *l, const sim_fpu *r)
+{
+  sim_fpu_un (is, l, r);
+
+  /* Invert result.  */
+  *is = !*is;
+  return 0;
+}
+
+INLINE_SIM_FPU(int)
+sim_fpu_classify (const sim_fpu *f)
+{
+  switch (f->class)
+    {
+    case sim_fpu_class_snan: return SIM_FPU_IS_SNAN;
+    case sim_fpu_class_qnan: return SIM_FPU_IS_QNAN;
+    case sim_fpu_class_infinity:
+      return f->sign ? SIM_FPU_IS_NINF : SIM_FPU_IS_PINF;
+    case sim_fpu_class_zero:
+      return f->sign ? SIM_FPU_IS_NZERO : SIM_FPU_IS_PZERO;
+    case sim_fpu_class_number:
+      return f->sign ? SIM_FPU_IS_NNUMBER : SIM_FPU_IS_PNUMBER;
+    case sim_fpu_class_denorm:
+      return f->sign ? SIM_FPU_IS_NDENORM : SIM_FPU_IS_PDENORM;
+    default:
+      fprintf (stderr, "Bad switch\n");
+      abort ();
+    }
+  return 0;
+}
 
 /* A number of useful constants */
 
-EXTERN_SIM_FPU (const sim_fpu) sim_fpu_zero = {
-  sim_fpu_class_zero,
+#if EXTERN_SIM_FPU_P
+sim_fpu_state _sim_fpu = {
+  .quiet_nan_inverted = false,
+  .current_mode = sim_fpu_ieee754_1985,
 };
-EXTERN_SIM_FPU (const sim_fpu) sim_fpu_qnan = {
-  sim_fpu_class_qnan,
+
+const sim_fpu sim_fpu_zero = {
+  sim_fpu_class_zero, 0, 0, 0
 };
-EXTERN_SIM_FPU (const sim_fpu) sim_fpu_one = {
-  sim_fpu_class_number, 0, IMPLICIT_1, 1
+const sim_fpu sim_fpu_qnan = {
+  sim_fpu_class_qnan, 0, 0, 0
 };
-EXTERN_SIM_FPU (const sim_fpu) sim_fpu_two = {
-  sim_fpu_class_number, 0, IMPLICIT_1, 2
+const sim_fpu sim_fpu_one = {
+  sim_fpu_class_number, 0, IMPLICIT_1, 0
 };
-EXTERN_SIM_FPU (const sim_fpu) sim_fpu_max32 = {
+const sim_fpu sim_fpu_two = {
+  sim_fpu_class_number, 0, IMPLICIT_1, 1
+};
+const sim_fpu sim_fpu_max32 = {
   sim_fpu_class_number, 0, LSMASK64 (NR_FRAC_GUARD, NR_GUARDS32), NORMAL_EXPMAX32
 };
-EXTERN_SIM_FPU (const sim_fpu) sim_fpu_max64 = {
+const sim_fpu sim_fpu_max64 = {
   sim_fpu_class_number, 0, LSMASK64 (NR_FRAC_GUARD, NR_GUARDS64), NORMAL_EXPMAX64
 };
+#endif
 
+/* Specification swapping behaviour */
+INLINE_SIM_FPU (bool)
+sim_fpu_is_ieee754_1985 (void)
+{
+  return (sim_fpu_current_mode == sim_fpu_ieee754_1985);
+}
+
+INLINE_SIM_FPU (bool)
+sim_fpu_is_ieee754_2008 (void)
+{
+  return (sim_fpu_current_mode == sim_fpu_ieee754_2008);
+}
+
+INLINE_SIM_FPU (void)
+sim_fpu_set_mode (const sim_fpu_mode m)
+{
+  sim_fpu_current_mode = m;
+}
 
 /* For debugging */
 
@@ -2427,18 +2515,27 @@ INLINE_SIM_FPU (void)
 sim_fpu_print_fpu (const sim_fpu *f,
                   sim_fpu_print_func *print,
                   void *arg)
+{
+  sim_fpu_printn_fpu (f, print, -1, arg);
+}
+
+INLINE_SIM_FPU (void)
+sim_fpu_printn_fpu (const sim_fpu *f,
+                  sim_fpu_print_func *print,
+                  int digits,
+                  void *arg)
 {
   print (arg, "%s", f->sign ? "-" : "+");
   switch (f->class)
     {
     case sim_fpu_class_qnan:
       print (arg, "0.");
-      print_bits (f->fraction, NR_FRAC_GUARD - 1, print, arg);
+      print_bits (f->fraction, NR_FRAC_GUARD - 1, digits, print, arg);
       print (arg, "*QuietNaN");
       break;
     case sim_fpu_class_snan:
       print (arg, "0.");
-      print_bits (f->fraction, NR_FRAC_GUARD - 1, print, arg);
+      print_bits (f->fraction, NR_FRAC_GUARD - 1, digits, print, arg);
       print (arg, "*SignalNaN");
       break;
     case sim_fpu_class_zero:
@@ -2450,8 +2547,8 @@ sim_fpu_print_fpu (const sim_fpu *f,
     case sim_fpu_class_number:
     case sim_fpu_class_denorm:
       print (arg, "1.");
-      print_bits (f->fraction, NR_FRAC_GUARD - 1, print, arg);
-      print (arg, "*2^%+-5d", f->normal_exp);
+      print_bits (f->fraction, NR_FRAC_GUARD - 1, digits, print, arg);
+      print (arg, "*2^%+d", f->normal_exp);
       ASSERT (f->fraction >= IMPLICIT_1);
       ASSERT (f->fraction < IMPLICIT_2);
     }
@@ -2464,7 +2561,7 @@ sim_fpu_print_status (int status,
                      void *arg)
 {
   int i = 1;
-  char *prefix = "";
+  const char *prefix = "";
   while (status >= i)
     {
       switch ((sim_fpu_status) (status & i))
@@ -2499,27 +2596,24 @@ sim_fpu_print_status (int status,
        case sim_fpu_status_invalid_sqrt:
          print (arg, "%sSQRT", prefix);
          break;
+       case sim_fpu_status_invalid_irx:
+         print (arg, "%sIRX", prefix);
          break;
        case sim_fpu_status_inexact:
          print (arg, "%sX", prefix);
          break;
-         break;
        case sim_fpu_status_overflow:
          print (arg, "%sO", prefix);
          break;
-         break;
        case sim_fpu_status_underflow:
          print (arg, "%sU", prefix);
          break;
-         break;
        case sim_fpu_status_invalid_div0:
          print (arg, "%s/", prefix);
          break;
-         break;
        case sim_fpu_status_rounded:
          print (arg, "%sR", prefix);
          break;
-         break;
        }
       i <<= 1;
       prefix = ",";