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
5 /* Copyright (C) 1994, 1995, 1996 Free Software Foundation, Inc.
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
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.)
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.
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, 59 Temple Place - Suite 330,
28 Boston, MA 02111-1307, USA. */
30 /* As a special exception, if you link this library with other files,
31 some of which are compiled with GCC, to produce an executable,
32 this library does not by itself cause the resulting executable
33 to be covered by the GNU General Public License.
34 This exception does not however invalidate any other reasons why
35 the executable file might be covered by the GNU General Public License. */
37 /* This implements IEEE 754 format arithmetic, but does not provide a
38 mechanism for setting the rounding mode, or for generating or handling
41 The original code by Steve Chamberlain, hacked by Mark Eichin and Jim
42 Wilson, all of Cygnus Support. */
44 /* The intended way to use this file is to make two copies, add `#define FLOAT'
45 to one copy, then compile both copies and add them to libgcc.a. */
47 /* The following macros can be defined to change the behaviour of this file:
48 FLOAT: Implement a `float', aka SFmode, fp library. If this is not
49 defined, then this file implements a `double', aka DFmode, fp library.
50 FLOAT_ONLY: Used with FLOAT, to implement a `float' only library, i.e.
51 don't include float->double conversion which requires the double library.
52 This is useful only for machines which can't support doubles, e.g. some
54 CMPtype: Specify the type that floating point compares should return.
55 This defaults to SItype, aka int.
56 US_SOFTWARE_GOFAST: This makes all entry points use the same names as the
57 US Software goFast library. If this is not defined, the entry points use
58 the same names as libgcc1.c.
59 _DEBUG_BITFLOAT: This makes debugging the code a little easier, by adding
60 two integers to the FLO_union_type.
61 NO_NANS: Disable nan and infinity handling
62 SMALL_MACHINE: Useful when operations on QIs and HIs are faster
65 /* We don't currently support extended floats (long doubles) on machines
66 without hardware to deal with them.
68 These stubs are just to keep the linker from complaining about unresolved
69 references which can be pulled in from libio & libstdc++, even if the
70 user isn't using long doubles. However, they may generate an unresolved
71 external to abort if abort is not used by the function, and the stubs
72 are referenced from within libc, since libgcc goes before and after the
75 #ifdef EXTENDED_FLOAT_STUBS
76 __truncxfsf2 (){ abort(); }
77 __extendsfxf2 (){ abort(); }
78 __addxf3 (){ abort(); }
79 __divxf3 (){ abort(); }
80 __eqxf2 (){ abort(); }
81 __extenddfxf2 (){ abort(); }
82 __gtxf2 (){ abort(); }
83 __lexf2 (){ abort(); }
84 __ltxf2 (){ abort(); }
85 __mulxf3 (){ abort(); }
86 __negxf2 (){ abort(); }
87 __nexf2 (){ abort(); }
88 __subxf3 (){ abort(); }
89 __truncxfdf2 (){ abort(); }
91 __trunctfsf2 (){ abort(); }
92 __extendsftf2 (){ abort(); }
93 __addtf3 (){ abort(); }
94 __divtf3 (){ abort(); }
95 __eqtf2 (){ abort(); }
96 __extenddftf2 (){ abort(); }
97 __gttf2 (){ abort(); }
98 __letf2 (){ abort(); }
99 __lttf2 (){ abort(); }
100 __multf3 (){ abort(); }
101 __negtf2 (){ abort(); }
102 __netf2 (){ abort(); }
103 __subtf3 (){ abort(); }
104 __trunctfdf2 (){ abort(); }
105 #else /* !EXTENDED_FLOAT_STUBS, rest of file */
108 typedef SFtype __attribute__ ((mode (SF)));
109 typedef DFtype __attribute__ ((mode (DF)));
111 typedef int HItype __attribute__ ((mode (HI)));
112 typedef int SItype __attribute__ ((mode (SI)));
113 typedef int DItype __attribute__ ((mode (DI)));
115 /* The type of the result of a fp compare */
117 #define CMPtype SItype
120 typedef unsigned int UHItype __attribute__ ((mode (HI)));
121 typedef unsigned int USItype __attribute__ ((mode (SI)));
122 typedef unsigned int UDItype __attribute__ ((mode (DI)));
124 #define MAX_SI_INT ((SItype) ((unsigned) (~0)>>1))
125 #define MAX_USI_INT ((USItype) ~0)
134 # define GARDROUND 0x3f
135 # define GARDMASK 0x7f
136 # define GARDMSB 0x40
140 # define EXPMAX (0xff)
141 # define QUIET_NAN 0x100000L
142 # define FRAC_NBITS 32
143 # define FRACHIGH 0x80000000L
144 # define FRACHIGH2 0xc0000000L
145 # define pack_d pack_f
146 # define unpack_d unpack_f
147 typedef USItype fractype;
148 typedef UHItype halffractype;
149 typedef SFtype FLO_type;
150 typedef SItype intfrac;
153 # define PREFIXFPDP dp
154 # define PREFIXSFDF df
156 # define GARDROUND 0x7f
157 # define GARDMASK 0xff
158 # define GARDMSB 0x80
160 # define EXPBIAS 1023
162 # define EXPMAX (0x7ff)
163 # define QUIET_NAN 0x8000000000000LL
164 # define FRAC_NBITS 64
165 # define FRACHIGH 0x8000000000000000LL
166 # define FRACHIGH2 0xc000000000000000LL
167 typedef UDItype fractype;
168 typedef USItype halffractype;
169 typedef DFtype FLO_type;
170 typedef DItype intfrac;
173 #ifdef US_SOFTWARE_GOFAST
177 # define multiply fpmul
178 # define divide fpdiv
179 # define compare fpcmp
180 # define si_to_float sitofp
181 # define float_to_si fptosi
182 # define float_to_usi fptoui
183 # define negate __negsf2
184 # define sf_to_df fptodp
185 # define dptofp dptofp
189 # define multiply dpmul
190 # define divide dpdiv
191 # define compare dpcmp
192 # define si_to_float litodp
193 # define float_to_si dptoli
194 # define float_to_usi dptoul
195 # define negate __negdf2
196 # define df_to_sf dptofp
200 # define add __addsf3
201 # define sub __subsf3
202 # define multiply __mulsf3
203 # define divide __divsf3
204 # define compare __cmpsf2
205 # define _eq_f2 __eqsf2
206 # define _ne_f2 __nesf2
207 # define _gt_f2 __gtsf2
208 # define _ge_f2 __gesf2
209 # define _lt_f2 __ltsf2
210 # define _le_f2 __lesf2
211 # define si_to_float __floatsisf
212 # define float_to_si __fixsfsi
213 # define float_to_usi __fixunssfsi
214 # define negate __negsf2
215 # define sf_to_df __extendsfdf2
217 # define add __adddf3
218 # define sub __subdf3
219 # define multiply __muldf3
220 # define divide __divdf3
221 # define compare __cmpdf2
222 # define _eq_f2 __eqdf2
223 # define _ne_f2 __nedf2
224 # define _gt_f2 __gtdf2
225 # define _ge_f2 __gedf2
226 # define _lt_f2 __ltdf2
227 # define _le_f2 __ledf2
228 # define si_to_float __floatsidf
229 # define float_to_si __fixdfsi
230 # define float_to_usi __fixunsdfsi
231 # define negate __negdf2
232 # define df_to_sf __truncdfsf2
237 #define INLINE __inline__
239 /* Preserve the sticky-bit when shifting fractions to the right. */
240 #define LSHIFT(a) { a = (a & 1) | (a >> 1); }
242 /* numeric parameters */
243 /* F_D_BITOFF is the number of bits offset between the MSB of the mantissa
244 of a float and of a double. Assumes there are only two float types.
245 (double::FRAC_BITS+double::NGARGS-(float::FRAC_BITS-float::NGARDS))
247 #define F_D_BITOFF (52+8-(23+7))
250 #define NORMAL_EXPMIN (-(EXPBIAS)+1)
251 #define IMPLICIT_1 (1LL<<(FRACBITS+NGARDS))
252 #define IMPLICIT_2 (1LL<<(FRACBITS+1+NGARDS))
290 halffractype words[2];
293 #ifdef FLOAT_BIT_ORDER_MISMATCH
296 fractype fraction:FRACBITS __attribute__ ((packed));
297 unsigned int exp:EXPBITS __attribute__ ((packed));
298 unsigned int sign:1 __attribute__ ((packed));
303 #ifdef _DEBUG_BITFLOAT
306 unsigned int sign:1 __attribute__ ((packed));
307 unsigned int exp:EXPBITS __attribute__ ((packed));
308 fractype fraction:FRACBITS __attribute__ ((packed));
314 fractype fraction:FRACBITS __attribute__ ((packed));
315 unsigned int exp:EXPBITS __attribute__ ((packed));
316 unsigned int sign:1 __attribute__ ((packed));
326 /* IEEE "special" number predicates */
336 static fp_number_type *
339 static fp_number_type thenan;
346 isnan ( fp_number_type * x)
348 return x->class == CLASS_SNAN || x->class == CLASS_QNAN;
353 isinf ( fp_number_type * x)
355 return x->class == CLASS_INFINITY;
362 iszero ( fp_number_type * x)
364 return x->class == CLASS_ZERO;
369 flip_sign ( fp_number_type * x)
375 pack_d ( fp_number_type * src)
378 fractype fraction = src->fraction.ll; /* wasn't unsigned before? */
379 int sign = src->sign;
385 if (src->class == CLASS_QNAN || 1)
387 fraction |= QUIET_NAN;
390 else if (isinf (src))
395 else if (iszero (src))
400 else if (fraction == 0)
407 if (src->normal_exp < NORMAL_EXPMIN)
409 /* This number's exponent is too low to fit into the bits
410 available in the number, so we'll store 0 in the exponent and
411 shift the fraction to the right to make up for it. */
413 int shift = NORMAL_EXPMIN - src->normal_exp;
417 if (shift > FRAC_NBITS - NGARDS)
419 /* No point shifting, since it's more that 64 out. */
424 /* Shift by the value */
429 else if (src->normal_exp > EXPBIAS)
436 exp = src->normal_exp + EXPBIAS;
437 /* IF the gard bits are the all zero, but the first, then we're
438 half way between two numbers, choose the one which makes the
439 lsb of the answer 0. */
440 if ((fraction & GARDMASK) == GARDMSB)
442 if (fraction & (1 << NGARDS))
443 fraction += GARDROUND + 1;
447 /* Add a one to the guards to round up */
448 fraction += GARDROUND;
450 if (fraction >= IMPLICIT_2)
459 /* We previously used bitfields to store the number, but this doesn't
460 handle little/big endian systems conviently, so use shifts and
462 #ifdef FLOAT_BIT_ORDER_MISMATCH
463 dst.bits.fraction = fraction;
465 dst.bits.sign = sign;
467 dst.value_raw = fraction & ((((fractype)1) << FRACBITS) - (fractype)1);
468 dst.value_raw |= ((fractype) (exp & ((1 << EXPBITS) - 1))) << FRACBITS;
469 dst.value_raw |= ((fractype) (sign & 1)) << (FRACBITS | EXPBITS);
472 #if defined(FLOAT_WORD_ORDER_MISMATCH) && !defined(FLOAT)
474 halffractype tmp = dst.words[0];
475 dst.words[0] = dst.words[1];
484 unpack_d (FLO_union_type * src, fp_number_type * dst)
486 /* We previously used bitfields to store the number, but this doesn't
487 handle little/big endian systems conviently, so use shifts and
493 #if defined(FLOAT_WORD_ORDER_MISMATCH) && !defined(FLOAT)
494 FLO_union_type swapped;
496 swapped.words[0] = src->words[1];
497 swapped.words[1] = src->words[0];
501 #ifdef FLOAT_BIT_ORDER_MISMATCH
502 fraction = src->bits.fraction;
504 sign = src->bits.sign;
506 fraction = src->value_raw & ((((fractype)1) << FRACBITS) - (fractype)1);
507 exp = ((int)(src->value_raw >> FRACBITS)) & ((1 << EXPBITS) - 1);
508 sign = ((int)(src->value_raw >> (FRACBITS + EXPBITS))) & 1;
514 /* Hmm. Looks like 0 */
517 /* tastes like zero */
518 dst->class = CLASS_ZERO;
522 /* Zero exponent with non zero fraction - it's denormalized,
523 so there isn't a leading implicit one - we'll shift it so
525 dst->normal_exp = exp - EXPBIAS + 1;
528 dst->class = CLASS_NUMBER;
530 while (fraction < IMPLICIT_1)
536 dst->fraction.ll = fraction;
539 else if (exp == EXPMAX)
544 /* Attached to a zero fraction - means infinity */
545 dst->class = CLASS_INFINITY;
549 /* Non zero fraction, means nan */
552 dst->class = CLASS_SNAN;
556 dst->class = CLASS_QNAN;
558 /* Keep the fraction part as the nan number */
559 dst->fraction.ll = fraction;
564 /* Nothing strange about this number */
565 dst->normal_exp = exp - EXPBIAS;
566 dst->class = CLASS_NUMBER;
567 dst->fraction.ll = (fraction << NGARDS) | IMPLICIT_1;
571 static fp_number_type *
572 _fpadd_parts (fp_number_type * a,
574 fp_number_type * tmp)
578 /* Put commonly used fields in local variables. */
594 /* Adding infinities with opposite signs yields a NaN. */
595 if (isinf (b) && a->sign != b->sign)
612 /* Got two numbers. shift the smaller and increment the exponent till
617 a_normal_exp = a->normal_exp;
618 b_normal_exp = b->normal_exp;
619 a_fraction = a->fraction.ll;
620 b_fraction = b->fraction.ll;
622 diff = a_normal_exp - b_normal_exp;
626 if (diff < FRAC_NBITS)
628 /* ??? This does shifts one bit at a time. Optimize. */
629 while (a_normal_exp > b_normal_exp)
634 while (b_normal_exp > a_normal_exp)
642 /* Somethings's up.. choose the biggest */
643 if (a_normal_exp > b_normal_exp)
645 b_normal_exp = a_normal_exp;
650 a_normal_exp = b_normal_exp;
656 if (a->sign != b->sign)
660 tfraction = -a_fraction + b_fraction;
664 tfraction = a_fraction - b_fraction;
669 tmp->normal_exp = a_normal_exp;
670 tmp->fraction.ll = tfraction;
675 tmp->normal_exp = a_normal_exp;
676 tmp->fraction.ll = -tfraction;
678 /* and renormalize it */
680 while (tmp->fraction.ll < IMPLICIT_1 && tmp->fraction.ll)
682 tmp->fraction.ll <<= 1;
689 tmp->normal_exp = a_normal_exp;
690 tmp->fraction.ll = a_fraction + b_fraction;
692 tmp->class = CLASS_NUMBER;
693 /* Now the fraction is added, we have to shift down to renormalize the
696 if (tmp->fraction.ll >= IMPLICIT_2)
698 LSHIFT (tmp->fraction.ll);
706 add (FLO_type arg_a, FLO_type arg_b)
713 unpack_d ((FLO_union_type *) & arg_a, &a);
714 unpack_d ((FLO_union_type *) & arg_b, &b);
716 res = _fpadd_parts (&a, &b, &tmp);
722 sub (FLO_type arg_a, FLO_type arg_b)
729 unpack_d ((FLO_union_type *) & arg_a, &a);
730 unpack_d ((FLO_union_type *) & arg_b, &b);
734 res = _fpadd_parts (&a, &b, &tmp);
739 static fp_number_type *
740 _fpmul_parts ( fp_number_type * a,
742 fp_number_type * tmp)
749 a->sign = a->sign != b->sign;
754 b->sign = a->sign != b->sign;
761 a->sign = a->sign != b->sign;
770 b->sign = a->sign != b->sign;
775 a->sign = a->sign != b->sign;
780 b->sign = a->sign != b->sign;
784 /* Calculate the mantissa by multiplying both 64bit numbers to get a
787 fractype x = a->fraction.ll;
788 fractype ylow = b->fraction.ll;
792 #if defined(NO_DI_MODE)
794 /* ??? This does multiplies one bit at a time. Optimize. */
795 for (bit = 0; bit < FRAC_NBITS; bit++)
801 carry = (low += ylow) < ylow;
802 high += yhigh + carry;
815 /* Multiplying two 32 bit numbers to get a 64 bit number on
816 a machine with DI, so we're safe */
818 DItype answer = (DItype)(a->fraction.ll) * (DItype)(b->fraction.ll);
824 /* Doing a 64*64 to 128 */
826 UDItype nl = a->fraction.ll & 0xffffffff;
827 UDItype nh = a->fraction.ll >> 32;
828 UDItype ml = b->fraction.ll & 0xffffffff;
829 UDItype mh = b->fraction.ll >>32;
830 UDItype pp_ll = ml * nl;
831 UDItype pp_hl = mh * nl;
832 UDItype pp_lh = ml * nh;
833 UDItype pp_hh = mh * nh;
836 UDItype ps_hh__ = pp_hl + pp_lh;
838 res2 += 0x100000000LL;
839 pp_hl = (ps_hh__ << 32) & 0xffffffff00000000LL;
840 res0 = pp_ll + pp_hl;
843 res2 += ((ps_hh__ >> 32) & 0xffffffffL) + pp_hh;
850 tmp->normal_exp = a->normal_exp + b->normal_exp;
851 tmp->sign = a->sign != b->sign;
853 tmp->normal_exp += 2; /* ??????????????? */
855 tmp->normal_exp += 4; /* ??????????????? */
857 while (high >= IMPLICIT_2)
867 while (high < IMPLICIT_1)
876 /* rounding is tricky. if we only round if it won't make us round later. */
880 if (((high & GARDMASK) != GARDMSB)
881 && (((high + 1) & GARDMASK) == GARDMSB))
883 /* don't round, it gets done again later. */
891 if ((high & GARDMASK) == GARDMSB)
893 if (high & (1 << NGARDS))
895 /* half way, so round to even */
896 high += GARDROUND + 1;
900 /* but we really weren't half way */
901 high += GARDROUND + 1;
904 tmp->fraction.ll = high;
905 tmp->class = CLASS_NUMBER;
910 multiply (FLO_type arg_a, FLO_type arg_b)
917 unpack_d ((FLO_union_type *) & arg_a, &a);
918 unpack_d ((FLO_union_type *) & arg_b, &b);
920 res = _fpmul_parts (&a, &b, &tmp);
925 static fp_number_type *
926 _fpdiv_parts (fp_number_type * a,
928 fp_number_type * tmp)
932 fractype r0, r1, y0, y1, bit;
935 fractype denominator;
948 a->sign = a->sign ^ b->sign;
950 if (isinf (a) || iszero (a))
952 if (a->class == b->class)
965 a->class = CLASS_INFINITY;
969 /* Calculate the mantissa by multiplying both 64bit numbers to get a
973 intfrac d0, d1; /* weren't unsigned before ??? */
976 ( numerator / denominator) * 2^(numerator exponent - denominator exponent)
979 a->normal_exp = a->normal_exp - b->normal_exp;
980 numerator = a->fraction.ll;
981 denominator = b->fraction.ll;
983 if (numerator < denominator)
985 /* Fraction will be less than 1.0 */
991 /* ??? Does divide one bit at a time. Optimize. */
994 if (numerator >= denominator)
997 numerator -= denominator;
1003 if ((quotient & GARDMASK) == GARDMSB)
1005 if (quotient & (1 << NGARDS))
1007 /* half way, so round to even */
1008 quotient += GARDROUND + 1;
1012 /* but we really weren't half way, more bits exist */
1013 quotient += GARDROUND + 1;
1017 a->fraction.ll = quotient;
1023 divide (FLO_type arg_a, FLO_type arg_b)
1028 fp_number_type *res;
1030 unpack_d ((FLO_union_type *) & arg_a, &a);
1031 unpack_d ((FLO_union_type *) & arg_b, &b);
1033 res = _fpdiv_parts (&a, &b, &tmp);
1035 return pack_d (res);
1038 /* according to the demo, fpcmp returns a comparison with 0... thus
1045 _fpcmp_parts (fp_number_type * a, fp_number_type * b)
1048 /* either nan -> unordered. Must be checked outside of this routine. */
1049 if (isnan (a) && isnan (b))
1051 return 1; /* still unordered! */
1055 if (isnan (a) || isnan (b))
1057 return 1; /* how to indicate unordered compare? */
1059 if (isinf (a) && isinf (b))
1061 /* +inf > -inf, but +inf != +inf */
1062 /* b \a| +inf(0)| -inf(1)
1063 ______\+--------+--------
1064 +inf(0)| a==b(0)| a<b(-1)
1065 -------+--------+--------
1066 -inf(1)| a>b(1) | a==b(0)
1067 -------+--------+--------
1068 So since unordered must be non zero, just line up the columns...
1070 return b->sign - a->sign;
1072 /* but not both... */
1075 return a->sign ? -1 : 1;
1079 return b->sign ? 1 : -1;
1081 if (iszero (a) && iszero (b))
1087 return b->sign ? 1 : -1;
1091 return a->sign ? -1 : 1;
1093 /* now both are "normal". */
1094 if (a->sign != b->sign)
1096 /* opposite signs */
1097 return a->sign ? -1 : 1;
1099 /* same sign; exponents? */
1100 if (a->normal_exp > b->normal_exp)
1102 return a->sign ? -1 : 1;
1104 if (a->normal_exp < b->normal_exp)
1106 return a->sign ? 1 : -1;
1108 /* same exponents; check size. */
1109 if (a->fraction.ll > b->fraction.ll)
1111 return a->sign ? -1 : 1;
1113 if (a->fraction.ll < b->fraction.ll)
1115 return a->sign ? 1 : -1;
1117 /* after all that, they're equal. */
1122 compare (FLO_type arg_a, FLO_type arg_b)
1127 unpack_d ((FLO_union_type *) & arg_a, &a);
1128 unpack_d ((FLO_union_type *) & arg_b, &b);
1130 return _fpcmp_parts (&a, &b);
1133 #ifndef US_SOFTWARE_GOFAST
1135 /* These should be optimized for their specific tasks someday. */
1138 _eq_f2 (FLO_type arg_a, FLO_type arg_b)
1143 unpack_d ((FLO_union_type *) & arg_a, &a);
1144 unpack_d ((FLO_union_type *) & arg_b, &b);
1146 if (isnan (&a) || isnan (&b))
1147 return 1; /* false, truth == 0 */
1149 return _fpcmp_parts (&a, &b) ;
1153 _ne_f2 (FLO_type arg_a, FLO_type arg_b)
1158 unpack_d ((FLO_union_type *) & arg_a, &a);
1159 unpack_d ((FLO_union_type *) & arg_b, &b);
1161 if (isnan (&a) || isnan (&b))
1162 return 1; /* true, truth != 0 */
1164 return _fpcmp_parts (&a, &b) ;
1168 _gt_f2 (FLO_type arg_a, FLO_type arg_b)
1173 unpack_d ((FLO_union_type *) & arg_a, &a);
1174 unpack_d ((FLO_union_type *) & arg_b, &b);
1176 if (isnan (&a) || isnan (&b))
1177 return -1; /* false, truth > 0 */
1179 return _fpcmp_parts (&a, &b);
1183 _ge_f2 (FLO_type arg_a, FLO_type arg_b)
1188 unpack_d ((FLO_union_type *) & arg_a, &a);
1189 unpack_d ((FLO_union_type *) & arg_b, &b);
1191 if (isnan (&a) || isnan (&b))
1192 return -1; /* false, truth >= 0 */
1193 return _fpcmp_parts (&a, &b) ;
1197 _lt_f2 (FLO_type arg_a, FLO_type arg_b)
1202 unpack_d ((FLO_union_type *) & arg_a, &a);
1203 unpack_d ((FLO_union_type *) & arg_b, &b);
1205 if (isnan (&a) || isnan (&b))
1206 return 1; /* false, truth < 0 */
1208 return _fpcmp_parts (&a, &b);
1212 _le_f2 (FLO_type arg_a, FLO_type arg_b)
1217 unpack_d ((FLO_union_type *) & arg_a, &a);
1218 unpack_d ((FLO_union_type *) & arg_b, &b);
1220 if (isnan (&a) || isnan (&b))
1221 return 1; /* false, truth <= 0 */
1223 return _fpcmp_parts (&a, &b) ;
1226 #endif /* ! US_SOFTWARE_GOFAST */
1229 si_to_float (SItype arg_a)
1233 in.class = CLASS_NUMBER;
1234 in.sign = arg_a < 0;
1237 in.class = CLASS_ZERO;
1241 in.normal_exp = FRACBITS + NGARDS;
1244 /* Special case for minint, since there is no +ve integer
1245 representation for it */
1246 if (arg_a == 0x80000000)
1248 return -2147483648.0;
1250 in.fraction.ll = (-arg_a);
1253 in.fraction.ll = arg_a;
1255 while (in.fraction.ll < (1LL << (FRACBITS + NGARDS)))
1257 in.fraction.ll <<= 1;
1261 return pack_d (&in);
1265 float_to_si (FLO_type arg_a)
1270 unpack_d ((FLO_union_type *) & arg_a, &a);
1275 /* get reasonable MAX_SI_INT... */
1277 return a.sign ? MAX_SI_INT : (-MAX_SI_INT)-1;
1278 /* it is a number, but a small one */
1279 if (a.normal_exp < 0)
1281 if (a.normal_exp > 30)
1282 return a.sign ? (-MAX_SI_INT)-1 : MAX_SI_INT;
1283 tmp = a.fraction.ll >> ((FRACBITS + NGARDS) - a.normal_exp);
1284 return a.sign ? (-tmp) : (tmp);
1287 #ifdef US_SOFTWARE_GOFAST
1288 /* While libgcc2.c defines its own __fixunssfsi and __fixunsdfsi routines,
1289 we also define them for GOFAST because the ones in libgcc2.c have the
1290 wrong names and I'd rather define these here and keep GOFAST CYG-LOC's
1291 out of libgcc2.c. We can't define these here if not GOFAST because then
1292 there'd be duplicate copies. */
1295 float_to_usi (FLO_type arg_a)
1299 unpack_d ((FLO_union_type *) & arg_a, &a);
1304 /* get reasonable MAX_USI_INT... */
1306 return a.sign ? MAX_USI_INT : 0;
1307 /* it is a negative number */
1310 /* it is a number, but a small one */
1311 if (a.normal_exp < 0)
1313 if (a.normal_exp > 31)
1315 else if (a.normal_exp > (FRACBITS + NGARDS))
1316 return a.fraction.ll << (a.normal_exp - (FRACBITS + NGARDS));
1318 return a.fraction.ll >> ((FRACBITS + NGARDS) - a.normal_exp);
1323 negate (FLO_type arg_a)
1327 unpack_d ((FLO_union_type *) & arg_a, &a);
1335 __make_fp(fp_class_type class,
1344 in.normal_exp = exp;
1345 in.fraction.ll = frac;
1346 return pack_d (&in);
1351 /* This enables one to build an fp library that supports float but not double.
1352 Otherwise, we would get an undefined reference to __make_dp.
1353 This is needed for some 8-bit ports that can't handle well values that
1354 are 8-bytes in size, so we just don't support double for them at all. */
1356 extern DFtype __make_dp (fp_class_type, unsigned int, int, UDItype frac);
1359 sf_to_df (SFtype arg_a)
1363 unpack_d ((FLO_union_type *) & arg_a, &in);
1364 return __make_dp (in.class, in.sign, in.normal_exp,
1365 ((UDItype) in.fraction.ll) << F_D_BITOFF);
1373 extern SFtype __make_fp (fp_class_type, unsigned int, int, USItype);
1376 __make_dp (fp_class_type class, unsigned int sign, int exp, UDItype frac)
1382 in.normal_exp = exp;
1383 in.fraction.ll = frac;
1384 return pack_d (&in);
1388 df_to_sf (DFtype arg_a)
1393 unpack_d ((FLO_union_type *) & arg_a, &in);
1395 sffrac = in.fraction.ll >> F_D_BITOFF;
1397 /* We set the lowest guard bit in SFFRAC if we discarded any non
1399 if ((in.fraction.ll & (((USItype) 1 << F_D_BITOFF) - 1)) != 0)
1402 return __make_fp (in.class, in.sign, in.normal_exp, sffrac);
1406 #endif /* !EXTENDED_FLOAT_STUBS */