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