gcc/config/fp-bit.c

   1 /* This is a software floating point library which can be used instead of
   2    the floating point routines in libgcc1.c for targets without hardware
   3    floating point.
   4    Copyright (C) 1994-1998, 2000 Free Software Foundation, Inc.
   5
   6 This file is free software; you can redistribute it and/or modify it
   7 under the terms of the GNU General Public License as published by the
   8 Free Software Foundation; either version 2, or (at your option) any
   9 later version.
  10
  11 In addition to the permissions in the GNU General Public License, the
  12 Free Software Foundation gives you unlimited permission to link the
  13 compiled version of this file with other programs, and to distribute
  14 those programs without any restriction coming from the use of this
  15 file.  (The General Public License restrictions do apply in other
  16 respects; for example, they cover modification of the file, and
  17 distribution when not linked into another program.)
  18
  19 This file is distributed in the hope that it will be useful, but
  20 WITHOUT ANY WARRANTY; without even the implied warranty of
  21 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  22 General Public License for more details.
  23
  24 You should have received a copy of the GNU General Public License
  25 along with this program; see the file COPYING.  If not, write to
  26 the Free Software Foundation, 59 Temple Place - Suite 330,
  27 Boston, MA 02111-1307, 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 /* The intended way to use this file is to make two copies, add `#define FLOAT'
  44    to one copy, then compile both copies and add them to libgcc.a.  */
  45
  46 /* Defining FINE_GRAINED_LIBRARIES allows one to select which routines
  47    from this file are compiled via additional -D options.
  48
  49    This avoids the need to pull in the entire fp emulation library
  50    when only a small number of functions are needed.
  51
  52    If FINE_GRAINED_LIBRARIES is not defined, then compile every
  53    suitable routine.  */
  54 #ifndef FINE_GRAINED_LIBRARIES
  55 #define L_pack_df
  56 #define L_unpack_df
  57 #define L_pack_sf
  58 #define L_unpack_sf
  59 #define L_addsub_sf
  60 #define L_addsub_df
  61 #define L_mul_sf
  62 #define L_mul_df
  63 #define L_div_sf
  64 #define L_div_df
  65 #define L_fpcmp_parts_sf
  66 #define L_fpcmp_parts_df
  67 #define L_compare_sf
  68 #define L_compare_df
  69 #define L_eq_sf
  70 #define L_eq_df
  71 #define L_ne_sf
  72 #define L_ne_df
  73 #define L_gt_sf
  74 #define L_gt_df
  75 #define L_ge_sf
  76 #define L_ge_df
  77 #define L_lt_sf
  78 #define L_lt_df
  79 #define L_le_sf
  80 #define L_le_df
  81 #define L_unord_sf
  82 #define L_unord_df
  83 #define L_si_to_sf
  84 #define L_si_to_df
  85 #define L_sf_to_si
  86 #define L_df_to_si
  87 #define L_f_to_usi
  88 #define L_df_to_usi
  89 #define L_negate_sf
  90 #define L_negate_df
  91 #define L_make_sf
  92 #define L_make_df
  93 #define L_sf_to_df
  94 #define L_df_to_sf
  95 #ifdef FLOAT
  96 #define L_thenan_sf
  97 #else
  98 #define L_thenan_df
  99 #endif
 100 #endif
 101
 102 /* The following macros can be defined to change the behaviour of this file:
 103    FLOAT: Implement a `float', aka SFmode, fp library.  If this is not
 104      defined, then this file implements a `double', aka DFmode, fp library.
 105    FLOAT_ONLY: Used with FLOAT, to implement a `float' only library, i.e.
 106      don't include float->double conversion which requires the double library.
 107      This is useful only for machines which can't support doubles, e.g. some
 108      8-bit processors.
 109    CMPtype: Specify the type that floating point compares should return.
 110      This defaults to SItype, aka int.
 111    US_SOFTWARE_GOFAST: This makes all entry points use the same names as the
 112      US Software goFast library.  If this is not defined, the entry points use
 113      the same names as libgcc1.c.
 114    _DEBUG_BITFLOAT: This makes debugging the code a little easier, by adding
 115      two integers to the FLO_union_type.
 116    NO_NANS: Disable nan and infinity handling
 117    SMALL_MACHINE: Useful when operations on QIs and HIs are faster
 118      than on an SI */
 119
 120 /* We don't currently support extended floats (long doubles) on machines
 121    without hardware to deal with them.
 122
 123    These stubs are just to keep the linker from complaining about unresolved
 124    references which can be pulled in from libio & libstdc++, even if the
 125    user isn't using long doubles.  However, they may generate an unresolved
 126    external to abort if abort is not used by the function, and the stubs
 127    are referenced from within libc, since libgcc goes before and after the
 128    system library.  */
 129
 130 #ifdef EXTENDED_FLOAT_STUBS
 131 __truncxfsf2 (){ abort(); }
 132 __extendsfxf2 (){ abort(); }
 133 __addxf3 (){ abort(); }
 134 __divxf3 (){ abort(); }
 135 __eqxf2 (){ abort(); }
 136 __extenddfxf2 (){ abort(); }
 137 __gtxf2 (){ abort(); }
 138 __lexf2 (){ abort(); }
 139 __ltxf2 (){ abort(); }
 140 __mulxf3 (){ abort(); }
 141 __negxf2 (){ abort(); }
 142 __nexf2 (){ abort(); }
 143 __subxf3 (){ abort(); }
 144 __truncxfdf2 (){ abort(); }
 145
 146 __trunctfsf2 (){ abort(); }
 147 __extendsftf2 (){ abort(); }
 148 __addtf3 (){ abort(); }
 149 __divtf3 (){ abort(); }
 150 __eqtf2 (){ abort(); }
 151 __extenddftf2 (){ abort(); }
 152 __gttf2 (){ abort(); }
 153 __letf2 (){ abort(); }
 154 __lttf2 (){ abort(); }
 155 __multf3 (){ abort(); }
 156 __negtf2 (){ abort(); }
 157 __netf2 (){ abort(); }
 158 __subtf3 (){ abort(); }
 159 __trunctfdf2 (){ abort(); }
 160 __gexf2 (){ abort(); }
 161 __fixxfsi (){ abort(); }
 162 __floatsixf (){ abort(); }
 163 #else   /* !EXTENDED_FLOAT_STUBS, rest of file */
 164
 165
 166 typedef float SFtype __attribute__ ((mode (SF)));
 167 typedef float DFtype __attribute__ ((mode (DF)));
 168
 169 typedef int HItype __attribute__ ((mode (HI)));
 170 typedef int SItype __attribute__ ((mode (SI)));
 171 typedef int DItype __attribute__ ((mode (DI)));
 172
 173 /* The type of the result of a fp compare */
 174 #ifndef CMPtype
 175 #define CMPtype SItype
 176 #endif
 177
 178 typedef unsigned int UHItype __attribute__ ((mode (HI)));
 179 typedef unsigned int USItype __attribute__ ((mode (SI)));
 180 typedef unsigned int UDItype __attribute__ ((mode (DI)));
 181
 182 #define MAX_SI_INT   ((SItype) ((unsigned) (~0)>>1))
 183 #define MAX_USI_INT  ((USItype) ~0)
 184
 185
 186 #ifdef FLOAT_ONLY
 187 #define NO_DI_MODE
 188 #endif
 189
 190 #ifdef FLOAT
 191 #       define NGARDS    7L
 192 #       define GARDROUND 0x3f
 193 #       define GARDMASK  0x7f
 194 #       define GARDMSB   0x40
 195 #       define EXPBITS 8
 196 #       define EXPBIAS 127
 197 #       define FRACBITS 23
 198 #       define EXPMAX (0xff)
 199 #       define QUIET_NAN 0x100000L
 200 #       define FRAC_NBITS 32
 201 #       define FRACHIGH  0x80000000L
 202 #       define FRACHIGH2 0xc0000000L
 203 #       define pack_d __pack_f
 204 #       define unpack_d __unpack_f
 205 #       define __fpcmp_parts __fpcmp_parts_f
 206         typedef USItype fractype;
 207         typedef UHItype halffractype;
 208         typedef SFtype FLO_type;
 209         typedef SItype intfrac;
 210
 211 #else
 212 #       define PREFIXFPDP dp
 213 #       define PREFIXSFDF df
 214 #       define NGARDS 8L
 215 #       define GARDROUND 0x7f
 216 #       define GARDMASK  0xff
 217 #       define GARDMSB   0x80
 218 #       define EXPBITS 11
 219 #       define EXPBIAS 1023
 220 #       define FRACBITS 52
 221 #       define EXPMAX (0x7ff)
 222 #       define QUIET_NAN 0x8000000000000LL
 223 #       define FRAC_NBITS 64
 224 #       define FRACHIGH  0x8000000000000000LL
 225 #       define FRACHIGH2 0xc000000000000000LL
 226 #       define pack_d __pack_d
 227 #       define unpack_d __unpack_d
 228 #       define __fpcmp_parts __fpcmp_parts_d
 229         typedef UDItype fractype;
 230         typedef USItype halffractype;
 231         typedef DFtype FLO_type;
 232         typedef DItype intfrac;
 233 #endif
 234
 235 #ifdef US_SOFTWARE_GOFAST
 236 #       ifdef FLOAT
 237 #               define add              fpadd
 238 #               define sub              fpsub
 239 #               define multiply         fpmul
 240 #               define divide           fpdiv
 241 #               define compare          fpcmp
 242 #               define si_to_float      sitofp
 243 #               define float_to_si      fptosi
 244 #               define float_to_usi     fptoui
 245 #               define negate           __negsf2
 246 #               define sf_to_df         fptodp
 247 #               define dptofp           dptofp
 248 #else
 249 #               define add              dpadd
 250 #               define sub              dpsub
 251 #               define multiply         dpmul
 252 #               define divide           dpdiv
 253 #               define compare          dpcmp
 254 #               define si_to_float      litodp
 255 #               define float_to_si      dptoli
 256 #               define float_to_usi     dptoul
 257 #               define negate           __negdf2
 258 #               define df_to_sf         dptofp
 259 #endif
 260 #else
 261 #       ifdef FLOAT
 262 #               define add              __addsf3
 263 #               define sub              __subsf3
 264 #               define multiply         __mulsf3
 265 #               define divide           __divsf3
 266 #               define compare          __cmpsf2
 267 #               define _eq_f2           __eqsf2
 268 #               define _ne_f2           __nesf2
 269 #               define _gt_f2           __gtsf2
 270 #               define _ge_f2           __gesf2
 271 #               define _lt_f2           __ltsf2
 272 #               define _le_f2           __lesf2
 273 #               define _unord_f2        __unordsf2
 274 #               define si_to_float      __floatsisf
 275 #               define float_to_si      __fixsfsi
 276 #               define float_to_usi     __fixunssfsi
 277 #               define negate           __negsf2
 278 #               define sf_to_df         __extendsfdf2
 279 #else
 280 #               define add              __adddf3
 281 #               define sub              __subdf3
 282 #               define multiply         __muldf3
 283 #               define divide           __divdf3
 284 #               define compare          __cmpdf2
 285 #               define _eq_f2           __eqdf2
 286 #               define _ne_f2           __nedf2
 287 #               define _gt_f2           __gtdf2
 288 #               define _ge_f2           __gedf2
 289 #               define _lt_f2           __ltdf2
 290 #               define _le_f2           __ledf2
 291 #               define _unord_f2        __unorddf2
 292 #               define si_to_float      __floatsidf
 293 #               define float_to_si      __fixdfsi
 294 #               define float_to_usi     __fixunsdfsi
 295 #               define negate           __negdf2
 296 #               define df_to_sf         __truncdfsf2
 297 #       endif
 298 #endif
 299
 300
 301 #ifndef INLINE
 302 #define INLINE __inline__
 303 #endif
 304
 305 /* Preserve the sticky-bit when shifting fractions to the right.  */
 306 #define LSHIFT(a) { a = (a & 1) | (a >> 1); }
 307
 308 /* numeric parameters */
 309 /* F_D_BITOFF is the number of bits offset between the MSB of the mantissa
 310    of a float and of a double. Assumes there are only two float types.
 311    (double::FRAC_BITS+double::NGARDS-(float::FRAC_BITS-float::NGARDS))
 312  */
 313 #define F_D_BITOFF (52+8-(23+7))
 314
 315
 316 #define NORMAL_EXPMIN (-(EXPBIAS)+1)
 317 #define IMPLICIT_1 (1LL<<(FRACBITS+NGARDS))
 318 #define IMPLICIT_2 (1LL<<(FRACBITS+1+NGARDS))
 319
 320 /* common types */
 321
 322 typedef enum
 323 {
 324   CLASS_SNAN,
 325   CLASS_QNAN,
 326   CLASS_ZERO,
 327   CLASS_NUMBER,
 328   CLASS_INFINITY
 329 } fp_class_type;
 330
 331 typedef struct
 332 {
 333 #ifdef SMALL_MACHINE
 334   char class;
 335   unsigned char sign;
 336   short normal_exp;
 337 #else
 338   fp_class_type class;
 339   unsigned int sign;
 340   int normal_exp;
 341 #endif
 342
 343   union
 344     {
 345       fractype ll;
 346       halffractype l[2];
 347     } fraction;
 348 } fp_number_type;
 349
 350 typedef union
 351 {
 352   FLO_type value;
 353   fractype value_raw;
 354
 355 #ifndef FLOAT
 356   halffractype words[2];
 357 #endif
 358
 359 #ifdef FLOAT_BIT_ORDER_MISMATCH
 360   struct
 361     {
 362       fractype fraction:FRACBITS __attribute__ ((packed));
 363       unsigned int exp:EXPBITS __attribute__ ((packed));
 364       unsigned int sign:1 __attribute__ ((packed));
 365     }
 366   bits;
 367 #endif
 368
 369 #ifdef _DEBUG_BITFLOAT
 370   struct
 371     {
 372       unsigned int sign:1 __attribute__ ((packed));
 373       unsigned int exp:EXPBITS __attribute__ ((packed));
 374       fractype fraction:FRACBITS __attribute__ ((packed));
 375     }
 376   bits_big_endian;
 377
 378   struct
 379     {
 380       fractype fraction:FRACBITS __attribute__ ((packed));
 381       unsigned int exp:EXPBITS __attribute__ ((packed));
 382       unsigned int sign:1 __attribute__ ((packed));
 383     }
 384   bits_little_endian;
 385 #endif
 386 }
 387 FLO_union_type;
 388
 389
 390 /* end of header */
 391
 392 /* IEEE "special" number predicates */
 393
 394 #ifdef NO_NANS
 395
 396 #define nan() 0
 397 #define isnan(x) 0
 398 #define isinf(x) 0
 399 #else
 400
 401 #if   defined L_thenan_sf
 402 const fp_number_type __thenan_sf = { CLASS_SNAN, 0, 0, (fractype) 0 };
 403 #elif defined L_thenan_df
 404 const fp_number_type __thenan_df = { CLASS_SNAN, 0, 0, (fractype) 0 };
 405 #elif defined FLOAT
 406 extern const fp_number_type __thenan_sf;
 407 #else
 408 extern const fp_number_type __thenan_df;
 409 #endif
 410
 411 INLINE
 412 static fp_number_type *
 413 nan ()
 414 {
 415   /* Discard the const qualifier... */
 416 #ifdef FLOAT
 417   return (fp_number_type *) (& __thenan_sf);
 418 #else
 419   return (fp_number_type *) (& __thenan_df);
 420 #endif
 421 }
 422
 423 INLINE
 424 static int
 425 isnan ( fp_number_type *  x)
 426 {
 427   return x->class == CLASS_SNAN || x->class == CLASS_QNAN;
 428 }
 429
 430 INLINE
 431 static int
 432 isinf ( fp_number_type *  x)
 433 {
 434   return x->class == CLASS_INFINITY;
 435 }
 436
 437 #endif
 438
 439 INLINE
 440 static int
 441 iszero ( fp_number_type *  x)
 442 {
 443   return x->class == CLASS_ZERO;
 444 }
 445
 446 INLINE
 447 static void
 448 flip_sign ( fp_number_type *  x)
 449 {
 450   x->sign = !x->sign;
 451 }
 452
 453 extern FLO_type pack_d ( fp_number_type * );
 454
 455 #if defined(L_pack_df) || defined(L_pack_sf)
 456 FLO_type
 457 pack_d ( fp_number_type *  src)
 458 {
 459   FLO_union_type dst;
 460   fractype fraction = src->fraction.ll; /* wasn't unsigned before? */
 461   int sign = src->sign;
 462   int exp = 0;
 463
 464   if (isnan (src))
 465     {
 466       exp = EXPMAX;
 467       if (src->class == CLASS_QNAN || 1)
 468         {
 469           fraction |= QUIET_NAN;
 470         }
 471     }
 472   else if (isinf (src))
 473     {
 474       exp = EXPMAX;
 475       fraction = 0;
 476     }
 477   else if (iszero (src))
 478     {
 479       exp = 0;
 480       fraction = 0;
 481     }
 482   else if (fraction == 0)
 483     {
 484       exp = 0;
 485     }
 486   else
 487     {
 488       if (src->normal_exp < NORMAL_EXPMIN)
 489         {
 490           /* This number's exponent is too low to fit into the bits
 491              available in the number, so we'll store 0 in the exponent and
 492              shift the fraction to the right to make up for it.  */
 493
 494           int shift = NORMAL_EXPMIN - src->normal_exp;
 495
 496           exp = 0;
 497
 498           if (shift > FRAC_NBITS - NGARDS)
 499             {
 500               /* No point shifting, since it's more that 64 out.  */
 501               fraction = 0;
 502             }
 503           else
 504             {
 505               /* Shift by the value */
 506               fraction >>= shift;
 507             }
 508           fraction >>= NGARDS;
 509         }
 510       else if (src->normal_exp > EXPBIAS)
 511         {
 512           exp = EXPMAX;
 513           fraction = 0;
 514         }
 515       else
 516         {
 517           exp = src->normal_exp + EXPBIAS;
 518           /* IF the gard bits are the all zero, but the first, then we're
 519              half way between two numbers, choose the one which makes the
 520              lsb of the answer 0.  */
 521           if ((fraction & GARDMASK) == GARDMSB)
 522             {
 523               if (fraction & (1 << NGARDS))
 524                 fraction += GARDROUND + 1;
 525             }
 526           else
 527             {
 528               /* Add a one to the guards to round up */
 529               fraction += GARDROUND;
 530             }
 531           if (fraction >= IMPLICIT_2)
 532             {
 533               fraction >>= 1;
 534               exp += 1;
 535             }
 536           fraction >>= NGARDS;
 537         }
 538     }
 539
 540   /* We previously used bitfields to store the number, but this doesn't
 541      handle little/big endian systems conveniently, so use shifts and
 542      masks */
 543 #ifdef FLOAT_BIT_ORDER_MISMATCH
 544   dst.bits.fraction = fraction;
 545   dst.bits.exp = exp;
 546   dst.bits.sign = sign;
 547 #else
 548   dst.value_raw = fraction & ((((fractype)1) << FRACBITS) - (fractype)1);
 549   dst.value_raw |= ((fractype) (exp & ((1 << EXPBITS) - 1))) << FRACBITS;
 550   dst.value_raw |= ((fractype) (sign & 1)) << (FRACBITS | EXPBITS);
 551 #endif
 552
 553 #if defined(FLOAT_WORD_ORDER_MISMATCH) && !defined(FLOAT)
 554   {
 555     halffractype tmp = dst.words[0];
 556     dst.words[0] = dst.words[1];
 557     dst.words[1] = tmp;
 558   }
 559 #endif
 560
 561   return dst.value;
 562 }
 563 #endif
 564
 565 extern void unpack_d (FLO_union_type *, fp_number_type *);
 566
 567 #if defined(L_unpack_df) || defined(L_unpack_sf)
 568 void
 569 unpack_d (FLO_union_type * src, fp_number_type * dst)
 570 {
 571   /* We previously used bitfields to store the number, but this doesn't
 572      handle little/big endian systems conveniently, so use shifts and
 573      masks */
 574   fractype fraction;
 575   int exp;
 576   int sign;
 577
 578 #if defined(FLOAT_WORD_ORDER_MISMATCH) && !defined(FLOAT)
 579   FLO_union_type swapped;
 580
 581   swapped.words[0] = src->words[1];
 582   swapped.words[1] = src->words[0];
 583   src = &swapped;
 584 #endif
 585
 586 #ifdef FLOAT_BIT_ORDER_MISMATCH
 587   fraction = src->bits.fraction;
 588   exp = src->bits.exp;
 589   sign = src->bits.sign;
 590 #else
 591   fraction = src->value_raw & ((((fractype)1) << FRACBITS) - (fractype)1);
 592   exp = ((int)(src->value_raw >> FRACBITS)) & ((1 << EXPBITS) - 1);
 593   sign = ((int)(src->value_raw >> (FRACBITS + EXPBITS))) & 1;
 594 #endif
 595
 596   dst->sign = sign;
 597   if (exp == 0)
 598     {
 599       /* Hmm.  Looks like 0 */
 600       if (fraction == 0)
 601         {
 602           /* tastes like zero */
 603           dst->class = CLASS_ZERO;
 604         }
 605       else
 606         {
 607           /* Zero exponent with non zero fraction - it's denormalized,
 608              so there isn't a leading implicit one - we'll shift it so
 609              it gets one.  */
 610           dst->normal_exp = exp - EXPBIAS + 1;
 611           fraction <<= NGARDS;
 612
 613           dst->class = CLASS_NUMBER;
 614 #if 1
 615           while (fraction < IMPLICIT_1)
 616             {
 617               fraction <<= 1;
 618               dst->normal_exp--;
 619             }
 620 #endif
 621           dst->fraction.ll = fraction;
 622         }
 623     }
 624   else if (exp == EXPMAX)
 625     {
 626       /* Huge exponent*/
 627       if (fraction == 0)
 628         {
 629           /* Attached to a zero fraction - means infinity */
 630           dst->class = CLASS_INFINITY;
 631         }
 632       else
 633         {
 634           /* Non zero fraction, means nan */
 635           if (fraction & QUIET_NAN)
 636             {
 637               dst->class = CLASS_QNAN;
 638             }
 639           else
 640             {
 641               dst->class = CLASS_SNAN;
 642             }
 643           /* Keep the fraction part as the nan number */
 644           dst->fraction.ll = fraction;
 645         }
 646     }
 647   else
 648     {
 649       /* Nothing strange about this number */
 650       dst->normal_exp = exp - EXPBIAS;
 651       dst->class = CLASS_NUMBER;
 652       dst->fraction.ll = (fraction << NGARDS) | IMPLICIT_1;
 653     }
 654 }
 655 #endif
 656
 657 #if defined(L_addsub_sf) || defined(L_addsub_df)
 658 static fp_number_type *
 659 _fpadd_parts (fp_number_type * a,
 660               fp_number_type * b,
 661               fp_number_type * tmp)
 662 {
 663   intfrac tfraction;
 664
 665   /* Put commonly used fields in local variables.  */
 666   int a_normal_exp;
 667   int b_normal_exp;
 668   fractype a_fraction;
 669   fractype b_fraction;
 670
 671   if (isnan (a))
 672     {
 673       return a;
 674     }
 675   if (isnan (b))
 676     {
 677       return b;
 678     }
 679   if (isinf (a))
 680     {
 681       /* Adding infinities with opposite signs yields a NaN.  */
 682       if (isinf (b) && a->sign != b->sign)
 683         return nan ();
 684       return a;
 685     }
 686   if (isinf (b))
 687     {
 688       return b;
 689     }
 690   if (iszero (b))
 691     {
 692       if (iszero (a))
 693         {
 694           *tmp = *a;
 695           tmp->sign = a->sign & b->sign;
 696           return tmp;
 697         }
 698       return a;
 699     }
 700   if (iszero (a))
 701     {
 702       return b;
 703     }
 704
 705   /* Got two numbers. shift the smaller and increment the exponent till
 706      they're the same */
 707   {
 708     int diff;
 709
 710     a_normal_exp = a->normal_exp;
 711     b_normal_exp = b->normal_exp;
 712     a_fraction = a->fraction.ll;
 713     b_fraction = b->fraction.ll;
 714
 715     diff = a_normal_exp - b_normal_exp;
 716
 717     if (diff < 0)
 718       diff = -diff;
 719     if (diff < FRAC_NBITS)
 720       {
 721         /* ??? This does shifts one bit at a time.  Optimize.  */
 722         while (a_normal_exp > b_normal_exp)
 723           {
 724             b_normal_exp++;
 725             LSHIFT (b_fraction);
 726           }
 727         while (b_normal_exp > a_normal_exp)
 728           {
 729             a_normal_exp++;
 730             LSHIFT (a_fraction);
 731           }
 732       }
 733     else
 734       {
 735         /* Somethings's up.. choose the biggest */
 736         if (a_normal_exp > b_normal_exp)
 737           {
 738             b_normal_exp = a_normal_exp;
 739             b_fraction = 0;
 740           }
 741         else
 742           {
 743             a_normal_exp = b_normal_exp;
 744             a_fraction = 0;
 745           }
 746       }
 747   }
 748
 749   if (a->sign != b->sign)
 750     {
 751       if (a->sign)
 752         {
 753           tfraction = -a_fraction + b_fraction;
 754         }
 755       else
 756         {
 757           tfraction = a_fraction - b_fraction;
 758         }
 759       if (tfraction >= 0)
 760         {
 761           tmp->sign = 0;
 762           tmp->normal_exp = a_normal_exp;
 763           tmp->fraction.ll = tfraction;
 764         }
 765       else
 766         {
 767           tmp->sign = 1;
 768           tmp->normal_exp = a_normal_exp;
 769           tmp->fraction.ll = -tfraction;
 770         }
 771       /* and renormalize it */
 772
 773       while (tmp->fraction.ll < IMPLICIT_1 && tmp->fraction.ll)
 774         {
 775           tmp->fraction.ll <<= 1;
 776           tmp->normal_exp--;
 777         }
 778     }
 779   else
 780     {
 781       tmp->sign = a->sign;
 782       tmp->normal_exp = a_normal_exp;
 783       tmp->fraction.ll = a_fraction + b_fraction;
 784     }
 785   tmp->class = CLASS_NUMBER;
 786   /* Now the fraction is added, we have to shift down to renormalize the
 787      number */
 788
 789   if (tmp->fraction.ll >= IMPLICIT_2)
 790     {
 791       LSHIFT (tmp->fraction.ll);
 792       tmp->normal_exp++;
 793     }
 794   return tmp;
 795
 796 }
 797
 798 FLO_type
 799 add (FLO_type arg_a, FLO_type arg_b)
 800 {
 801   fp_number_type a;
 802   fp_number_type b;
 803   fp_number_type tmp;
 804   fp_number_type *res;
 805   FLO_union_type au, bu;
 806
 807   au.value = arg_a;
 808   bu.value = arg_b;
 809
 810   unpack_d (&au, &a);
 811   unpack_d (&bu, &b);
 812
 813   res = _fpadd_parts (&a, &b, &tmp);
 814
 815   return pack_d (res);
 816 }
 817
 818 FLO_type
 819 sub (FLO_type arg_a, FLO_type arg_b)
 820 {
 821   fp_number_type a;
 822   fp_number_type b;
 823   fp_number_type tmp;
 824   fp_number_type *res;
 825   FLO_union_type au, bu;
 826
 827   au.value = arg_a;
 828   bu.value = arg_b;
 829
 830   unpack_d (&au, &a);
 831   unpack_d (&bu, &b);
 832
 833   b.sign ^= 1;
 834
 835   res = _fpadd_parts (&a, &b, &tmp);
 836
 837   return pack_d (res);
 838 }
 839 #endif
 840
 841 #if defined(L_mul_sf) || defined(L_mul_df)
 842 static INLINE fp_number_type *
 843 _fpmul_parts ( fp_number_type *  a,
 844                fp_number_type *  b,
 845                fp_number_type * tmp)
 846 {
 847   fractype low = 0;
 848   fractype high = 0;
 849
 850   if (isnan (a))
 851     {
 852       a->sign = a->sign != b->sign;
 853       return a;
 854     }
 855   if (isnan (b))
 856     {
 857       b->sign = a->sign != b->sign;
 858       return b;
 859     }
 860   if (isinf (a))
 861     {
 862       if (iszero (b))
 863         return nan ();
 864       a->sign = a->sign != b->sign;
 865       return a;
 866     }
 867   if (isinf (b))
 868     {
 869       if (iszero (a))
 870         {
 871           return nan ();
 872         }
 873       b->sign = a->sign != b->sign;
 874       return b;
 875     }
 876   if (iszero (a))
 877     {
 878       a->sign = a->sign != b->sign;
 879       return a;
 880     }
 881   if (iszero (b))
 882     {
 883       b->sign = a->sign != b->sign;
 884       return b;
 885     }
 886
 887   /* Calculate the mantissa by multiplying both 64bit numbers to get a
 888      128 bit number */
 889   {
 890 #if defined(NO_DI_MODE)
 891     {
 892       fractype x = a->fraction.ll;
 893       fractype ylow = b->fraction.ll;
 894       fractype yhigh = 0;
 895       int bit;
 896
 897       /* ??? This does multiplies one bit at a time.  Optimize.  */
 898       for (bit = 0; bit < FRAC_NBITS; bit++)
 899         {
 900           int carry;
 901
 902           if (x & 1)
 903             {
 904               carry = (low += ylow) < ylow;
 905               high += yhigh + carry;
 906             }
 907           yhigh <<= 1;
 908           if (ylow & FRACHIGH)
 909             {
 910               yhigh |= 1;
 911             }
 912           ylow <<= 1;
 913           x >>= 1;
 914         }
 915     }
 916 #elif defined(FLOAT)
 917     {
 918       /* Multiplying two 32 bit numbers to get a 64 bit number  on
 919         a machine with DI, so we're safe */
 920
 921       DItype answer = (DItype)(a->fraction.ll) * (DItype)(b->fraction.ll);
 922
 923       high = answer >> 32;
 924       low = answer;
 925     }
 926 #else
 927     /* Doing a 64*64 to 128 */
 928     {
 929       UDItype nl = a->fraction.ll & 0xffffffff;
 930       UDItype nh = a->fraction.ll >> 32;
 931       UDItype ml = b->fraction.ll & 0xffffffff;
 932       UDItype mh = b->fraction.ll >>32;
 933       UDItype pp_ll = ml * nl;
 934       UDItype pp_hl = mh * nl;
 935       UDItype pp_lh = ml * nh;
 936       UDItype pp_hh = mh * nh;
 937       UDItype res2 = 0;
 938       UDItype res0 = 0;
 939       UDItype ps_hh__ = pp_hl + pp_lh;
 940       if (ps_hh__ < pp_hl)
 941         res2 += 0x100000000LL;
 942       pp_hl = (ps_hh__ << 32) & 0xffffffff00000000LL;
 943       res0 = pp_ll + pp_hl;
 944       if (res0 < pp_ll)
 945         res2++;
 946       res2 += ((ps_hh__ >> 32) & 0xffffffffL) + pp_hh;
 947       high = res2;
 948       low = res0;
 949     }
 950 #endif
 951   }
 952
 953   tmp->normal_exp = a->normal_exp + b->normal_exp;
 954   tmp->sign = a->sign != b->sign;
 955 #ifdef FLOAT
 956   tmp->normal_exp += 2;         /* ??????????????? */
 957 #else
 958   tmp->normal_exp += 4;         /* ??????????????? */
 959 #endif
 960   while (high >= IMPLICIT_2)
 961     {
 962       tmp->normal_exp++;
 963       if (high & 1)
 964         {
 965           low >>= 1;
 966           low |= FRACHIGH;
 967         }
 968       high >>= 1;
 969     }
 970   while (high < IMPLICIT_1)
 971     {
 972       tmp->normal_exp--;
 973
 974       high <<= 1;
 975       if (low & FRACHIGH)
 976         high |= 1;
 977       low <<= 1;
 978     }
 979   /* rounding is tricky. if we only round if it won't make us round later. */
 980 #if 0
 981   if (low & FRACHIGH2)
 982     {
 983       if (((high & GARDMASK) != GARDMSB)
 984           && (((high + 1) & GARDMASK) == GARDMSB))
 985         {
 986           /* don't round, it gets done again later. */
 987         }
 988       else
 989         {
 990           high++;
 991         }
 992     }
 993 #endif
 994   if ((high & GARDMASK) == GARDMSB)
 995     {
 996       if (high & (1 << NGARDS))
 997         {
 998           /* half way, so round to even */
 999           high += GARDROUND + 1;
1000         }
1001       else if (low)
1002         {
1003           /* but we really weren't half way */
1004           high += GARDROUND + 1;
1005         }
1006     }
1007   tmp->fraction.ll = high;
1008   tmp->class = CLASS_NUMBER;
1009   return tmp;
1010 }
1011
1012 FLO_type
1013 multiply (FLO_type arg_a, FLO_type arg_b)
1014 {
1015   fp_number_type a;
1016   fp_number_type b;
1017   fp_number_type tmp;
1018   fp_number_type *res;
1019   FLO_union_type au, bu;
1020
1021   au.value = arg_a;
1022   bu.value = arg_b;
1023
1024   unpack_d (&au, &a);
1025   unpack_d (&bu, &b);
1026
1027   res = _fpmul_parts (&a, &b, &tmp);
1028
1029   return pack_d (res);
1030 }
1031 #endif
1032
1033 #if defined(L_div_sf) || defined(L_div_df)
1034 static INLINE fp_number_type *
1035 _fpdiv_parts (fp_number_type * a,
1036               fp_number_type * b)
1037 {
1038   fractype bit;
1039   fractype numerator;
1040   fractype denominator;
1041   fractype quotient;
1042
1043   if (isnan (a))
1044     {
1045       return a;
1046     }
1047   if (isnan (b))
1048     {
1049       return b;
1050     }
1051
1052   a->sign = a->sign ^ b->sign;
1053
1054   if (isinf (a) || iszero (a))
1055     {
1056       if (a->class == b->class)
1057         return nan ();
1058       return a;
1059     }
1060
1061   if (isinf (b))
1062     {
1063       a->fraction.ll = 0;
1064       a->normal_exp = 0;
1065       return a;
1066     }
1067   if (iszero (b))
1068     {
1069       a->class = CLASS_INFINITY;
1070       return a;
1071     }
1072
1073   /* Calculate the mantissa by multiplying both 64bit numbers to get a
1074      128 bit number */
1075   {
1076     /* quotient =
1077        ( numerator / denominator) * 2^(numerator exponent -  denominator exponent)
1078      */
1079
1080     a->normal_exp = a->normal_exp - b->normal_exp;
1081     numerator = a->fraction.ll;
1082     denominator = b->fraction.ll;
1083
1084     if (numerator < denominator)
1085       {
1086         /* Fraction will be less than 1.0 */
1087         numerator *= 2;
1088         a->normal_exp--;
1089       }
1090     bit = IMPLICIT_1;
1091     quotient = 0;
1092     /* ??? Does divide one bit at a time.  Optimize.  */
1093     while (bit)
1094       {
1095         if (numerator >= denominator)
1096           {
1097             quotient |= bit;
1098             numerator -= denominator;
1099           }
1100         bit >>= 1;
1101         numerator *= 2;
1102       }
1103
1104     if ((quotient & GARDMASK) == GARDMSB)
1105       {
1106         if (quotient & (1 << NGARDS))
1107           {
1108             /* half way, so round to even */
1109             quotient += GARDROUND + 1;
1110           }
1111         else if (numerator)
1112           {
1113             /* but we really weren't half way, more bits exist */
1114             quotient += GARDROUND + 1;
1115           }
1116       }
1117
1118     a->fraction.ll = quotient;
1119     return (a);
1120   }
1121 }
1122
1123 FLO_type
1124 divide (FLO_type arg_a, FLO_type arg_b)
1125 {
1126   fp_number_type a;
1127   fp_number_type b;
1128   fp_number_type *res;
1129   FLO_union_type au, bu;
1130
1131   au.value = arg_a;
1132   bu.value = arg_b;
1133
1134   unpack_d (&au, &a);
1135   unpack_d (&bu, &b);
1136
1137   res = _fpdiv_parts (&a, &b);
1138
1139   return pack_d (res);
1140 }
1141 #endif
1142
1143 int __fpcmp_parts (fp_number_type * a, fp_number_type *b);
1144
1145 #if defined(L_fpcmp_parts_sf) || defined(L_fpcmp_parts_df)
1146 /* according to the demo, fpcmp returns a comparison with 0... thus
1147    a<b -> -1
1148    a==b -> 0
1149    a>b -> +1
1150  */
1151
1152 int
1153 __fpcmp_parts (fp_number_type * a, fp_number_type * b)
1154 {
1155 #if 0
1156   /* either nan -> unordered. Must be checked outside of this routine. */
1157   if (isnan (a) && isnan (b))
1158     {
1159       return 1;                 /* still unordered! */
1160     }
1161 #endif
1162
1163   if (isnan (a) || isnan (b))
1164     {
1165       return 1;                 /* how to indicate unordered compare? */
1166     }
1167   if (isinf (a) && isinf (b))
1168     {
1169       /* +inf > -inf, but +inf != +inf */
1170       /* b    \a| +inf(0)| -inf(1)
1171        ______\+--------+--------
1172        +inf(0)| a==b(0)| a<b(-1)
1173        -------+--------+--------
1174        -inf(1)| a>b(1) | a==b(0)
1175        -------+--------+--------
1176        So since unordered must be non zero, just line up the columns...
1177        */
1178       return b->sign - a->sign;
1179     }
1180   /* but not both... */
1181   if (isinf (a))
1182     {
1183       return a->sign ? -1 : 1;
1184     }
1185   if (isinf (b))
1186     {
1187       return b->sign ? 1 : -1;
1188     }
1189   if (iszero (a) && iszero (b))
1190     {
1191       return 0;
1192     }
1193   if (iszero (a))
1194     {
1195       return b->sign ? 1 : -1;
1196     }
1197   if (iszero (b))
1198     {
1199       return a->sign ? -1 : 1;
1200     }
1201   /* now both are "normal". */
1202   if (a->sign != b->sign)
1203     {
1204       /* opposite signs */
1205       return a->sign ? -1 : 1;
1206     }
1207   /* same sign; exponents? */
1208   if (a->normal_exp > b->normal_exp)
1209     {
1210       return a->sign ? -1 : 1;
1211     }
1212   if (a->normal_exp < b->normal_exp)
1213     {
1214       return a->sign ? 1 : -1;
1215     }
1216   /* same exponents; check size. */
1217   if (a->fraction.ll > b->fraction.ll)
1218     {
1219       return a->sign ? -1 : 1;
1220     }
1221   if (a->fraction.ll < b->fraction.ll)
1222     {
1223       return a->sign ? 1 : -1;
1224     }
1225   /* after all that, they're equal. */
1226   return 0;
1227 }
1228 #endif
1229
1230 #if defined(L_compare_sf) || defined(L_compare_df)
1231 CMPtype
1232 compare (FLO_type arg_a, FLO_type arg_b)
1233 {
1234   fp_number_type a;
1235   fp_number_type b;
1236   FLO_union_type au, bu;
1237
1238   au.value = arg_a;
1239   bu.value = arg_b;
1240
1241   unpack_d (&au, &a);
1242   unpack_d (&bu, &b);
1243
1244   return __fpcmp_parts (&a, &b);
1245 }
1246 #endif
1247
1248 #ifndef US_SOFTWARE_GOFAST
1249
1250 /* These should be optimized for their specific tasks someday.  */
1251
1252 #if defined(L_eq_sf) || defined(L_eq_df)
1253 CMPtype
1254 _eq_f2 (FLO_type arg_a, FLO_type arg_b)
1255 {
1256   fp_number_type a;
1257   fp_number_type b;
1258   FLO_union_type au, bu;
1259
1260   au.value = arg_a;
1261   bu.value = arg_b;
1262
1263   unpack_d (&au, &a);
1264   unpack_d (&bu, &b);
1265
1266   if (isnan (&a) || isnan (&b))
1267     return 1;                   /* false, truth == 0 */
1268
1269   return __fpcmp_parts (&a, &b) ;
1270 }
1271 #endif
1272
1273 #if defined(L_ne_sf) || defined(L_ne_df)
1274 CMPtype
1275 _ne_f2 (FLO_type arg_a, FLO_type arg_b)
1276 {
1277   fp_number_type a;
1278   fp_number_type b;
1279   FLO_union_type au, bu;
1280
1281   au.value = arg_a;
1282   bu.value = arg_b;
1283
1284   unpack_d (&au, &a);
1285   unpack_d (&bu, &b);
1286
1287   if (isnan (&a) || isnan (&b))
1288     return 1;                   /* true, truth != 0 */
1289
1290   return  __fpcmp_parts (&a, &b) ;
1291 }
1292 #endif
1293
1294 #if defined(L_gt_sf) || defined(L_gt_df)
1295 CMPtype
1296 _gt_f2 (FLO_type arg_a, FLO_type arg_b)
1297 {
1298   fp_number_type a;
1299   fp_number_type b;
1300   FLO_union_type au, bu;
1301
1302   au.value = arg_a;
1303   bu.value = arg_b;
1304
1305   unpack_d (&au, &a);
1306   unpack_d (&bu, &b);
1307
1308   if (isnan (&a) || isnan (&b))
1309     return -1;                  /* false, truth > 0 */
1310
1311   return __fpcmp_parts (&a, &b);
1312 }
1313 #endif
1314
1315 #if defined(L_ge_sf) || defined(L_ge_df)
1316 CMPtype
1317 _ge_f2 (FLO_type arg_a, FLO_type arg_b)
1318 {
1319   fp_number_type a;
1320   fp_number_type b;
1321   FLO_union_type au, bu;
1322
1323   au.value = arg_a;
1324   bu.value = arg_b;
1325
1326   unpack_d (&au, &a);
1327   unpack_d (&bu, &b);
1328
1329   if (isnan (&a) || isnan (&b))
1330     return -1;                  /* false, truth >= 0 */
1331   return __fpcmp_parts (&a, &b) ;
1332 }
1333 #endif
1334
1335 #if defined(L_lt_sf) || defined(L_lt_df)
1336 CMPtype
1337 _lt_f2 (FLO_type arg_a, FLO_type arg_b)
1338 {
1339   fp_number_type a;
1340   fp_number_type b;
1341   FLO_union_type au, bu;
1342
1343   au.value = arg_a;
1344   bu.value = arg_b;
1345
1346   unpack_d (&au, &a);
1347   unpack_d (&bu, &b);
1348
1349   if (isnan (&a) || isnan (&b))
1350     return 1;                   /* false, truth < 0 */
1351
1352   return __fpcmp_parts (&a, &b);
1353 }
1354 #endif
1355
1356 #if defined(L_le_sf) || defined(L_le_df)
1357 CMPtype
1358 _le_f2 (FLO_type arg_a, FLO_type arg_b)
1359 {
1360   fp_number_type a;
1361   fp_number_type b;
1362   FLO_union_type au, bu;
1363
1364   au.value = arg_a;
1365   bu.value = arg_b;
1366
1367   unpack_d (&au, &a);
1368   unpack_d (&bu, &b);
1369
1370   if (isnan (&a) || isnan (&b))
1371     return 1;                   /* false, truth <= 0 */
1372
1373   return __fpcmp_parts (&a, &b) ;
1374 }
1375 #endif
1376
1377 #if defined(L_unord_sf) || defined(L_unord_df)
1378 CMPtype
1379 _unord_f2 (FLO_type arg_a, FLO_type arg_b)
1380 {
1381   fp_number_type a;
1382   fp_number_type b;
1383   FLO_union_type au, bu;
1384
1385   au.value = arg_a;
1386   bu.value = arg_b;
1387
1388   unpack_d (&au, &a);
1389   unpack_d (&bu, &b);
1390
1391   return (isnan (&a) || isnan (&b);
1392 }
1393 #endif
1394
1395 #endif /* ! US_SOFTWARE_GOFAST */
1396
1397 #if defined(L_si_to_sf) || defined(L_si_to_df)
1398 FLO_type
1399 si_to_float (SItype arg_a)
1400 {
1401   fp_number_type in;
1402
1403   in.class = CLASS_NUMBER;
1404   in.sign = arg_a < 0;
1405   if (!arg_a)
1406     {
1407       in.class = CLASS_ZERO;
1408     }
1409   else
1410     {
1411       in.normal_exp = FRACBITS + NGARDS;
1412       if (in.sign)
1413         {
1414           /* Special case for minint, since there is no +ve integer
1415              representation for it */
1416           if (arg_a == (SItype) 0x80000000)
1417             {
1418               return -2147483648.0;
1419             }
1420           in.fraction.ll = (-arg_a);
1421         }
1422       else
1423         in.fraction.ll = arg_a;
1424
1425       while (in.fraction.ll < (1LL << (FRACBITS + NGARDS)))
1426         {
1427           in.fraction.ll <<= 1;
1428           in.normal_exp -= 1;
1429         }
1430     }
1431   return pack_d (&in);
1432 }
1433 #endif
1434
1435 #if defined(L_sf_to_si) || defined(L_df_to_si)
1436 SItype
1437 float_to_si (FLO_type arg_a)
1438 {
1439   fp_number_type a;
1440   SItype tmp;
1441   FLO_union_type au;
1442
1443   au.value = arg_a;
1444   unpack_d (&au, &a);
1445
1446   if (iszero (&a))
1447     return 0;
1448   if (isnan (&a))
1449     return 0;
1450   /* get reasonable MAX_SI_INT... */
1451   if (isinf (&a))
1452     return a.sign ? (-MAX_SI_INT)-1 : MAX_SI_INT;
1453   /* it is a number, but a small one */
1454   if (a.normal_exp < 0)
1455     return 0;
1456   if (a.normal_exp > 30)
1457     return a.sign ? (-MAX_SI_INT)-1 : MAX_SI_INT;
1458   tmp = a.fraction.ll >> ((FRACBITS + NGARDS) - a.normal_exp);
1459   return a.sign ? (-tmp) : (tmp);
1460 }
1461 #endif
1462
1463 #if defined(L_sf_to_usi) || defined(L_df_to_usi)
1464 #ifdef US_SOFTWARE_GOFAST
1465 /* While libgcc2.c defines its own __fixunssfsi and __fixunsdfsi routines,
1466    we also define them for GOFAST because the ones in libgcc2.c have the
1467    wrong names and I'd rather define these here and keep GOFAST CYG-LOC's
1468    out of libgcc2.c.  We can't define these here if not GOFAST because then
1469    there'd be duplicate copies.  */
1470
1471 USItype
1472 float_to_usi (FLO_type arg_a)
1473 {
1474   fp_number_type a;
1475   FLO_union_type au;
1476
1477   au.value = arg_a;
1478   unpack_d (&au, &a);
1479
1480   if (iszero (&a))
1481     return 0;
1482   if (isnan (&a))
1483     return 0;
1484   /* it is a negative number */
1485   if (a.sign)
1486     return 0;
1487   /* get reasonable MAX_USI_INT... */
1488   if (isinf (&a))
1489     return MAX_USI_INT;
1490   /* it is a number, but a small one */
1491   if (a.normal_exp < 0)
1492     return 0;
1493   if (a.normal_exp > 31)
1494     return MAX_USI_INT;
1495   else if (a.normal_exp > (FRACBITS + NGARDS))
1496     return a.fraction.ll << (a.normal_exp - (FRACBITS + NGARDS));
1497   else
1498     return a.fraction.ll >> ((FRACBITS + NGARDS) - a.normal_exp);
1499 }
1500 #endif
1501 #endif
1502
1503 #if defined(L_negate_sf) || defined(L_negate_df)
1504 FLO_type
1505 negate (FLO_type arg_a)
1506 {
1507   fp_number_type a;
1508   FLO_union_type au;
1509
1510   au.value = arg_a;
1511   unpack_d (&au, &a);
1512
1513   flip_sign (&a);
1514   return pack_d (&a);
1515 }
1516 #endif
1517
1518 #ifdef FLOAT
1519
1520 #if defined(L_make_sf)
1521 SFtype
1522 __make_fp(fp_class_type class,
1523              unsigned int sign,
1524              int exp,
1525              USItype frac)
1526 {
1527   fp_number_type in;
1528
1529   in.class = class;
1530   in.sign = sign;
1531   in.normal_exp = exp;
1532   in.fraction.ll = frac;
1533   return pack_d (&in);
1534 }
1535 #endif
1536
1537 #ifndef FLOAT_ONLY
1538
1539 /* This enables one to build an fp library that supports float but not double.
1540    Otherwise, we would get an undefined reference to __make_dp.
1541    This is needed for some 8-bit ports that can't handle well values that
1542    are 8-bytes in size, so we just don't support double for them at all.  */
1543
1544 extern DFtype __make_dp (fp_class_type, unsigned int, int, UDItype frac);
1545
1546 #if defined(L_sf_to_df)
1547 DFtype
1548 sf_to_df (SFtype arg_a)
1549 {
1550   fp_number_type in;
1551   FLO_union_type au;
1552
1553   au.value = arg_a;
1554   unpack_d (&au, &in);
1555
1556   return __make_dp (in.class, in.sign, in.normal_exp,
1557                     ((UDItype) in.fraction.ll) << F_D_BITOFF);
1558 }
1559 #endif
1560
1561 #endif
1562 #endif
1563
1564 #ifndef FLOAT
1565
1566 extern SFtype __make_fp (fp_class_type, unsigned int, int, USItype);
1567
1568 #if defined(L_make_df)
1569 DFtype
1570 __make_dp (fp_class_type class, unsigned int sign, int exp, UDItype frac)
1571 {
1572   fp_number_type in;
1573
1574   in.class = class;
1575   in.sign = sign;
1576   in.normal_exp = exp;
1577   in.fraction.ll = frac;
1578   return pack_d (&in);
1579 }
1580 #endif
1581
1582 #if defined(L_df_to_sf)
1583 SFtype
1584 df_to_sf (DFtype arg_a)
1585 {
1586   fp_number_type in;
1587   USItype sffrac;
1588   FLO_union_type au;
1589
1590   au.value = arg_a;
1591   unpack_d (&au, &in);
1592
1593   sffrac = in.fraction.ll >> F_D_BITOFF;
1594
1595   /* We set the lowest guard bit in SFFRAC if we discarded any non
1596      zero bits.  */
1597   if ((in.fraction.ll & (((USItype) 1 << F_D_BITOFF) - 1)) != 0)
1598     sffrac |= 1;
1599
1600   return __make_fp (in.class, in.sign, in.normal_exp, sffrac);
1601 }
1602 #endif
1603
1604 #endif
1605 #endif /* !EXTENDED_FLOAT_STUBS */