libgcc/config/m68k/fpgnulib.c

   1 /* This is a stripped down version of floatlib.c.  It supplies only those
   2    functions which exist in libgcc, but for which there is not assembly
   3    language versions in m68k/lb1sf68.S.
   4
   5    It also includes simplistic support for extended floats (by working in
   6    double precision).  You must compile this file again with -DEXTFLOAT
   7    to get this support.  */
   8
   9 /*
  10 ** gnulib support for software floating point.
  11 ** Copyright (C) 1991 by Pipeline Associates, Inc.  All rights reserved.
  12 ** Permission is granted to do *anything* you want with this file,
  13 ** commercial or otherwise, provided this message remains intact.  So there!
  14 ** I would appreciate receiving any updates/patches/changes that anyone
  15 ** makes, and am willing to be the repository for said changes (am I
  16 ** making a big mistake?).
  17 **
  18 ** Pat Wood
  19 ** Pipeline Associates, Inc.
  20 ** pipeline!phw@motown.com or
  21 ** sun!pipeline!phw or
  22 ** uunet!motown!pipeline!phw
  23 **
  24 ** 05/01/91 -- V1.0 -- first release to gcc mailing lists
  25 ** 05/04/91 -- V1.1 -- added float and double prototypes and return values
  26 **                  -- fixed problems with adding and subtracting zero
  27 **                  -- fixed rounding in truncdfsf2
  28 **                  -- fixed SWAP define and tested on 386
  29 */
  30
  31 /*
  32 ** The following are routines that replace the gnulib soft floating point
  33 ** routines that are called automatically when -msoft-float is selected.
  34 ** The support single and double precision IEEE format, with provisions
  35 ** for byte-swapped machines (tested on 386).  Some of the double-precision
  36 ** routines work at full precision, but most of the hard ones simply punt
  37 ** and call the single precision routines, producing a loss of accuracy.
  38 ** long long support is not assumed or included.
  39 ** Overall accuracy is close to IEEE (actually 68882) for single-precision
  40 ** arithmetic.  I think there may still be a 1 in 1000 chance of a bit
  41 ** being rounded the wrong way during a multiply.  I'm not fussy enough to
  42 ** bother with it, but if anyone is, knock yourself out.
  43 **
  44 ** Efficiency has only been addressed where it was obvious that something
  45 ** would make a big difference.  Anyone who wants to do this right for
  46 ** best speed should go in and rewrite in assembler.
  47 **
  48 ** I have tested this only on a 68030 workstation and 386/ix integrated
  49 ** in with -msoft-float.
  50 */
  51
  52 /* the following deal with IEEE single-precision numbers */
  53 #define EXCESS          126L
  54 #define SIGNBIT         0x80000000L
  55 #define HIDDEN          (1L << 23L)
  56 #define SIGN(fp)        ((fp) & SIGNBIT)
  57 #define EXPMASK         0xFFL
  58 #define EXP(fp)         (((fp) >> 23L) & 0xFF)
  59 #define MANT(fp)        (((fp) & 0x7FFFFFL) | HIDDEN)
  60 #define PACK(s,e,m)     ((s) | ((e) << 23L) | (m))
  61
  62 /* the following deal with IEEE double-precision numbers */
  63 #define EXCESSD         1022L
  64 #define HIDDEND         (1L << 20L)
  65 #define EXPDBITS        11
  66 #define EXPDMASK        0x7FFL
  67 #define EXPD(fp)        (((fp.l.upper) >> 20L) & 0x7FFL)
  68 #define SIGND(fp)       ((fp.l.upper) & SIGNBIT)
  69 #define MANTD(fp)       (((((fp.l.upper) & 0xFFFFF) | HIDDEND) << 10) | \
  70                                 (fp.l.lower >> 22))
  71 #define MANTDMASK       0xFFFFFL /* mask of upper part */
  72
  73 /* the following deal with IEEE extended-precision numbers */
  74 #define EXCESSX         16382L
  75 #define HIDDENX         (1L << 31L)
  76 #define EXPXBITS        15
  77 #define EXPXMASK        0x7FFF
  78 #define EXPX(fp)        (((fp.l.upper) >> 16) & EXPXMASK)
  79 #define SIGNX(fp)       ((fp.l.upper) & SIGNBIT)
  80 #define MANTXMASK       0x7FFFFFFFL /* mask of upper part */
  81
  82 union double_long
  83 {
  84   double d;
  85   struct {
  86       long upper;
  87       unsigned long lower;
  88     } l;
  89 };
  90
  91 union float_long {
  92   float f;
  93   long l;
  94 };
  95
  96 union long_double_long
  97 {
  98   long double ld;
  99   struct
 100     {
 101       long upper;
 102       unsigned long middle;
 103       unsigned long lower;
 104     } l;
 105 };
 106 \f
 107 #ifndef EXTFLOAT
 108
 109 int
 110 __unordsf2(float a, float b)
 111 {
 112   union float_long fl;
 113
 114   fl.f = a;
 115   if (EXP(fl.l) == EXP(~0u) && (MANT(fl.l) & ~HIDDEN) != 0)
 116     return 1;
 117   fl.f = b;
 118   if (EXP(fl.l) == EXP(~0u) && (MANT(fl.l) & ~HIDDEN) != 0)
 119     return 1;
 120   return 0;
 121 }
 122
 123 int
 124 __unorddf2(double a, double b)
 125 {
 126   union double_long dl;
 127
 128   dl.d = a;
 129   if (EXPD(dl) == EXPDMASK
 130       && ((dl.l.upper & MANTDMASK) != 0 || dl.l.lower != 0))
 131     return 1;
 132   dl.d = b;
 133   if (EXPD(dl) == EXPDMASK
 134       && ((dl.l.upper & MANTDMASK) != 0 || dl.l.lower != 0))
 135     return 1;
 136   return 0;
 137 }
 138
 139 /* convert unsigned int to double */
 140 double
 141 __floatunsidf (unsigned long a1)
 142 {
 143   long exp = 32 + EXCESSD;
 144   union double_long dl;
 145
 146   if (!a1)
 147     {
 148       dl.l.upper = dl.l.lower = 0;
 149       return dl.d;
 150     }
 151
 152   while (a1 < 0x2000000L)
 153     {
 154       a1 <<= 4;
 155       exp -= 4;
 156     }
 157
 158   while (a1 < 0x80000000L)
 159     {
 160       a1 <<= 1;
 161       exp--;
 162     }
 163
 164   /* pack up and go home */
 165   dl.l.upper = exp << 20L;
 166   dl.l.upper |= (a1 >> 11L) & ~HIDDEND;
 167   dl.l.lower = a1 << 21L;
 168
 169   return dl.d;
 170 }
 171
 172 /* convert int to double */
 173 double
 174 __floatsidf (long a1)
 175 {
 176   long sign = 0, exp = 31 + EXCESSD;
 177   union double_long dl;
 178
 179   if (!a1)
 180     {
 181       dl.l.upper = dl.l.lower = 0;
 182       return dl.d;
 183     }
 184
 185   if (a1 < 0)
 186     {
 187       sign = SIGNBIT;
 188       a1 = (long)-(unsigned long)a1;
 189       if (a1 < 0)
 190         {
 191           dl.l.upper = SIGNBIT | ((32 + EXCESSD) << 20L);
 192           dl.l.lower = 0;
 193           return dl.d;
 194         }
 195     }
 196
 197   while (a1 < 0x1000000L)
 198     {
 199       a1 <<= 4;
 200       exp -= 4;
 201     }
 202
 203   while (a1 < 0x40000000L)
 204     {
 205       a1 <<= 1;
 206       exp--;
 207     }
 208
 209   /* pack up and go home */
 210   dl.l.upper = sign;
 211   dl.l.upper |= exp << 20L;
 212   dl.l.upper |= (a1 >> 10L) & ~HIDDEND;
 213   dl.l.lower = a1 << 22L;
 214
 215   return dl.d;
 216 }
 217
 218 /* convert unsigned int to float */
 219 float
 220 __floatunsisf (unsigned long l)
 221 {
 222   double foo = __floatunsidf (l);
 223   return foo;
 224 }
 225
 226 /* convert int to float */
 227 float
 228 __floatsisf (long l)
 229 {
 230   double foo = __floatsidf (l);
 231   return foo;
 232 }
 233
 234 /* convert float to double */
 235 double
 236 __extendsfdf2 (float a1)
 237 {
 238   register union float_long fl1;
 239   register union double_long dl;
 240   register long exp;
 241   register long mant;
 242
 243   fl1.f = a1;
 244
 245   dl.l.upper = SIGN (fl1.l);
 246   if ((fl1.l & ~SIGNBIT) == 0)
 247     {
 248       dl.l.lower = 0;
 249       return dl.d;
 250     }
 251
 252   exp = EXP(fl1.l);
 253   mant = MANT (fl1.l) & ~HIDDEN;
 254   if (exp == 0)
 255     {
 256       /* Denormal.  */
 257       exp = 1;
 258       while (!(mant & HIDDEN))
 259         {
 260           mant <<= 1;
 261           exp--;
 262         }
 263       mant &= ~HIDDEN;
 264     }
 265   exp = exp - EXCESS + EXCESSD;
 266   /* Handle inf and NaN */
 267   if (exp == EXPMASK - EXCESS + EXCESSD)
 268     exp = EXPDMASK;
 269   dl.l.upper |= exp << 20;
 270   dl.l.upper |= mant >> 3;
 271   dl.l.lower = mant << 29;
 272
 273   return dl.d;
 274 }
 275
 276 /* convert double to float */
 277 float
 278 __truncdfsf2 (double a1)
 279 {
 280   register long exp;
 281   register long mant;
 282   register union float_long fl;
 283   register union double_long dl1;
 284   int sticky;
 285   int shift;
 286
 287   dl1.d = a1;
 288
 289   if ((dl1.l.upper & ~SIGNBIT) == 0 && !dl1.l.lower)
 290     {
 291       fl.l = SIGND(dl1);
 292       return fl.f;
 293     }
 294
 295   exp = EXPD (dl1) - EXCESSD + EXCESS;
 296
 297   sticky = dl1.l.lower & ((1 << 22) - 1);
 298   mant = MANTD (dl1);
 299   /* shift double mantissa 6 bits so we can round */
 300   sticky |= mant & ((1 << 6) - 1);
 301   mant >>= 6;
 302   if (exp == EXPDMASK - EXCESSD + EXCESS)
 303     {
 304       exp = EXPMASK;
 305       mant = (mant >> 1) | (mant & 1) | (!!sticky);
 306     }
 307   else
 308     {
 309       /* Check for underflow and denormals.  */
 310       if (exp <= 0)
 311         {
 312           if (exp < -24)
 313             {
 314               sticky |= mant;
 315               mant = 0;
 316             }
 317           else
 318             {
 319               sticky |= mant & ((1 << (1 - exp)) - 1);
 320               mant >>= 1 - exp;
 321             }
 322           exp = 0;
 323         }
 324
 325       /* now round */
 326       shift = 1;
 327       if ((mant & 1) && (sticky || (mant & 2)))
 328         {
 329           int rounding = exp ? 2 : 1;
 330
 331           mant += 1;
 332
 333           /* did the round overflow? */
 334           if (mant >= (HIDDEN << rounding))
 335             {
 336               exp++;
 337               shift = rounding;
 338             }
 339         }
 340       /* shift down */
 341       mant >>= shift;
 342       if (exp >= EXPMASK)
 343         {
 344           exp = EXPMASK;
 345           mant = 0;
 346         }
 347     }
 348
 349   mant &= ~HIDDEN;
 350
 351   /* pack up and go home */
 352   fl.l = PACK (SIGND (dl1), exp, mant);
 353   return (fl.f);
 354 }
 355
 356 /* convert double to int */
 357 long
 358 __fixdfsi (double a1)
 359 {
 360   register union double_long dl1;
 361   register long exp;
 362   register long l;
 363
 364   dl1.d = a1;
 365
 366   if (!dl1.l.upper && !dl1.l.lower)
 367     return 0;
 368
 369   exp = EXPD (dl1) - EXCESSD - 31;
 370   l = MANTD (dl1);
 371
 372   if (exp > 0)
 373     {
 374       /* Return largest integer.  */
 375       return SIGND (dl1) ? 0x80000000L : 0x7fffffffL;
 376     }
 377
 378   if (exp <= -32)
 379     return 0;
 380
 381   /* shift down until exp = 0 */
 382   if (exp < 0)
 383     l >>= -exp;
 384
 385   return (SIGND (dl1) ? -l : l);
 386 }
 387
 388 /* convert float to int */
 389 long
 390 __fixsfsi (float a1)
 391 {
 392   double foo = a1;
 393   return __fixdfsi (foo);
 394 }
 395 \f
 396 #else /* EXTFLOAT */
 397
 398 /* We do not need these routines for coldfire, as it has no extended
 399    float format. */
 400 #if !defined (__mcoldfire__)
 401
 402 /* Primitive extended precision floating point support.
 403
 404    We assume all numbers are normalized, don't do any rounding, etc.  */
 405
 406 /* Prototypes for the above in case we use them.  */
 407 double __floatunsidf (unsigned long);
 408 double __floatsidf (long);
 409 float __floatsisf (long);
 410 double __extendsfdf2 (float);
 411 float __truncdfsf2 (double);
 412 long __fixdfsi (double);
 413 long __fixsfsi (float);
 414 int __cmpdf2 (double, double);
 415
 416 int
 417 __unordxf2(long double a, long double b)
 418 {
 419   union long_double_long ldl;
 420
 421   ldl.ld = a;
 422   if (EXPX(ldl) == EXPXMASK
 423       && ((ldl.l.middle & MANTXMASK) != 0 || ldl.l.lower != 0))
 424     return 1;
 425   ldl.ld = b;
 426   if (EXPX(ldl) == EXPXMASK
 427       && ((ldl.l.middle & MANTXMASK) != 0 || ldl.l.lower != 0))
 428     return 1;
 429   return 0;
 430 }
 431
 432 /* convert double to long double */
 433 long double
 434 __extenddfxf2 (double d)
 435 {
 436   register union double_long dl;
 437   register union long_double_long ldl;
 438   register long exp;
 439
 440   dl.d = d;
 441   /*printf ("dfxf in: %g\n", d);*/
 442
 443   ldl.l.upper = SIGND (dl);
 444   if ((dl.l.upper & ~SIGNBIT) == 0 && !dl.l.lower)
 445     {
 446       ldl.l.middle = 0;
 447       ldl.l.lower = 0;
 448       return ldl.ld;
 449     }
 450
 451   exp = EXPD (dl) - EXCESSD + EXCESSX;
 452   /* Check for underflow and denormals. */
 453   if (exp < 0)
 454     {
 455       if (exp < -53)
 456         {
 457           ldl.l.middle = 0;
 458           ldl.l.lower = 0;
 459         }
 460       else if (exp < -30)
 461         {
 462           ldl.l.lower = (ldl.l.middle & MANTXMASK) >> ((1 - exp) - 32);
 463           ldl.l.middle &= ~MANTXMASK;
 464         }
 465       else
 466         {
 467           ldl.l.lower >>= 1 - exp;
 468           ldl.l.lower |= (ldl.l.middle & MANTXMASK) << (32 - (1 - exp));
 469           ldl.l.middle = (ldl.l.middle & ~MANTXMASK) | (ldl.l.middle & MANTXMASK >> (1 - exp));
 470         }
 471       exp = 0;
 472     }
 473   /* Handle inf and NaN */
 474   if (exp == EXPDMASK - EXCESSD + EXCESSX)
 475     exp = EXPXMASK;
 476   ldl.l.upper |= exp << 16;
 477   ldl.l.middle = HIDDENX;
 478   /* 31-20: # mantissa bits in ldl.l.middle - # mantissa bits in dl.l.upper */
 479   ldl.l.middle |= (dl.l.upper & MANTDMASK) << (31 - 20);
 480   /* 1+20: explicit-integer-bit + # mantissa bits in dl.l.upper */
 481   ldl.l.middle |= dl.l.lower >> (1 + 20);
 482   /* 32 - 21: # bits of dl.l.lower in ldl.l.middle */
 483   ldl.l.lower = dl.l.lower << (32 - 21);
 484
 485   /*printf ("dfxf out: %s\n", dumpxf (ldl.ld));*/
 486   return ldl.ld;
 487 }
 488
 489 /* convert long double to double */
 490 double
 491 __truncxfdf2 (long double ld)
 492 {
 493   register long exp;
 494   register union double_long dl;
 495   register union long_double_long ldl;
 496
 497   ldl.ld = ld;
 498   /*printf ("xfdf in: %s\n", dumpxf (ld));*/
 499
 500   dl.l.upper = SIGNX (ldl);
 501   if ((ldl.l.upper & ~SIGNBIT) == 0 && !ldl.l.middle && !ldl.l.lower)
 502     {
 503       dl.l.lower = 0;
 504       return dl.d;
 505     }
 506
 507   exp = EXPX (ldl) - EXCESSX + EXCESSD;
 508   /* Check for underflow and denormals. */
 509   if (exp <= 0)
 510     {
 511       if (exp < -53)
 512         {
 513           ldl.l.middle = 0;
 514           ldl.l.lower = 0;
 515         }
 516       else if (exp < -30)
 517         {
 518           ldl.l.lower = (ldl.l.middle & MANTXMASK) >> ((1 - exp) - 32);
 519           ldl.l.middle &= ~MANTXMASK;
 520         }
 521       else
 522         {
 523           ldl.l.lower >>= 1 - exp;
 524           ldl.l.lower |= (ldl.l.middle & MANTXMASK) << (32 - (1 - exp));
 525           ldl.l.middle = (ldl.l.middle & ~MANTXMASK) | (ldl.l.middle & MANTXMASK >> (1 - exp));
 526         }
 527       exp = 0;
 528     }
 529   else if (exp == EXPXMASK - EXCESSX + EXCESSD)
 530     {
 531       exp = EXPDMASK;
 532       ldl.l.middle |= ldl.l.lower;
 533     }
 534   else if (exp >= EXPDMASK)
 535     {
 536       exp = EXPDMASK;
 537       ldl.l.middle = 0;
 538       ldl.l.lower = 0;
 539     }
 540   dl.l.upper |= exp << (32 - (EXPDBITS + 1));
 541   /* +1-1: add one for sign bit, but take one off for explicit-integer-bit */
 542   dl.l.upper |= (ldl.l.middle & MANTXMASK) >> (EXPDBITS + 1 - 1);
 543   dl.l.lower = (ldl.l.middle & MANTXMASK) << (32 - (EXPDBITS + 1 - 1));
 544   dl.l.lower |= ldl.l.lower >> (EXPDBITS + 1 - 1);
 545
 546   /*printf ("xfdf out: %g\n", dl.d);*/
 547   return dl.d;
 548 }
 549
 550 /* convert a float to a long double */
 551 long double
 552 __extendsfxf2 (float f)
 553 {
 554   long double foo = __extenddfxf2 (__extendsfdf2 (f));
 555   return foo;
 556 }
 557
 558 /* convert a long double to a float */
 559 float
 560 __truncxfsf2 (long double ld)
 561 {
 562   float foo = __truncdfsf2 (__truncxfdf2 (ld));
 563   return foo;
 564 }
 565
 566 /* convert an int to a long double */
 567 long double
 568 __floatsixf (long l)
 569 {
 570   double foo = __floatsidf (l);
 571   return foo;
 572 }
 573
 574 /* convert an unsigned int to a long double */
 575 long double
 576 __floatunsixf (unsigned long l)
 577 {
 578   double foo = __floatunsidf (l);
 579   return foo;
 580 }
 581
 582 /* convert a long double to an int */
 583 long
 584 __fixxfsi (long double a)
 585 {
 586   union long_double_long ldl;
 587   long exp;
 588   long l;
 589
 590   ldl.ld = a;
 591
 592   exp = EXPX (ldl);
 593   if (exp == 0 && ldl.l.middle == 0 && ldl.l.lower == 0)
 594     return 0;
 595
 596   exp = exp - EXCESSX - 63;
 597
 598   if (exp > 0)
 599     {
 600       /* Return largest integer.  */
 601       return SIGNX (ldl) ? 0x80000000L : 0x7fffffffL;
 602     }
 603
 604   if (exp <= -64)
 605     return 0;
 606
 607   if (exp <= -32)
 608     {
 609       ldl.l.lower = ldl.l.middle >> (-exp - 32);
 610     }
 611   else if (exp < 0)
 612     {
 613       ldl.l.lower = ldl.l.lower >> -exp;
 614       ldl.l.lower |= ldl.l.middle << (32 + exp);
 615     }
 616
 617   return SIGNX (ldl) ? -ldl.l.lower : ldl.l.lower;
 618 }
 619
 620 /* The remaining provide crude math support by working in double precision.  */
 621
 622 long double
 623 __addxf3 (long double x1, long double x2)
 624 {
 625   return (double) x1 + (double) x2;
 626 }
 627
 628 long double
 629 __subxf3 (long double x1, long double x2)
 630 {
 631   return (double) x1 - (double) x2;
 632 }
 633
 634 long double
 635 __mulxf3 (long double x1, long double x2)
 636 {
 637   return (double) x1 * (double) x2;
 638 }
 639
 640 long double
 641 __divxf3 (long double x1, long double x2)
 642 {
 643   return (double) x1 / (double) x2;
 644 }
 645
 646 long double
 647 __negxf2 (long double x1)
 648 {
 649   return - (double) x1;
 650 }
 651
 652 long
 653 __cmpxf2 (long double x1, long double x2)
 654 {
 655   return __cmpdf2 ((double) x1, (double) x2);
 656 }
 657
 658 long
 659 __eqxf2 (long double x1, long double x2)
 660 {
 661   return __cmpdf2 ((double) x1, (double) x2);
 662 }
 663
 664 long
 665 __nexf2 (long double x1, long double x2)
 666 {
 667   return __cmpdf2 ((double) x1, (double) x2);
 668 }
 669
 670 long
 671 __ltxf2 (long double x1, long double x2)
 672 {
 673   return __cmpdf2 ((double) x1, (double) x2);
 674 }
 675
 676 long
 677 __lexf2 (long double x1, long double x2)
 678 {
 679   return __cmpdf2 ((double) x1, (double) x2);
 680 }
 681
 682 long
 683 __gtxf2 (long double x1, long double x2)
 684 {
 685   return __cmpdf2 ((double) x1, (double) x2);
 686 }
 687
 688 long
 689 __gexf2 (long double x1, long double x2)
 690 {
 691   return __cmpdf2 ((double) x1, (double) x2);
 692 }
 693
 694 #endif /* !__mcoldfire__ */
 695 #endif /* EXTFLOAT */