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 EXP(fp) (((fp) >> 23L) & 0xFF)
58 #define MANT(fp) (((fp) & 0x7FFFFFL) | HIDDEN)
59 #define PACK(s,e,m) ((s) | ((e) << 23L) | (m))
60
61 /* the following deal with IEEE double-precision numbers */
62 #define EXCESSD 1022L
63 #define HIDDEND (1L << 20L)
64 #define EXPDBITS 11
65 #define EXPDMASK 0x7FFL
66 #define EXPD(fp) (((fp.l.upper) >> 20L) & 0x7FFL)
67 #define SIGND(fp) ((fp.l.upper) & SIGNBIT)
68 #define MANTD(fp) (((((fp.l.upper) & 0xFFFFF) | HIDDEND) << 10) | \
69 (fp.l.lower >> 22))
70 #define MANTDMASK 0xFFFFFL /* mask of upper part */
71
72 /* the following deal with IEEE extended-precision numbers */
73 #define EXCESSX 16382L
74 #define HIDDENX (1L << 31L)
75 #define EXPXBITS 15
76 #define EXPXMASK 0x7FFF
77 #define EXPX(fp) (((fp.l.upper) >> 16) & EXPXMASK)
78 #define SIGNX(fp) ((fp.l.upper) & SIGNBIT)
79 #define MANTXMASK 0x7FFFFFFFL /* mask of upper part */
80
81 union double_long
82 {
83 double d;
84 struct {
85 long upper;
86 unsigned long lower;
87 } l;
88 };
89
90 union float_long {
91 float f;
92 long l;
93 };
94
95 union long_double_long
96 {
97 long double ld;
98 struct
99 {
100 long upper;
101 unsigned long middle;
102 unsigned long lower;
103 } l;
104 };
105
106 #ifndef EXTFLOAT
107
108 int
__unordsf2(float a,float b)109 __unordsf2(float a, float b)
110 {
111 union float_long fl;
112
113 fl.f = a;
114 if (EXP(fl.l) == EXP(~0u) && (MANT(fl.l) & ~HIDDEN) != 0)
115 return 1;
116 fl.f = b;
117 if (EXP(fl.l) == EXP(~0u) && (MANT(fl.l) & ~HIDDEN) != 0)
118 return 1;
119 return 0;
120 }
121
122 int
__unorddf2(double a,double b)123 __unorddf2(double a, double b)
124 {
125 union double_long dl;
126
127 dl.d = a;
128 if (EXPD(dl) == EXPDMASK
129 && ((dl.l.upper & MANTDMASK) != 0 || dl.l.lower != 0))
130 return 1;
131 dl.d = b;
132 if (EXPD(dl) == EXPDMASK
133 && ((dl.l.upper & MANTDMASK) != 0 || dl.l.lower != 0))
134 return 1;
135 return 0;
136 }
137
138 /* convert unsigned int to double */
139 double
__floatunsidf(unsigned long a1)140 __floatunsidf (unsigned long a1)
141 {
142 long exp = 32 + EXCESSD;
143 union double_long dl;
144
145 if (!a1)
146 {
147 dl.l.upper = dl.l.lower = 0;
148 return dl.d;
149 }
150
151 while (a1 < 0x2000000L)
152 {
153 a1 <<= 4;
154 exp -= 4;
155 }
156
157 while (a1 < 0x80000000L)
158 {
159 a1 <<= 1;
160 exp--;
161 }
162
163 /* pack up and go home */
164 dl.l.upper = exp << 20L;
165 dl.l.upper |= (a1 >> 11L) & ~HIDDEND;
166 dl.l.lower = a1 << 21L;
167
168 return dl.d;
169 }
170
171 /* convert int to double */
172 double
__floatsidf(long a1)173 __floatsidf (long a1)
174 {
175 long sign = 0, exp = 31 + EXCESSD;
176 union double_long dl;
177
178 if (!a1)
179 {
180 dl.l.upper = dl.l.lower = 0;
181 return dl.d;
182 }
183
184 if (a1 < 0)
185 {
186 sign = SIGNBIT;
187 a1 = (long)-(unsigned long)a1;
188 if (a1 < 0)
189 {
190 dl.l.upper = SIGNBIT | ((32 + EXCESSD) << 20L);
191 dl.l.lower = 0;
192 return dl.d;
193 }
194 }
195
196 while (a1 < 0x1000000L)
197 {
198 a1 <<= 4;
199 exp -= 4;
200 }
201
202 while (a1 < 0x40000000L)
203 {
204 a1 <<= 1;
205 exp--;
206 }
207
208 /* pack up and go home */
209 dl.l.upper = sign;
210 dl.l.upper |= exp << 20L;
211 dl.l.upper |= (a1 >> 10L) & ~HIDDEND;
212 dl.l.lower = a1 << 22L;
213
214 return dl.d;
215 }
216
217 /* convert unsigned int to float */
218 float
__floatunsisf(unsigned long l)219 __floatunsisf (unsigned long l)
220 {
221 double foo = __floatunsidf (l);
222 return foo;
223 }
224
225 /* convert int to float */
226 float
__floatsisf(long l)227 __floatsisf (long l)
228 {
229 double foo = __floatsidf (l);
230 return foo;
231 }
232
233 /* convert float to double */
234 double
__extendsfdf2(float a1)235 __extendsfdf2 (float a1)
236 {
237 register union float_long fl1;
238 register union double_long dl;
239 register long exp;
240 register long mant;
241
242 fl1.f = a1;
243
244 dl.l.upper = SIGN (fl1.l);
245 if ((fl1.l & ~SIGNBIT) == 0)
246 {
247 dl.l.lower = 0;
248 return dl.d;
249 }
250
251 exp = EXP(fl1.l);
252 mant = MANT (fl1.l) & ~HIDDEN;
253 if (exp == 0)
254 {
255 /* Denormal. */
256 exp = 1;
257 while (!(mant & HIDDEN))
258 {
259 mant <<= 1;
260 exp--;
261 }
262 mant &= ~HIDDEN;
263 }
264 exp = exp - EXCESS + EXCESSD;
265 dl.l.upper |= exp << 20;
266 dl.l.upper |= mant >> 3;
267 dl.l.lower = mant << 29;
268
269 return dl.d;
270 }
271
272 /* convert double to float */
273 float
__truncdfsf2(double a1)274 __truncdfsf2 (double a1)
275 {
276 register long exp;
277 register long mant;
278 register union float_long fl;
279 register union double_long dl1;
280 int sticky;
281 int shift;
282
283 dl1.d = a1;
284
285 if ((dl1.l.upper & ~SIGNBIT) == 0 && !dl1.l.lower)
286 {
287 fl.l = SIGND(dl1);
288 return fl.f;
289 }
290
291 exp = EXPD (dl1) - EXCESSD + EXCESS;
292
293 sticky = dl1.l.lower & ((1 << 22) - 1);
294 mant = MANTD (dl1);
295 /* shift double mantissa 6 bits so we can round */
296 sticky |= mant & ((1 << 6) - 1);
297 mant >>= 6;
298
299 /* Check for underflow and denormals. */
300 if (exp <= 0)
301 {
302 if (exp < -24)
303 {
304 sticky |= mant;
305 mant = 0;
306 }
307 else
308 {
309 sticky |= mant & ((1 << (1 - exp)) - 1);
310 mant >>= 1 - exp;
311 }
312 exp = 0;
313 }
314
315 /* now round */
316 shift = 1;
317 if ((mant & 1) && (sticky || (mant & 2)))
318 {
319 int rounding = exp ? 2 : 1;
320
321 mant += 1;
322
323 /* did the round overflow? */
324 if (mant >= (HIDDEN << rounding))
325 {
326 exp++;
327 shift = rounding;
328 }
329 }
330 /* shift down */
331 mant >>= shift;
332
333 mant &= ~HIDDEN;
334
335 /* pack up and go home */
336 fl.l = PACK (SIGND (dl1), exp, mant);
337 return (fl.f);
338 }
339
340 /* convert double to int */
341 long
__fixdfsi(double a1)342 __fixdfsi (double a1)
343 {
344 register union double_long dl1;
345 register long exp;
346 register long l;
347
348 dl1.d = a1;
349
350 if (!dl1.l.upper && !dl1.l.lower)
351 return 0;
352
353 exp = EXPD (dl1) - EXCESSD - 31;
354 l = MANTD (dl1);
355
356 if (exp > 0)
357 {
358 /* Return largest integer. */
359 return SIGND (dl1) ? 0x80000000L : 0x7fffffffL;
360 }
361
362 if (exp <= -32)
363 return 0;
364
365 /* shift down until exp = 0 */
366 if (exp < 0)
367 l >>= -exp;
368
369 return (SIGND (dl1) ? -l : l);
370 }
371
372 /* convert float to int */
373 long
__fixsfsi(float a1)374 __fixsfsi (float a1)
375 {
376 double foo = a1;
377 return __fixdfsi (foo);
378 }
379
380 #else /* EXTFLOAT */
381
382 /* We do not need these routines for coldfire, as it has no extended
383 float format. */
384 #if !defined (__mcoldfire__)
385
386 /* Primitive extended precision floating point support.
387
388 We assume all numbers are normalized, don't do any rounding, etc. */
389
390 /* Prototypes for the above in case we use them. */
391 double __floatunsidf (unsigned long);
392 double __floatsidf (long);
393 float __floatsisf (long);
394 double __extendsfdf2 (float);
395 float __truncdfsf2 (double);
396 long __fixdfsi (double);
397 long __fixsfsi (float);
398 long __cmpdf2 (double, double);
399
400 int
__unordxf2(long double a,long double b)401 __unordxf2(long double a, long double b)
402 {
403 union long_double_long ldl;
404
405 ldl.ld = a;
406 if (EXPX(ldl) == EXPXMASK
407 && ((ldl.l.middle & MANTXMASK) != 0 || ldl.l.lower != 0))
408 return 1;
409 ldl.ld = b;
410 if (EXPX(ldl) == EXPXMASK
411 && ((ldl.l.middle & MANTXMASK) != 0 || ldl.l.lower != 0))
412 return 1;
413 return 0;
414 }
415
416 /* convert double to long double */
417 long double
__extenddfxf2(double d)418 __extenddfxf2 (double d)
419 {
420 register union double_long dl;
421 register union long_double_long ldl;
422 register long exp;
423
424 dl.d = d;
425 /*printf ("dfxf in: %g\n", d);*/
426
427 ldl.l.upper = SIGND (dl);
428 if ((dl.l.upper & ~SIGNBIT) == 0 && !dl.l.lower)
429 {
430 ldl.l.middle = 0;
431 ldl.l.lower = 0;
432 return ldl.ld;
433 }
434
435 exp = EXPD (dl) - EXCESSD + EXCESSX;
436 ldl.l.upper |= exp << 16;
437 ldl.l.middle = HIDDENX;
438 /* 31-20: # mantissa bits in ldl.l.middle - # mantissa bits in dl.l.upper */
439 ldl.l.middle |= (dl.l.upper & MANTDMASK) << (31 - 20);
440 /* 1+20: explicit-integer-bit + # mantissa bits in dl.l.upper */
441 ldl.l.middle |= dl.l.lower >> (1 + 20);
442 /* 32 - 21: # bits of dl.l.lower in ldl.l.middle */
443 ldl.l.lower = dl.l.lower << (32 - 21);
444
445 /*printf ("dfxf out: %s\n", dumpxf (ldl.ld));*/
446 return ldl.ld;
447 }
448
449 /* convert long double to double */
450 double
__truncxfdf2(long double ld)451 __truncxfdf2 (long double ld)
452 {
453 register long exp;
454 register union double_long dl;
455 register union long_double_long ldl;
456
457 ldl.ld = ld;
458 /*printf ("xfdf in: %s\n", dumpxf (ld));*/
459
460 dl.l.upper = SIGNX (ldl);
461 if ((ldl.l.upper & ~SIGNBIT) == 0 && !ldl.l.middle && !ldl.l.lower)
462 {
463 dl.l.lower = 0;
464 return dl.d;
465 }
466
467 exp = EXPX (ldl) - EXCESSX + EXCESSD;
468 /* ??? quick and dirty: keep `exp' sane */
469 if (exp >= EXPDMASK)
470 exp = EXPDMASK - 1;
471 dl.l.upper |= exp << (32 - (EXPDBITS + 1));
472 /* +1-1: add one for sign bit, but take one off for explicit-integer-bit */
473 dl.l.upper |= (ldl.l.middle & MANTXMASK) >> (EXPDBITS + 1 - 1);
474 dl.l.lower = (ldl.l.middle & MANTXMASK) << (32 - (EXPDBITS + 1 - 1));
475 dl.l.lower |= ldl.l.lower >> (EXPDBITS + 1 - 1);
476
477 /*printf ("xfdf out: %g\n", dl.d);*/
478 return dl.d;
479 }
480
481 /* convert a float to a long double */
482 long double
__extendsfxf2(float f)483 __extendsfxf2 (float f)
484 {
485 long double foo = __extenddfxf2 (__extendsfdf2 (f));
486 return foo;
487 }
488
489 /* convert a long double to a float */
490 float
__truncxfsf2(long double ld)491 __truncxfsf2 (long double ld)
492 {
493 float foo = __truncdfsf2 (__truncxfdf2 (ld));
494 return foo;
495 }
496
497 /* convert an int to a long double */
498 long double
__floatsixf(long l)499 __floatsixf (long l)
500 {
501 double foo = __floatsidf (l);
502 return foo;
503 }
504
505 /* convert an unsigned int to a long double */
506 long double
__floatunsixf(unsigned long l)507 __floatunsixf (unsigned long l)
508 {
509 double foo = __floatunsidf (l);
510 return foo;
511 }
512
513 /* convert a long double to an int */
514 long
__fixxfsi(long double ld)515 __fixxfsi (long double ld)
516 {
517 long foo = __fixdfsi ((double) ld);
518 return foo;
519 }
520
521 /* The remaining provide crude math support by working in double precision. */
522
523 long double
__addxf3(long double x1,long double x2)524 __addxf3 (long double x1, long double x2)
525 {
526 return (double) x1 + (double) x2;
527 }
528
529 long double
__subxf3(long double x1,long double x2)530 __subxf3 (long double x1, long double x2)
531 {
532 return (double) x1 - (double) x2;
533 }
534
535 long double
__mulxf3(long double x1,long double x2)536 __mulxf3 (long double x1, long double x2)
537 {
538 return (double) x1 * (double) x2;
539 }
540
541 long double
__divxf3(long double x1,long double x2)542 __divxf3 (long double x1, long double x2)
543 {
544 return (double) x1 / (double) x2;
545 }
546
547 long double
__negxf2(long double x1)548 __negxf2 (long double x1)
549 {
550 return - (double) x1;
551 }
552
553 long
__cmpxf2(long double x1,long double x2)554 __cmpxf2 (long double x1, long double x2)
555 {
556 return __cmpdf2 ((double) x1, (double) x2);
557 }
558
559 long
__eqxf2(long double x1,long double x2)560 __eqxf2 (long double x1, long double x2)
561 {
562 return __cmpdf2 ((double) x1, (double) x2);
563 }
564
565 long
__nexf2(long double x1,long double x2)566 __nexf2 (long double x1, long double x2)
567 {
568 return __cmpdf2 ((double) x1, (double) x2);
569 }
570
571 long
__ltxf2(long double x1,long double x2)572 __ltxf2 (long double x1, long double x2)
573 {
574 return __cmpdf2 ((double) x1, (double) x2);
575 }
576
577 long
__lexf2(long double x1,long double x2)578 __lexf2 (long double x1, long double x2)
579 {
580 return __cmpdf2 ((double) x1, (double) x2);
581 }
582
583 long
__gtxf2(long double x1,long double x2)584 __gtxf2 (long double x1, long double x2)
585 {
586 return __cmpdf2 ((double) x1, (double) x2);
587 }
588
589 long
__gexf2(long double x1,long double x2)590 __gexf2 (long double x1, long double x2)
591 {
592 return __cmpdf2 ((double) x1, (double) x2);
593 }
594
595 #endif /* !__mcoldfire__ */
596 #endif /* EXTFLOAT */
597