]>
gcc.gnu.org Git - gcc.git/blob - gcc/real.c
082cfd03af130c477df65df0dcb4b044dbbb1bec
1 /* real.c - implementation of REAL_ARITHMETIC, REAL_VALUE_ATOF,
2 and support for XFmode IEEE extended real floating point arithmetic.
3 Copyright (C) 1993, 1994, 1995 Free Software Foundation, Inc.
4 Contributed by Stephen L. Moshier (moshier@world.std.com).
6 This file is part of GNU CC.
8 GNU CC is free software; you can redistribute it and/or modify
9 it under the terms of the GNU General Public License as published by
10 the Free Software Foundation; either version 2, or (at your option)
13 GNU CC is distributed in the hope that it will be useful,
14 but WITHOUT ANY WARRANTY; without even the implied warranty of
15 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16 GNU General Public License for more details.
18 You should have received a copy of the GNU General Public License
19 along with GNU CC; see the file COPYING. If not, write to
20 the Free Software Foundation, 59 Temple Place - Suite 330,
21 Boston, MA 02111-1307, USA. */
32 /* To enable support of XFmode extended real floating point, define
33 LONG_DOUBLE_TYPE_SIZE 96 in the tm.h file (m68k.h or i386.h).
35 To support cross compilation between IEEE, VAX and IBM floating
36 point formats, define REAL_ARITHMETIC in the tm.h file.
38 In either case the machine files (tm.h) must not contain any code
39 that tries to use host floating point arithmetic to convert
40 REAL_VALUE_TYPEs from `double' to `float', pass them to fprintf,
41 etc. In cross-compile situations a REAL_VALUE_TYPE may not
42 be intelligible to the host computer's native arithmetic.
44 The emulator defaults to the host's floating point format so that
45 its decimal conversion functions can be used if desired (see
48 The first part of this file interfaces gcc to a floating point
49 arithmetic suite that was not written with gcc in mind. Avoid
50 changing the low-level arithmetic routines unless you have suitable
51 test programs available. A special version of the PARANOIA floating
52 point arithmetic tester, modified for this purpose, can be found on
53 usc.edu: /pub/C-numanal/ieeetest.zoo. Other tests, and libraries of
54 XFmode and TFmode transcendental functions, can be obtained by ftp from
55 netlib.att.com: netlib/cephes. */
57 /* Type of computer arithmetic.
58 Only one of DEC, IBM, IEEE, or UNK should get defined.
60 `IEEE', when REAL_WORDS_BIG_ENDIAN is non-zero, refers generically
61 to big-endian IEEE floating-point data structure. This definition
62 should work in SFmode `float' type and DFmode `double' type on
63 virtually all big-endian IEEE machines. If LONG_DOUBLE_TYPE_SIZE
64 has been defined to be 96, then IEEE also invokes the particular
65 XFmode (`long double' type) data structure used by the Motorola
66 680x0 series processors.
68 `IEEE', when REAL_WORDS_BIG_ENDIAN is zero, refers generally to
69 little-endian IEEE machines. In this case, if LONG_DOUBLE_TYPE_SIZE
70 has been defined to be 96, then IEEE also invokes the particular
71 XFmode `long double' data structure used by the Intel 80x86 series
74 `DEC' refers specifically to the Digital Equipment Corp PDP-11
75 and VAX floating point data structure. This model currently
76 supports no type wider than DFmode.
78 `IBM' refers specifically to the IBM System/370 and compatible
79 floating point data structure. This model currently supports
80 no type wider than DFmode. The IBM conversions were contributed by
81 frank@atom.ansto.gov.au (Frank Crawford).
83 If LONG_DOUBLE_TYPE_SIZE = 64 (the default, unless tm.h defines it)
84 then `long double' and `double' are both implemented, but they
85 both mean DFmode. In this case, the software floating-point
86 support available here is activated by writing
87 #define REAL_ARITHMETIC
90 The case LONG_DOUBLE_TYPE_SIZE = 128 activates TFmode support
91 and may deactivate XFmode since `long double' is used to refer
94 The macros FLOAT_WORDS_BIG_ENDIAN, HOST_FLOAT_WORDS_BIG_ENDIAN,
95 contributed by Richard Earnshaw <Richard.Earnshaw@cl.cam.ac.uk>,
96 separate the floating point unit's endian-ness from that of
97 the integer addressing. This permits one to define a big-endian
98 FPU on a little-endian machine (e.g., ARM). An extension to
99 BYTES_BIG_ENDIAN may be required for some machines in the future.
100 These optional macros may be defined in tm.h. In real.h, they
101 default to WORDS_BIG_ENDIAN, etc., so there is no need to define
102 them for any normal host or target machine on which the floats
103 and the integers have the same endian-ness. */
106 /* The following converts gcc macros into the ones used by this file. */
108 /* REAL_ARITHMETIC defined means that macros in real.h are
109 defined to call emulator functions. */
110 #ifdef REAL_ARITHMETIC
112 #if TARGET_FLOAT_FORMAT == VAX_FLOAT_FORMAT
113 /* PDP-11, Pro350, VAX: */
115 #else /* it's not VAX */
116 #if TARGET_FLOAT_FORMAT == IBM_FLOAT_FORMAT
117 /* IBM System/370 style */
119 #else /* it's also not an IBM */
120 #if TARGET_FLOAT_FORMAT == IEEE_FLOAT_FORMAT
122 #else /* it's not IEEE either */
123 /* UNKnown arithmetic. We don't support this and can't go on. */
124 unknown arithmetic type
126 #endif /* not IEEE */
130 #define REAL_WORDS_BIG_ENDIAN FLOAT_WORDS_BIG_ENDIAN
133 /* REAL_ARITHMETIC not defined means that the *host's* data
134 structure will be used. It may differ by endian-ness from the
135 target machine's structure and will get its ends swapped
136 accordingly (but not here). Probably only the decimal <-> binary
137 functions in this file will actually be used in this case. */
139 #if HOST_FLOAT_FORMAT == VAX_FLOAT_FORMAT
141 #else /* it's not VAX */
142 #if HOST_FLOAT_FORMAT == IBM_FLOAT_FORMAT
143 /* IBM System/370 style */
145 #else /* it's also not an IBM */
146 #if HOST_FLOAT_FORMAT == IEEE_FLOAT_FORMAT
148 #else /* it's not IEEE either */
149 unknown arithmetic type
151 #endif /* not IEEE */
155 #define REAL_WORDS_BIG_ENDIAN HOST_FLOAT_WORDS_BIG_ENDIAN
157 #endif /* REAL_ARITHMETIC not defined */
159 /* Define INFINITY for support of infinity.
160 Define NANS for support of Not-a-Number's (NaN's). */
161 #if !defined(DEC) && !defined(IBM)
166 /* Support of NaNs requires support of infinity. */
173 /* Find a host integer type that is at least 16 bits wide,
174 and another type at least twice whatever that size is. */
176 #if HOST_BITS_PER_CHAR >= 16
177 #define EMUSHORT char
178 #define EMUSHORT_SIZE HOST_BITS_PER_CHAR
179 #define EMULONG_SIZE (2 * HOST_BITS_PER_CHAR)
181 #if HOST_BITS_PER_SHORT >= 16
182 #define EMUSHORT short
183 #define EMUSHORT_SIZE HOST_BITS_PER_SHORT
184 #define EMULONG_SIZE (2 * HOST_BITS_PER_SHORT)
186 #if HOST_BITS_PER_INT >= 16
188 #define EMUSHORT_SIZE HOST_BITS_PER_INT
189 #define EMULONG_SIZE (2 * HOST_BITS_PER_INT)
191 #if HOST_BITS_PER_LONG >= 16
192 #define EMUSHORT long
193 #define EMUSHORT_SIZE HOST_BITS_PER_LONG
194 #define EMULONG_SIZE (2 * HOST_BITS_PER_LONG)
196 /* You will have to modify this program to have a smaller unit size. */
197 #define EMU_NON_COMPILE
203 #if HOST_BITS_PER_SHORT >= EMULONG_SIZE
204 #define EMULONG short
206 #if HOST_BITS_PER_INT >= EMULONG_SIZE
209 #if HOST_BITS_PER_LONG >= EMULONG_SIZE
212 #if HOST_BITS_PER_LONG_LONG >= EMULONG_SIZE
213 #define EMULONG long long int
215 /* You will have to modify this program to have a smaller unit size. */
216 #define EMU_NON_COMPILE
223 /* The host interface doesn't work if no 16-bit size exists. */
224 #if EMUSHORT_SIZE != 16
225 #define EMU_NON_COMPILE
228 /* OK to continue compilation. */
229 #ifndef EMU_NON_COMPILE
231 /* Construct macros to translate between REAL_VALUE_TYPE and e type.
232 In GET_REAL and PUT_REAL, r and e are pointers.
233 A REAL_VALUE_TYPE is guaranteed to occupy contiguous locations
234 in memory, with no holes. */
236 #if LONG_DOUBLE_TYPE_SIZE == 96
237 /* Number of 16 bit words in external e type format */
239 #define MAXDECEXP 4932
240 #define MINDECEXP -4956
241 #define GET_REAL(r,e) bcopy ((char *) r, (char *) e, 2*NE)
242 #define PUT_REAL(e,r) bcopy ((char *) e, (char *) r, 2*NE)
243 #else /* no XFmode */
244 #if LONG_DOUBLE_TYPE_SIZE == 128
246 #define MAXDECEXP 4932
247 #define MINDECEXP -4977
248 #define GET_REAL(r,e) bcopy ((char *) r, (char *) e, 2*NE)
249 #define PUT_REAL(e,r) bcopy ((char *) e, (char *) r, 2*NE)
252 #define MAXDECEXP 4932
253 #define MINDECEXP -4956
254 #ifdef REAL_ARITHMETIC
255 /* Emulator uses target format internally
256 but host stores it in host endian-ness. */
258 #define GET_REAL(r,e) \
260 if (HOST_FLOAT_WORDS_BIG_ENDIAN == REAL_WORDS_BIG_ENDIAN) \
261 e53toe ((unsigned EMUSHORT*) (r), (e)); \
264 unsigned EMUSHORT w[4]; \
265 w[3] = ((EMUSHORT *) r)[0]; \
266 w[2] = ((EMUSHORT *) r)[1]; \
267 w[1] = ((EMUSHORT *) r)[2]; \
268 w[0] = ((EMUSHORT *) r)[3]; \
273 #define PUT_REAL(e,r) \
275 if (HOST_FLOAT_WORDS_BIG_ENDIAN == REAL_WORDS_BIG_ENDIAN) \
276 etoe53 ((e), (unsigned EMUSHORT *) (r)); \
279 unsigned EMUSHORT w[4]; \
281 *((EMUSHORT *) r) = w[3]; \
282 *((EMUSHORT *) r + 1) = w[2]; \
283 *((EMUSHORT *) r + 2) = w[1]; \
284 *((EMUSHORT *) r + 3) = w[0]; \
288 #else /* not REAL_ARITHMETIC */
290 /* emulator uses host format */
291 #define GET_REAL(r,e) e53toe ((unsigned EMUSHORT *) (r), (e))
292 #define PUT_REAL(e,r) etoe53 ((e), (unsigned EMUSHORT *) (r))
294 #endif /* not REAL_ARITHMETIC */
295 #endif /* not TFmode */
296 #endif /* no XFmode */
299 /* Number of 16 bit words in internal format */
302 /* Array offset to exponent */
305 /* Array offset to high guard word */
308 /* Number of bits of precision */
309 #define NBITS ((NI-4)*16)
311 /* Maximum number of decimal digits in ASCII conversion
314 #define NDEC (NBITS*8/27)
316 /* The exponent of 1.0 */
317 #define EXONE (0x3fff)
319 extern int extra_warnings
;
320 extern unsigned EMUSHORT ezero
[], ehalf
[], eone
[], etwo
[];
321 extern unsigned EMUSHORT elog2
[], esqrt2
[];
323 static void endian
PROTO((unsigned EMUSHORT
*, long *,
325 static void eclear
PROTO((unsigned EMUSHORT
*));
326 static void emov
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
327 static void eabs
PROTO((unsigned EMUSHORT
*));
328 static void eneg
PROTO((unsigned EMUSHORT
*));
329 static int eisneg
PROTO((unsigned EMUSHORT
*));
330 static int eisinf
PROTO((unsigned EMUSHORT
*));
331 static int eisnan
PROTO((unsigned EMUSHORT
*));
332 static void einfin
PROTO((unsigned EMUSHORT
*));
333 static void enan
PROTO((unsigned EMUSHORT
*, int));
334 static void emovi
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
335 static void emovo
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
336 static void ecleaz
PROTO((unsigned EMUSHORT
*));
337 static void ecleazs
PROTO((unsigned EMUSHORT
*));
338 static void emovz
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
339 static void einan
PROTO((unsigned EMUSHORT
*));
340 static int eiisnan
PROTO((unsigned EMUSHORT
*));
341 static int eiisneg
PROTO((unsigned EMUSHORT
*));
342 static void eiinfin
PROTO((unsigned EMUSHORT
*));
343 static int eiisinf
PROTO((unsigned EMUSHORT
*));
344 static int ecmpm
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
345 static void eshdn1
PROTO((unsigned EMUSHORT
*));
346 static void eshup1
PROTO((unsigned EMUSHORT
*));
347 static void eshdn8
PROTO((unsigned EMUSHORT
*));
348 static void eshup8
PROTO((unsigned EMUSHORT
*));
349 static void eshup6
PROTO((unsigned EMUSHORT
*));
350 static void eshdn6
PROTO((unsigned EMUSHORT
*));
351 static void eaddm
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));\f
352 static void esubm
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
353 static void m16m
PROTO((unsigned int, unsigned short *,
355 static int edivm
PROTO((unsigned short *, unsigned short *));
356 static int emulm
PROTO((unsigned short *, unsigned short *));
357 static void emdnorm
PROTO((unsigned EMUSHORT
*, int, int, EMULONG
, int));
358 static void esub
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*,
359 unsigned EMUSHORT
*));
360 static void eadd
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*,
361 unsigned EMUSHORT
*));
362 static void eadd1
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*,
363 unsigned EMUSHORT
*));
364 static void ediv
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*,
365 unsigned EMUSHORT
*));
366 static void emul
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*,
367 unsigned EMUSHORT
*));
368 static void e53toe
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
369 static void e64toe
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
370 static void e113toe
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
371 static void e24toe
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
372 static void etoe113
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
373 static void toe113
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
374 static void etoe64
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
375 static void toe64
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
376 static void etoe53
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
377 static void toe53
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
378 static void etoe24
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
379 static void toe24
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
380 static int ecmp
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
381 static void eround
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
382 static void ltoe
PROTO((HOST_WIDE_INT
*, unsigned EMUSHORT
*));
383 static void ultoe
PROTO((unsigned HOST_WIDE_INT
*, unsigned EMUSHORT
*));
384 static void eifrac
PROTO((unsigned EMUSHORT
*, HOST_WIDE_INT
*,
385 unsigned EMUSHORT
*));
386 static void euifrac
PROTO((unsigned EMUSHORT
*, unsigned HOST_WIDE_INT
*,
387 unsigned EMUSHORT
*));
388 static int eshift
PROTO((unsigned EMUSHORT
*, int));
389 static int enormlz
PROTO((unsigned EMUSHORT
*));
390 static void e24toasc
PROTO((unsigned EMUSHORT
*, char *, int));
391 static void e53toasc
PROTO((unsigned EMUSHORT
*, char *, int));
392 static void e64toasc
PROTO((unsigned EMUSHORT
*, char *, int));
393 static void e113toasc
PROTO((unsigned EMUSHORT
*, char *, int));
394 static void etoasc
PROTO((unsigned EMUSHORT
*, char *, int));
395 static void asctoe24
PROTO((char *, unsigned EMUSHORT
*));
396 static void asctoe53
PROTO((char *, unsigned EMUSHORT
*));
397 static void asctoe64
PROTO((char *, unsigned EMUSHORT
*));
398 static void asctoe113
PROTO((char *, unsigned EMUSHORT
*));
399 static void asctoe
PROTO((char *, unsigned EMUSHORT
*));
400 static void asctoeg
PROTO((char *, unsigned EMUSHORT
*, int));
401 static void efloor
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
402 static void efrexp
PROTO((unsigned EMUSHORT
*, int *,
403 unsigned EMUSHORT
*));
404 static void eldexp
PROTO((unsigned EMUSHORT
*, int, unsigned EMUSHORT
*));
405 static void eremain
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*,
406 unsigned EMUSHORT
*));
407 static void eiremain
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
408 static void mtherr
PROTO((char *, int));
409 static void dectoe
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
410 static void etodec
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
411 static void todec
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
412 static void ibmtoe
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*,
414 static void etoibm
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*,
416 static void toibm
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*,
418 static void make_nan
PROTO((unsigned EMUSHORT
*, int, enum machine_mode
));
419 static void uditoe
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
420 static void ditoe
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
421 static void etoudi
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
422 static void etodi
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
423 static void esqrt
PROTO((unsigned EMUSHORT
*, unsigned EMUSHORT
*));
425 /* Copy 32-bit numbers obtained from array containing 16-bit numbers,
426 swapping ends if required, into output array of longs. The
427 result is normally passed to fprintf by the ASM_OUTPUT_ macros. */
431 unsigned EMUSHORT e
[];
433 enum machine_mode mode
;
437 if (REAL_WORDS_BIG_ENDIAN
)
443 /* Swap halfwords in the fourth long. */
444 th
= (unsigned long) e
[6] & 0xffff;
445 t
= (unsigned long) e
[7] & 0xffff;
451 /* Swap halfwords in the third long. */
452 th
= (unsigned long) e
[4] & 0xffff;
453 t
= (unsigned long) e
[5] & 0xffff;
456 /* fall into the double case */
460 /* swap halfwords in the second word */
461 th
= (unsigned long) e
[2] & 0xffff;
462 t
= (unsigned long) e
[3] & 0xffff;
465 /* fall into the float case */
470 /* swap halfwords in the first word */
471 th
= (unsigned long) e
[0] & 0xffff;
472 t
= (unsigned long) e
[1] & 0xffff;
483 /* Pack the output array without swapping. */
490 /* Pack the fourth long. */
491 th
= (unsigned long) e
[7] & 0xffff;
492 t
= (unsigned long) e
[6] & 0xffff;
498 /* Pack the third long.
499 Each element of the input REAL_VALUE_TYPE array has 16 useful bits
501 th
= (unsigned long) e
[5] & 0xffff;
502 t
= (unsigned long) e
[4] & 0xffff;
505 /* fall into the double case */
509 /* pack the second long */
510 th
= (unsigned long) e
[3] & 0xffff;
511 t
= (unsigned long) e
[2] & 0xffff;
514 /* fall into the float case */
519 /* pack the first long */
520 th
= (unsigned long) e
[1] & 0xffff;
521 t
= (unsigned long) e
[0] & 0xffff;
533 /* This is the implementation of the REAL_ARITHMETIC macro. */
536 earith (value
, icode
, r1
, r2
)
537 REAL_VALUE_TYPE
*value
;
542 unsigned EMUSHORT d1
[NE
], d2
[NE
], v
[NE
];
548 /* Return NaN input back to the caller. */
551 PUT_REAL (d1
, value
);
556 PUT_REAL (d2
, value
);
560 code
= (enum tree_code
) icode
;
568 esub (d2
, d1
, v
); /* d1 - d2 */
576 #ifndef REAL_INFINITY
577 if (ecmp (d2
, ezero
) == 0)
580 enan (v
, eisneg (d1
) ^ eisneg (d2
));
587 ediv (d2
, d1
, v
); /* d1/d2 */
590 case MIN_EXPR
: /* min (d1,d2) */
591 if (ecmp (d1
, d2
) < 0)
597 case MAX_EXPR
: /* max (d1,d2) */
598 if (ecmp (d1
, d2
) > 0)
611 /* Truncate REAL_VALUE_TYPE toward zero to signed HOST_WIDE_INT.
612 implements REAL_VALUE_RNDZINT (x) (etrunci (x)). */
618 unsigned EMUSHORT f
[NE
], g
[NE
];
634 /* Truncate REAL_VALUE_TYPE toward zero to unsigned HOST_WIDE_INT;
635 implements REAL_VALUE_UNSIGNED_RNDZINT (x) (etruncui (x)). */
641 unsigned EMUSHORT f
[NE
], g
[NE
];
643 unsigned HOST_WIDE_INT l
;
657 /* This is the REAL_VALUE_ATOF function. It converts a decimal string to
658 binary, rounding off as indicated by the machine_mode argument. Then it
659 promotes the rounded value to REAL_VALUE_TYPE. */
666 unsigned EMUSHORT tem
[NE
], e
[NE
];
696 /* Expansion of REAL_NEGATE. */
702 unsigned EMUSHORT e
[NE
];
712 /* Round real toward zero to HOST_WIDE_INT;
713 implements REAL_VALUE_FIX (x). */
719 unsigned EMUSHORT f
[NE
], g
[NE
];
726 warning ("conversion from NaN to int");
734 /* Round real toward zero to unsigned HOST_WIDE_INT
735 implements REAL_VALUE_UNSIGNED_FIX (x).
736 Negative input returns zero. */
738 unsigned HOST_WIDE_INT
742 unsigned EMUSHORT f
[NE
], g
[NE
];
743 unsigned HOST_WIDE_INT l
;
749 warning ("conversion from NaN to unsigned int");
758 /* REAL_VALUE_FROM_INT macro. */
761 ereal_from_int (d
, i
, j
)
765 unsigned EMUSHORT df
[NE
], dg
[NE
];
766 HOST_WIDE_INT low
, high
;
774 /* complement and add 1 */
781 eldexp (eone
, HOST_BITS_PER_WIDE_INT
, df
);
782 ultoe ((unsigned HOST_WIDE_INT
*) &high
, dg
);
784 ultoe ((unsigned HOST_WIDE_INT
*) &low
, df
);
792 /* REAL_VALUE_FROM_UNSIGNED_INT macro. */
795 ereal_from_uint (d
, i
, j
)
797 unsigned HOST_WIDE_INT i
, j
;
799 unsigned EMUSHORT df
[NE
], dg
[NE
];
800 unsigned HOST_WIDE_INT low
, high
;
804 eldexp (eone
, HOST_BITS_PER_WIDE_INT
, df
);
813 /* REAL_VALUE_TO_INT macro. */
816 ereal_to_int (low
, high
, rr
)
817 HOST_WIDE_INT
*low
, *high
;
820 unsigned EMUSHORT d
[NE
], df
[NE
], dg
[NE
], dh
[NE
];
827 warning ("conversion from NaN to int");
833 /* convert positive value */
840 eldexp (eone
, HOST_BITS_PER_WIDE_INT
, df
);
841 ediv (df
, d
, dg
); /* dg = d / 2^32 is the high word */
842 euifrac (dg
, (unsigned HOST_WIDE_INT
*) high
, dh
);
843 emul (df
, dh
, dg
); /* fractional part is the low word */
844 euifrac (dg
, (unsigned HOST_WIDE_INT
*)low
, dh
);
847 /* complement and add 1 */
857 /* REAL_VALUE_LDEXP macro. */
864 unsigned EMUSHORT e
[NE
], y
[NE
];
877 /* These routines are conditionally compiled because functions
878 of the same names may be defined in fold-const.c. */
880 #ifdef REAL_ARITHMETIC
882 /* Check for infinity in a REAL_VALUE_TYPE. */
888 unsigned EMUSHORT e
[NE
];
898 /* Check whether a REAL_VALUE_TYPE item is a NaN. */
904 unsigned EMUSHORT e
[NE
];
915 /* Check for a negative REAL_VALUE_TYPE number.
916 This just checks the sign bit, so that -0 counts as negative. */
922 return ereal_isneg (x
);
925 /* Expansion of REAL_VALUE_TRUNCATE.
926 The result is in floating point, rounded to nearest or even. */
929 real_value_truncate (mode
, arg
)
930 enum machine_mode mode
;
933 unsigned EMUSHORT e
[NE
], t
[NE
];
969 /* If an unsupported type was requested, presume that
970 the machine files know something useful to do with
971 the unmodified value. */
980 #endif /* REAL_ARITHMETIC defined */
982 /* Used for debugging--print the value of R in human-readable format
991 REAL_VALUE_TO_DECIMAL (r
, "%.20g", dstr
);
992 fprintf (stderr
, "%s", dstr
);
996 /* The following routines convert REAL_VALUE_TYPE to the various floating
997 point formats that are meaningful to supported computers.
999 The results are returned in 32-bit pieces, each piece stored in a `long'.
1000 This is so they can be printed by statements like
1002 fprintf (file, "%lx, %lx", L[0], L[1]);
1004 that will work on both narrow- and wide-word host computers. */
1006 /* Convert R to a 128-bit long double precision value. The output array L
1007 contains four 32-bit pieces of the result, in the order they would appear
1015 unsigned EMUSHORT e
[NE
];
1019 endian (e
, l
, TFmode
);
1022 /* Convert R to a double extended precision value. The output array L
1023 contains three 32-bit pieces of the result, in the order they would
1024 appear in memory. */
1031 unsigned EMUSHORT e
[NE
];
1035 endian (e
, l
, XFmode
);
1038 /* Convert R to a double precision value. The output array L contains two
1039 32-bit pieces of the result, in the order they would appear in memory. */
1046 unsigned EMUSHORT e
[NE
];
1050 endian (e
, l
, DFmode
);
1053 /* Convert R to a single precision float value stored in the least-significant
1054 bits of a `long'. */
1060 unsigned EMUSHORT e
[NE
];
1065 endian (e
, &l
, SFmode
);
1069 /* Convert X to a decimal ASCII string S for output to an assembly
1070 language file. Note, there is no standard way to spell infinity or
1071 a NaN, so these values may require special treatment in the tm.h
1075 ereal_to_decimal (x
, s
)
1079 unsigned EMUSHORT e
[NE
];
1085 /* Compare X and Y. Return 1 if X > Y, 0 if X == Y, -1 if X < Y,
1086 or -2 if either is a NaN. */
1090 REAL_VALUE_TYPE x
, y
;
1092 unsigned EMUSHORT ex
[NE
], ey
[NE
];
1096 return (ecmp (ex
, ey
));
1099 /* Return 1 if the sign bit of X is set, else return 0. */
1105 unsigned EMUSHORT ex
[NE
];
1108 return (eisneg (ex
));
1111 /* End of REAL_ARITHMETIC interface */
1114 Extended precision IEEE binary floating point arithmetic routines
1116 Numbers are stored in C language as arrays of 16-bit unsigned
1117 short integers. The arguments of the routines are pointers to
1120 External e type data structure, similar to Intel 8087 chip
1121 temporary real format but possibly with a larger significand:
1123 NE-1 significand words (least significant word first,
1124 most significant bit is normally set)
1125 exponent (value = EXONE for 1.0,
1126 top bit is the sign)
1129 Internal exploded e-type data structure of a number (a "word" is 16 bits):
1131 ei[0] sign word (0 for positive, 0xffff for negative)
1132 ei[1] biased exponent (value = EXONE for the number 1.0)
1133 ei[2] high guard word (always zero after normalization)
1135 to ei[NI-2] significand (NI-4 significand words,
1136 most significant word first,
1137 most significant bit is set)
1138 ei[NI-1] low guard word (0x8000 bit is rounding place)
1142 Routines for external format e-type numbers
1144 asctoe (string, e) ASCII string to extended double e type
1145 asctoe64 (string, &d) ASCII string to long double
1146 asctoe53 (string, &d) ASCII string to double
1147 asctoe24 (string, &f) ASCII string to single
1148 asctoeg (string, e, prec) ASCII string to specified precision
1149 e24toe (&f, e) IEEE single precision to e type
1150 e53toe (&d, e) IEEE double precision to e type
1151 e64toe (&d, e) IEEE long double precision to e type
1152 e113toe (&d, e) 128-bit long double precision to e type
1153 eabs (e) absolute value
1154 eadd (a, b, c) c = b + a
1156 ecmp (a, b) Returns 1 if a > b, 0 if a == b,
1157 -1 if a < b, -2 if either a or b is a NaN.
1158 ediv (a, b, c) c = b / a
1159 efloor (a, b) truncate to integer, toward -infinity
1160 efrexp (a, exp, s) extract exponent and significand
1161 eifrac (e, &l, frac) e to HOST_WIDE_INT and e type fraction
1162 euifrac (e, &l, frac) e to unsigned HOST_WIDE_INT and e type fraction
1163 einfin (e) set e to infinity, leaving its sign alone
1164 eldexp (a, n, b) multiply by 2**n
1166 emul (a, b, c) c = b * a
1168 eround (a, b) b = nearest integer value to a
1169 esub (a, b, c) c = b - a
1170 e24toasc (&f, str, n) single to ASCII string, n digits after decimal
1171 e53toasc (&d, str, n) double to ASCII string, n digits after decimal
1172 e64toasc (&d, str, n) 80-bit long double to ASCII string
1173 e113toasc (&d, str, n) 128-bit long double to ASCII string
1174 etoasc (e, str, n) e to ASCII string, n digits after decimal
1175 etoe24 (e, &f) convert e type to IEEE single precision
1176 etoe53 (e, &d) convert e type to IEEE double precision
1177 etoe64 (e, &d) convert e type to IEEE long double precision
1178 ltoe (&l, e) HOST_WIDE_INT to e type
1179 ultoe (&l, e) unsigned HOST_WIDE_INT to e type
1180 eisneg (e) 1 if sign bit of e != 0, else 0
1181 eisinf (e) 1 if e has maximum exponent (non-IEEE)
1182 or is infinite (IEEE)
1183 eisnan (e) 1 if e is a NaN
1186 Routines for internal format exploded e-type numbers
1188 eaddm (ai, bi) add significands, bi = bi + ai
1190 ecleazs (ei) set ei = 0 but leave its sign alone
1191 ecmpm (ai, bi) compare significands, return 1, 0, or -1
1192 edivm (ai, bi) divide significands, bi = bi / ai
1193 emdnorm (ai,l,s,exp) normalize and round off
1194 emovi (a, ai) convert external a to internal ai
1195 emovo (ai, a) convert internal ai to external a
1196 emovz (ai, bi) bi = ai, low guard word of bi = 0
1197 emulm (ai, bi) multiply significands, bi = bi * ai
1198 enormlz (ei) left-justify the significand
1199 eshdn1 (ai) shift significand and guards down 1 bit
1200 eshdn8 (ai) shift down 8 bits
1201 eshdn6 (ai) shift down 16 bits
1202 eshift (ai, n) shift ai n bits up (or down if n < 0)
1203 eshup1 (ai) shift significand and guards up 1 bit
1204 eshup8 (ai) shift up 8 bits
1205 eshup6 (ai) shift up 16 bits
1206 esubm (ai, bi) subtract significands, bi = bi - ai
1207 eiisinf (ai) 1 if infinite
1208 eiisnan (ai) 1 if a NaN
1209 eiisneg (ai) 1 if sign bit of ai != 0, else 0
1210 einan (ai) set ai = NaN
1211 eiinfin (ai) set ai = infinity
1213 The result is always normalized and rounded to NI-4 word precision
1214 after each arithmetic operation.
1216 Exception flags are NOT fully supported.
1218 Signaling NaN's are NOT supported; they are treated the same
1221 Define INFINITY for support of infinity; otherwise a
1222 saturation arithmetic is implemented.
1224 Define NANS for support of Not-a-Number items; otherwise the
1225 arithmetic will never produce a NaN output, and might be confused
1227 If NaN's are supported, the output of `ecmp (a,b)' is -2 if
1228 either a or b is a NaN. This means asking `if (ecmp (a,b) < 0)'
1229 may not be legitimate. Use `if (ecmp (a,b) == -1)' for `less than'
1232 Denormals are always supported here where appropriate (e.g., not
1233 for conversion to DEC numbers). */
1235 /* Definitions for error codes that are passed to the common error handling
1238 For Digital Equipment PDP-11 and VAX computers, certain
1239 IBM systems, and others that use numbers with a 56-bit
1240 significand, the symbol DEC should be defined. In this
1241 mode, most floating point constants are given as arrays
1242 of octal integers to eliminate decimal to binary conversion
1243 errors that might be introduced by the compiler.
1245 For computers, such as IBM PC, that follow the IEEE
1246 Standard for Binary Floating Point Arithmetic (ANSI/IEEE
1247 Std 754-1985), the symbol IEEE should be defined.
1248 These numbers have 53-bit significands. In this mode, constants
1249 are provided as arrays of hexadecimal 16 bit integers.
1250 The endian-ness of generated values is controlled by
1251 REAL_WORDS_BIG_ENDIAN.
1253 To accommodate other types of computer arithmetic, all
1254 constants are also provided in a normal decimal radix
1255 which one can hope are correctly converted to a suitable
1256 format by the available C language compiler. To invoke
1257 this mode, the symbol UNK is defined.
1259 An important difference among these modes is a predefined
1260 set of machine arithmetic constants for each. The numbers
1261 MACHEP (the machine roundoff error), MAXNUM (largest number
1262 represented), and several other parameters are preset by
1263 the configuration symbol. Check the file const.c to
1264 ensure that these values are correct for your computer.
1266 For ANSI C compatibility, define ANSIC equal to 1. Currently
1267 this affects only the atan2 function and others that use it. */
1269 /* Constant definitions for math error conditions. */
1271 #define DOMAIN 1 /* argument domain error */
1272 #define SING 2 /* argument singularity */
1273 #define OVERFLOW 3 /* overflow range error */
1274 #define UNDERFLOW 4 /* underflow range error */
1275 #define TLOSS 5 /* total loss of precision */
1276 #define PLOSS 6 /* partial loss of precision */
1277 #define INVALID 7 /* NaN-producing operation */
1279 /* e type constants used by high precision check routines */
1281 #if LONG_DOUBLE_TYPE_SIZE == 128
1283 unsigned EMUSHORT ezero
[NE
] =
1284 {0x0000, 0x0000, 0x0000, 0x0000,
1285 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000,};
1286 extern unsigned EMUSHORT ezero
[];
1289 unsigned EMUSHORT ehalf
[NE
] =
1290 {0x0000, 0x0000, 0x0000, 0x0000,
1291 0x0000, 0x0000, 0x0000, 0x0000, 0x8000, 0x3ffe,};
1292 extern unsigned EMUSHORT ehalf
[];
1295 unsigned EMUSHORT eone
[NE
] =
1296 {0x0000, 0x0000, 0x0000, 0x0000,
1297 0x0000, 0x0000, 0x0000, 0x0000, 0x8000, 0x3fff,};
1298 extern unsigned EMUSHORT eone
[];
1301 unsigned EMUSHORT etwo
[NE
] =
1302 {0x0000, 0x0000, 0x0000, 0x0000,
1303 0x0000, 0x0000, 0x0000, 0x0000, 0x8000, 0x4000,};
1304 extern unsigned EMUSHORT etwo
[];
1307 unsigned EMUSHORT e32
[NE
] =
1308 {0x0000, 0x0000, 0x0000, 0x0000,
1309 0x0000, 0x0000, 0x0000, 0x0000, 0x8000, 0x4004,};
1310 extern unsigned EMUSHORT e32
[];
1312 /* 6.93147180559945309417232121458176568075500134360255E-1 */
1313 unsigned EMUSHORT elog2
[NE
] =
1314 {0x40f3, 0xf6af, 0x03f2, 0xb398,
1315 0xc9e3, 0x79ab, 0150717, 0013767, 0130562, 0x3ffe,};
1316 extern unsigned EMUSHORT elog2
[];
1318 /* 1.41421356237309504880168872420969807856967187537695E0 */
1319 unsigned EMUSHORT esqrt2
[NE
] =
1320 {0x1d6f, 0xbe9f, 0x754a, 0x89b3,
1321 0x597d, 0x6484, 0174736, 0171463, 0132404, 0x3fff,};
1322 extern unsigned EMUSHORT esqrt2
[];
1324 /* 3.14159265358979323846264338327950288419716939937511E0 */
1325 unsigned EMUSHORT epi
[NE
] =
1326 {0x2902, 0x1cd1, 0x80dc, 0x628b,
1327 0xc4c6, 0xc234, 0020550, 0155242, 0144417, 0040000,};
1328 extern unsigned EMUSHORT epi
[];
1331 /* LONG_DOUBLE_TYPE_SIZE is other than 128 */
1332 unsigned EMUSHORT ezero
[NE
] =
1333 {0, 0000000, 0000000, 0000000, 0000000, 0000000,};
1334 unsigned EMUSHORT ehalf
[NE
] =
1335 {0, 0000000, 0000000, 0000000, 0100000, 0x3ffe,};
1336 unsigned EMUSHORT eone
[NE
] =
1337 {0, 0000000, 0000000, 0000000, 0100000, 0x3fff,};
1338 unsigned EMUSHORT etwo
[NE
] =
1339 {0, 0000000, 0000000, 0000000, 0100000, 0040000,};
1340 unsigned EMUSHORT e32
[NE
] =
1341 {0, 0000000, 0000000, 0000000, 0100000, 0040004,};
1342 unsigned EMUSHORT elog2
[NE
] =
1343 {0xc9e4, 0x79ab, 0150717, 0013767, 0130562, 0x3ffe,};
1344 unsigned EMUSHORT esqrt2
[NE
] =
1345 {0x597e, 0x6484, 0174736, 0171463, 0132404, 0x3fff,};
1346 unsigned EMUSHORT epi
[NE
] =
1347 {0xc4c6, 0xc234, 0020550, 0155242, 0144417, 0040000,};
1350 /* Control register for rounding precision.
1351 This can be set to 113 (if NE=10), 80 (if NE=6), 64, 56, 53, or 24 bits. */
1356 /* Clear out entire e-type number X. */
1360 register unsigned EMUSHORT
*x
;
1364 for (i
= 0; i
< NE
; i
++)
1368 /* Move e-type number from A to B. */
1372 register unsigned EMUSHORT
*a
, *b
;
1376 for (i
= 0; i
< NE
; i
++)
1381 /* Absolute value of e-type X. */
1385 unsigned EMUSHORT x
[];
1387 /* sign is top bit of last word of external format */
1388 x
[NE
- 1] &= 0x7fff;
1391 /* Negate the e-type number X. */
1395 unsigned EMUSHORT x
[];
1398 x
[NE
- 1] ^= 0x8000; /* Toggle the sign bit */
1401 /* Return 1 if sign bit of e-type number X is nonzero, else zero. */
1405 unsigned EMUSHORT x
[];
1408 if (x
[NE
- 1] & 0x8000)
1414 /* Return 1 if e-type number X is infinity, else return zero. */
1418 unsigned EMUSHORT x
[];
1425 if ((x
[NE
- 1] & 0x7fff) == 0x7fff)
1431 /* Check if e-type number is not a number. The bit pattern is one that we
1432 defined, so we know for sure how to detect it. */
1436 unsigned EMUSHORT x
[];
1441 /* NaN has maximum exponent */
1442 if ((x
[NE
- 1] & 0x7fff) != 0x7fff)
1444 /* ... and non-zero significand field. */
1445 for (i
= 0; i
< NE
- 1; i
++)
1455 /* Fill e-type number X with infinity pattern (IEEE)
1456 or largest possible number (non-IEEE). */
1460 register unsigned EMUSHORT
*x
;
1465 for (i
= 0; i
< NE
- 1; i
++)
1469 for (i
= 0; i
< NE
- 1; i
++)
1497 /* Output an e-type NaN.
1498 This generates Intel's quiet NaN pattern for extended real.
1499 The exponent is 7fff, the leading mantissa word is c000. */
1503 register unsigned EMUSHORT
*x
;
1508 for (i
= 0; i
< NE
- 2; i
++)
1511 *x
= (sign
<< 15) | 0x7fff;
1514 /* Move in an e-type number A, converting it to exploded e-type B. */
1518 unsigned EMUSHORT
*a
, *b
;
1520 register unsigned EMUSHORT
*p
, *q
;
1524 p
= a
+ (NE
- 1); /* point to last word of external number */
1525 /* get the sign bit */
1530 /* get the exponent */
1532 *q
++ &= 0x7fff; /* delete the sign bit */
1534 if ((*(q
- 1) & 0x7fff) == 0x7fff)
1540 for (i
= 3; i
< NI
; i
++)
1546 for (i
= 2; i
< NI
; i
++)
1552 /* clear high guard word */
1554 /* move in the significand */
1555 for (i
= 0; i
< NE
- 1; i
++)
1557 /* clear low guard word */
1561 /* Move out exploded e-type number A, converting it to e type B. */
1565 unsigned EMUSHORT
*a
, *b
;
1567 register unsigned EMUSHORT
*p
, *q
;
1568 unsigned EMUSHORT i
;
1572 q
= b
+ (NE
- 1); /* point to output exponent */
1573 /* combine sign and exponent */
1576 *q
-- = *p
++ | 0x8000;
1580 if (*(p
- 1) == 0x7fff)
1585 enan (b
, eiisneg (a
));
1593 /* skip over guard word */
1595 /* move the significand */
1596 for (j
= 0; j
< NE
- 1; j
++)
1600 /* Clear out exploded e-type number XI. */
1604 register unsigned EMUSHORT
*xi
;
1608 for (i
= 0; i
< NI
; i
++)
1612 /* Clear out exploded e-type XI, but don't touch the sign. */
1616 register unsigned EMUSHORT
*xi
;
1621 for (i
= 0; i
< NI
- 1; i
++)
1625 /* Move exploded e-type number from A to B. */
1629 register unsigned EMUSHORT
*a
, *b
;
1633 for (i
= 0; i
< NI
- 1; i
++)
1635 /* clear low guard word */
1639 /* Generate exploded e-type NaN.
1640 The explicit pattern for this is maximum exponent and
1641 top two significant bits set. */
1645 unsigned EMUSHORT x
[];
1653 /* Return nonzero if exploded e-type X is a NaN. */
1657 unsigned EMUSHORT x
[];
1661 if ((x
[E
] & 0x7fff) == 0x7fff)
1663 for (i
= M
+ 1; i
< NI
; i
++)
1672 /* Return nonzero if sign of exploded e-type X is nonzero. */
1676 unsigned EMUSHORT x
[];
1682 /* Fill exploded e-type X with infinity pattern.
1683 This has maximum exponent and significand all zeros. */
1687 unsigned EMUSHORT x
[];
1694 /* Return nonzero if exploded e-type X is infinite. */
1698 unsigned EMUSHORT x
[];
1705 if ((x
[E
] & 0x7fff) == 0x7fff)
1711 /* Compare significands of numbers in internal exploded e-type format.
1712 Guard words are included in the comparison.
1720 register unsigned EMUSHORT
*a
, *b
;
1724 a
+= M
; /* skip up to significand area */
1726 for (i
= M
; i
< NI
; i
++)
1734 if (*(--a
) > *(--b
))
1740 /* Shift significand of exploded e-type X down by 1 bit. */
1744 register unsigned EMUSHORT
*x
;
1746 register unsigned EMUSHORT bits
;
1749 x
+= M
; /* point to significand area */
1752 for (i
= M
; i
< NI
; i
++)
1764 /* Shift significand of exploded e-type X up by 1 bit. */
1768 register unsigned EMUSHORT
*x
;
1770 register unsigned EMUSHORT bits
;
1776 for (i
= M
; i
< NI
; i
++)
1789 /* Shift significand of exploded e-type X down by 8 bits. */
1793 register unsigned EMUSHORT
*x
;
1795 register unsigned EMUSHORT newbyt
, oldbyt
;
1800 for (i
= M
; i
< NI
; i
++)
1810 /* Shift significand of exploded e-type X up by 8 bits. */
1814 register unsigned EMUSHORT
*x
;
1817 register unsigned EMUSHORT newbyt
, oldbyt
;
1822 for (i
= M
; i
< NI
; i
++)
1832 /* Shift significand of exploded e-type X up by 16 bits. */
1836 register unsigned EMUSHORT
*x
;
1839 register unsigned EMUSHORT
*p
;
1844 for (i
= M
; i
< NI
- 1; i
++)
1850 /* Shift significand of exploded e-type X down by 16 bits. */
1854 register unsigned EMUSHORT
*x
;
1857 register unsigned EMUSHORT
*p
;
1862 for (i
= M
; i
< NI
- 1; i
++)
1868 /* Add significands of exploded e-type X and Y. X + Y replaces Y. */
1872 unsigned EMUSHORT
*x
, *y
;
1874 register unsigned EMULONG a
;
1881 for (i
= M
; i
< NI
; i
++)
1883 a
= (unsigned EMULONG
) (*x
) + (unsigned EMULONG
) (*y
) + carry
;
1888 *y
= (unsigned EMUSHORT
) a
;
1894 /* Subtract significands of exploded e-type X and Y. Y - X replaces Y. */
1898 unsigned EMUSHORT
*x
, *y
;
1907 for (i
= M
; i
< NI
; i
++)
1909 a
= (unsigned EMULONG
) (*y
) - (unsigned EMULONG
) (*x
) - carry
;
1914 *y
= (unsigned EMUSHORT
) a
;
1921 static unsigned EMUSHORT equot
[NI
];
1925 /* Radix 2 shift-and-add versions of multiply and divide */
1928 /* Divide significands */
1932 unsigned EMUSHORT den
[], num
[];
1935 register unsigned EMUSHORT
*p
, *q
;
1936 unsigned EMUSHORT j
;
1942 for (i
= M
; i
< NI
; i
++)
1947 /* Use faster compare and subtraction if denominator has only 15 bits of
1953 for (i
= M
+ 3; i
< NI
; i
++)
1958 if ((den
[M
+ 1] & 1) != 0)
1966 for (i
= 0; i
< NBITS
+ 2; i
++)
1984 /* The number of quotient bits to calculate is NBITS + 1 scaling guard
1985 bit + 1 roundoff bit. */
1990 for (i
= 0; i
< NBITS
+ 2; i
++)
1992 if (ecmpm (den
, num
) <= 0)
1995 j
= 1; /* quotient bit = 1 */
2009 /* test for nonzero remainder after roundoff bit */
2012 for (i
= M
; i
< NI
; i
++)
2020 for (i
= 0; i
< NI
; i
++)
2026 /* Multiply significands */
2029 unsigned EMUSHORT a
[], b
[];
2031 unsigned EMUSHORT
*p
, *q
;
2036 for (i
= M
; i
< NI
; i
++)
2041 while (*p
== 0) /* significand is not supposed to be zero */
2046 if ((*p
& 0xff) == 0)
2054 for (i
= 0; i
< k
; i
++)
2058 /* remember if there were any nonzero bits shifted out */
2065 for (i
= 0; i
< NI
; i
++)
2068 /* return flag for lost nonzero bits */
2074 /* Radix 65536 versions of multiply and divide. */
2076 /* Multiply significand of e-type number B
2077 by 16-bit quantity A, return e-type result to C. */
2082 unsigned EMUSHORT b
[], c
[];
2084 register unsigned EMUSHORT
*pp
;
2085 register unsigned EMULONG carry
;
2086 unsigned EMUSHORT
*ps
;
2087 unsigned EMUSHORT p
[NI
];
2088 unsigned EMULONG aa
, m
;
2097 for (i
=M
+1; i
<NI
; i
++)
2107 m
= (unsigned EMULONG
) aa
* *ps
--;
2108 carry
= (m
& 0xffff) + *pp
;
2109 *pp
-- = (unsigned EMUSHORT
)carry
;
2110 carry
= (carry
>> 16) + (m
>> 16) + *pp
;
2111 *pp
= (unsigned EMUSHORT
)carry
;
2112 *(pp
-1) = carry
>> 16;
2115 for (i
=M
; i
<NI
; i
++)
2119 /* Divide significands of exploded e-types NUM / DEN. Neither the
2120 numerator NUM nor the denominator DEN is permitted to have its high guard
2125 unsigned EMUSHORT den
[], num
[];
2128 register unsigned EMUSHORT
*p
;
2129 unsigned EMULONG tnum
;
2130 unsigned EMUSHORT j
, tdenm
, tquot
;
2131 unsigned EMUSHORT tprod
[NI
+1];
2137 for (i
=M
; i
<NI
; i
++)
2143 for (i
=M
; i
<NI
; i
++)
2145 /* Find trial quotient digit (the radix is 65536). */
2146 tnum
= (((unsigned EMULONG
) num
[M
]) << 16) + num
[M
+1];
2148 /* Do not execute the divide instruction if it will overflow. */
2149 if ((tdenm
* 0xffffL
) < tnum
)
2152 tquot
= tnum
/ tdenm
;
2153 /* Multiply denominator by trial quotient digit. */
2154 m16m ((unsigned int)tquot
, den
, tprod
);
2155 /* The quotient digit may have been overestimated. */
2156 if (ecmpm (tprod
, num
) > 0)
2160 if (ecmpm (tprod
, num
) > 0)
2170 /* test for nonzero remainder after roundoff bit */
2173 for (i
=M
; i
<NI
; i
++)
2180 for (i
=0; i
<NI
; i
++)
2186 /* Multiply significands of exploded e-type A and B, result in B. */
2190 unsigned EMUSHORT a
[], b
[];
2192 unsigned EMUSHORT
*p
, *q
;
2193 unsigned EMUSHORT pprod
[NI
];
2194 unsigned EMUSHORT j
;
2199 for (i
=M
; i
<NI
; i
++)
2205 for (i
=M
+1; i
<NI
; i
++)
2213 m16m ((unsigned int) *p
--, b
, pprod
);
2214 eaddm(pprod
, equot
);
2220 for (i
=0; i
<NI
; i
++)
2223 /* return flag for lost nonzero bits */
2229 /* Normalize and round off.
2231 The internal format number to be rounded is S.
2232 Input LOST is 0 if the value is exact. This is the so-called sticky bit.
2234 Input SUBFLG indicates whether the number was obtained
2235 by a subtraction operation. In that case if LOST is nonzero
2236 then the number is slightly smaller than indicated.
2238 Input EXP is the biased exponent, which may be negative.
2239 the exponent field of S is ignored but is replaced by
2240 EXP as adjusted by normalization and rounding.
2242 Input RCNTRL is the rounding control. If it is nonzero, the
2243 returned value will be rounded to RNDPRC bits.
2245 For future reference: In order for emdnorm to round off denormal
2246 significands at the right point, the input exponent must be
2247 adjusted to be the actual value it would have after conversion to
2248 the final floating point type. This adjustment has been
2249 implemented for all type conversions (etoe53, etc.) and decimal
2250 conversions, but not for the arithmetic functions (eadd, etc.).
2251 Data types having standard 15-bit exponents are not affected by
2252 this, but SFmode and DFmode are affected. For example, ediv with
2253 rndprc = 24 will not round correctly to 24-bit precision if the
2254 result is denormal. */
2256 static int rlast
= -1;
2258 static unsigned EMUSHORT rmsk
= 0;
2259 static unsigned EMUSHORT rmbit
= 0;
2260 static unsigned EMUSHORT rebit
= 0;
2262 static unsigned EMUSHORT rbit
[NI
];
2265 emdnorm (s
, lost
, subflg
, exp
, rcntrl
)
2266 unsigned EMUSHORT s
[];
2273 unsigned EMUSHORT r
;
2278 /* a blank significand could mean either zero or infinity. */
2291 if ((j
> NBITS
) && (exp
< 32767))
2299 if (exp
> (EMULONG
) (-NBITS
- 1))
2312 /* Round off, unless told not to by rcntrl. */
2315 /* Set up rounding parameters if the control register changed. */
2316 if (rndprc
!= rlast
)
2323 rw
= NI
- 1; /* low guard word */
2343 /* For DEC or IBM arithmetic */
2370 /* Shift down 1 temporarily if the data structure has an implied
2371 most significant bit and the number is denormal.
2372 Intel long double denormals also lose one bit of precision. */
2373 if ((exp
<= 0) && (rndprc
!= NBITS
)
2374 && ((rndprc
!= 64) || ((rndprc
== 64) && ! REAL_WORDS_BIG_ENDIAN
)))
2376 lost
|= s
[NI
- 1] & 1;
2379 /* Clear out all bits below the rounding bit,
2380 remembering in r if any were nonzero. */
2394 if ((r
& rmbit
) != 0)
2399 { /* round to even */
2400 if ((s
[re
] & rebit
) == 0)
2412 /* Undo the temporary shift for denormal values. */
2413 if ((exp
<= 0) && (rndprc
!= NBITS
)
2414 && ((rndprc
!= 64) || ((rndprc
== 64) && ! REAL_WORDS_BIG_ENDIAN
)))
2419 { /* overflow on roundoff */
2432 for (i
= 2; i
< NI
- 1; i
++)
2435 warning ("floating point overflow");
2439 for (i
= M
+ 1; i
< NI
- 1; i
++)
2442 if ((rndprc
< 64) || (rndprc
== 113))
2457 s
[1] = (unsigned EMUSHORT
) exp
;
2460 /* Subtract. C = B - A, all e type numbers. */
2462 static int subflg
= 0;
2466 unsigned EMUSHORT
*a
, *b
, *c
;
2480 /* Infinity minus infinity is a NaN.
2481 Test for subtracting infinities of the same sign. */
2482 if (eisinf (a
) && eisinf (b
)
2483 && ((eisneg (a
) ^ eisneg (b
)) == 0))
2485 mtherr ("esub", INVALID
);
2494 /* Add. C = A + B, all e type. */
2498 unsigned EMUSHORT
*a
, *b
, *c
;
2502 /* NaN plus anything is a NaN. */
2513 /* Infinity minus infinity is a NaN.
2514 Test for adding infinities of opposite signs. */
2515 if (eisinf (a
) && eisinf (b
)
2516 && ((eisneg (a
) ^ eisneg (b
)) != 0))
2518 mtherr ("esub", INVALID
);
2527 /* Arithmetic common to both addition and subtraction. */
2531 unsigned EMUSHORT
*a
, *b
, *c
;
2533 unsigned EMUSHORT ai
[NI
], bi
[NI
], ci
[NI
];
2535 EMULONG lt
, lta
, ltb
;
2556 /* compare exponents */
2561 { /* put the larger number in bi */
2571 if (lt
< (EMULONG
) (-NBITS
- 1))
2572 goto done
; /* answer same as larger addend */
2574 lost
= eshift (ai
, k
); /* shift the smaller number down */
2578 /* exponents were the same, so must compare significands */
2581 { /* the numbers are identical in magnitude */
2582 /* if different signs, result is zero */
2588 /* if same sign, result is double */
2589 /* double denormalized tiny number */
2590 if ((bi
[E
] == 0) && ((bi
[3] & 0x8000) == 0))
2595 /* add 1 to exponent unless both are zero! */
2596 for (j
= 1; j
< NI
- 1; j
++)
2600 /* This could overflow, but let emovo take care of that. */
2605 bi
[E
] = (unsigned EMUSHORT
) ltb
;
2609 { /* put the larger number in bi */
2625 emdnorm (bi
, lost
, subflg
, ltb
, 64);
2631 /* Divide: C = B/A, all e type. */
2635 unsigned EMUSHORT
*a
, *b
, *c
;
2637 unsigned EMUSHORT ai
[NI
], bi
[NI
];
2639 EMULONG lt
, lta
, ltb
;
2641 /* IEEE says if result is not a NaN, the sign is "-" if and only if
2642 operands have opposite signs -- but flush -0 to 0 later if not IEEE. */
2643 sign
= eisneg(a
) ^ eisneg(b
);
2646 /* Return any NaN input. */
2657 /* Zero over zero, or infinity over infinity, is a NaN. */
2658 if (((ecmp (a
, ezero
) == 0) && (ecmp (b
, ezero
) == 0))
2659 || (eisinf (a
) && eisinf (b
)))
2661 mtherr ("ediv", INVALID
);
2666 /* Infinity over anything else is infinity. */
2673 /* Anything else over infinity is zero. */
2685 { /* See if numerator is zero. */
2686 for (i
= 1; i
< NI
- 1; i
++)
2690 ltb
-= enormlz (bi
);
2700 { /* possible divide by zero */
2701 for (i
= 1; i
< NI
- 1; i
++)
2705 lta
-= enormlz (ai
);
2709 /* Divide by zero is not an invalid operation.
2710 It is a divide-by-zero operation! */
2712 mtherr ("ediv", SING
);
2718 /* calculate exponent */
2719 lt
= ltb
- lta
+ EXONE
;
2720 emdnorm (bi
, i
, 0, lt
, 64);
2727 && (ecmp (c
, ezero
) != 0)
2730 *(c
+(NE
-1)) |= 0x8000;
2732 *(c
+(NE
-1)) &= ~0x8000;
2735 /* Multiply e-types A and B, return e-type product C. */
2739 unsigned EMUSHORT
*a
, *b
, *c
;
2741 unsigned EMUSHORT ai
[NI
], bi
[NI
];
2743 EMULONG lt
, lta
, ltb
;
2745 /* IEEE says if result is not a NaN, the sign is "-" if and only if
2746 operands have opposite signs -- but flush -0 to 0 later if not IEEE. */
2747 sign
= eisneg(a
) ^ eisneg(b
);
2750 /* NaN times anything is the same NaN. */
2761 /* Zero times infinity is a NaN. */
2762 if ((eisinf (a
) && (ecmp (b
, ezero
) == 0))
2763 || (eisinf (b
) && (ecmp (a
, ezero
) == 0)))
2765 mtherr ("emul", INVALID
);
2770 /* Infinity times anything else is infinity. */
2772 if (eisinf (a
) || eisinf (b
))
2784 for (i
= 1; i
< NI
- 1; i
++)
2788 lta
-= enormlz (ai
);
2799 for (i
= 1; i
< NI
- 1; i
++)
2803 ltb
-= enormlz (bi
);
2812 /* Multiply significands */
2814 /* calculate exponent */
2815 lt
= lta
+ ltb
- (EXONE
- 1);
2816 emdnorm (bi
, j
, 0, lt
, 64);
2823 && (ecmp (c
, ezero
) != 0)
2826 *(c
+(NE
-1)) |= 0x8000;
2828 *(c
+(NE
-1)) &= ~0x8000;
2831 /* Convert double precision PE to e-type Y. */
2835 unsigned EMUSHORT
*pe
, *y
;
2844 ibmtoe (pe
, y
, DFmode
);
2847 register unsigned EMUSHORT r
;
2848 register unsigned EMUSHORT
*e
, *p
;
2849 unsigned EMUSHORT yy
[NI
];
2853 denorm
= 0; /* flag if denormalized number */
2855 if (! REAL_WORDS_BIG_ENDIAN
)
2861 yy
[M
] = (r
& 0x0f) | 0x10;
2862 r
&= ~0x800f; /* strip sign and 4 significand bits */
2867 if (! REAL_WORDS_BIG_ENDIAN
)
2869 if (((pe
[3] & 0xf) != 0) || (pe
[2] != 0)
2870 || (pe
[1] != 0) || (pe
[0] != 0))
2872 enan (y
, yy
[0] != 0);
2878 if (((pe
[0] & 0xf) != 0) || (pe
[1] != 0)
2879 || (pe
[2] != 0) || (pe
[3] != 0))
2881 enan (y
, yy
[0] != 0);
2892 #endif /* INFINITY */
2894 /* If zero exponent, then the significand is denormalized.
2895 So take back the understood high significand bit. */
2906 if (! REAL_WORDS_BIG_ENDIAN
)
2922 { /* if zero exponent, then normalize the significand */
2923 if ((k
= enormlz (yy
)) > NBITS
)
2926 yy
[E
] -= (unsigned EMUSHORT
) (k
- 1);
2929 #endif /* not IBM */
2930 #endif /* not DEC */
2933 /* Convert double extended precision float PE to e type Y. */
2937 unsigned EMUSHORT
*pe
, *y
;
2939 unsigned EMUSHORT yy
[NI
];
2940 unsigned EMUSHORT
*e
, *p
, *q
;
2945 for (i
= 0; i
< NE
- 5; i
++)
2947 /* This precision is not ordinarily supported on DEC or IBM. */
2949 for (i
= 0; i
< 5; i
++)
2953 p
= &yy
[0] + (NE
- 1);
2956 for (i
= 0; i
< 5; i
++)
2960 if (! REAL_WORDS_BIG_ENDIAN
)
2962 for (i
= 0; i
< 5; i
++)
2965 /* For denormal long double Intel format, shift significand up one
2966 -- but only if the top significand bit is zero. A top bit of 1
2967 is "pseudodenormal" when the exponent is zero. */
2968 if((yy
[NE
-1] & 0x7fff) == 0 && (yy
[NE
-2] & 0x8000) == 0)
2970 unsigned EMUSHORT temp
[NI
];
2980 p
= &yy
[0] + (NE
- 1);
2983 for (i
= 0; i
< 4; i
++)
2988 /* Point to the exponent field and check max exponent cases. */
2993 if (! REAL_WORDS_BIG_ENDIAN
)
2995 for (i
= 0; i
< 4; i
++)
2997 if ((i
!= 3 && pe
[i
] != 0)
2998 /* Anything but 0x8000 here, including 0, is a NaN. */
2999 || (i
== 3 && pe
[i
] != 0x8000))
3001 enan (y
, (*p
& 0x8000) != 0);
3008 for (i
= 1; i
<= 4; i
++)
3012 enan (y
, (*p
& 0x8000) != 0);
3024 #endif /* INFINITY */
3027 for (i
= 0; i
< NE
; i
++)
3031 /* Convert 128-bit long double precision float PE to e type Y. */
3035 unsigned EMUSHORT
*pe
, *y
;
3037 register unsigned EMUSHORT r
;
3038 unsigned EMUSHORT
*e
, *p
;
3039 unsigned EMUSHORT yy
[NI
];
3046 if (! REAL_WORDS_BIG_ENDIAN
)
3058 if (! REAL_WORDS_BIG_ENDIAN
)
3060 for (i
= 0; i
< 7; i
++)
3064 enan (y
, yy
[0] != 0);
3071 for (i
= 1; i
< 8; i
++)
3075 enan (y
, yy
[0] != 0);
3087 #endif /* INFINITY */
3091 if (! REAL_WORDS_BIG_ENDIAN
)
3093 for (i
= 0; i
< 7; i
++)
3099 for (i
= 0; i
< 7; i
++)
3103 /* If denormal, remove the implied bit; else shift down 1. */
3116 /* Convert single precision float PE to e type Y. */
3120 unsigned EMUSHORT
*pe
, *y
;
3124 ibmtoe (pe
, y
, SFmode
);
3127 register unsigned EMUSHORT r
;
3128 register unsigned EMUSHORT
*e
, *p
;
3129 unsigned EMUSHORT yy
[NI
];
3133 denorm
= 0; /* flag if denormalized number */
3136 if (! REAL_WORDS_BIG_ENDIAN
)
3146 yy
[M
] = (r
& 0x7f) | 0200;
3147 r
&= ~0x807f; /* strip sign and 7 significand bits */
3152 if (REAL_WORDS_BIG_ENDIAN
)
3154 if (((pe
[0] & 0x7f) != 0) || (pe
[1] != 0))
3156 enan (y
, yy
[0] != 0);
3162 if (((pe
[1] & 0x7f) != 0) || (pe
[0] != 0))
3164 enan (y
, yy
[0] != 0);
3175 #endif /* INFINITY */
3177 /* If zero exponent, then the significand is denormalized.
3178 So take back the understood high significand bit. */
3191 if (! REAL_WORDS_BIG_ENDIAN
)
3201 { /* if zero exponent, then normalize the significand */
3202 if ((k
= enormlz (yy
)) > NBITS
)
3205 yy
[E
] -= (unsigned EMUSHORT
) (k
- 1);
3208 #endif /* not IBM */
3211 /* Convert e-type X to IEEE 128-bit long double format E. */
3215 unsigned EMUSHORT
*x
, *e
;
3217 unsigned EMUSHORT xi
[NI
];
3224 make_nan (e
, eisneg (x
), TFmode
);
3229 exp
= (EMULONG
) xi
[E
];
3234 /* round off to nearest or even */
3237 emdnorm (xi
, 0, 0, exp
, 64);
3243 /* Convert exploded e-type X, that has already been rounded to
3244 113-bit precision, to IEEE 128-bit long double format Y. */
3248 unsigned EMUSHORT
*a
, *b
;
3250 register unsigned EMUSHORT
*p
, *q
;
3251 unsigned EMUSHORT i
;
3256 make_nan (b
, eiisneg (a
), TFmode
);
3261 if (REAL_WORDS_BIG_ENDIAN
)
3264 q
= b
+ 7; /* point to output exponent */
3266 /* If not denormal, delete the implied bit. */
3271 /* combine sign and exponent */
3273 if (REAL_WORDS_BIG_ENDIAN
)
3276 *q
++ = *p
++ | 0x8000;
3283 *q
-- = *p
++ | 0x8000;
3287 /* skip over guard word */
3289 /* move the significand */
3290 if (REAL_WORDS_BIG_ENDIAN
)
3292 for (i
= 0; i
< 7; i
++)
3297 for (i
= 0; i
< 7; i
++)
3302 /* Convert e-type X to IEEE double extended format E. */
3306 unsigned EMUSHORT
*x
, *e
;
3308 unsigned EMUSHORT xi
[NI
];
3315 make_nan (e
, eisneg (x
), XFmode
);
3320 /* adjust exponent for offset */
3321 exp
= (EMULONG
) xi
[E
];
3326 /* round off to nearest or even */
3329 emdnorm (xi
, 0, 0, exp
, 64);
3335 /* Convert exploded e-type X, that has already been rounded to
3336 64-bit precision, to IEEE double extended format Y. */
3340 unsigned EMUSHORT
*a
, *b
;
3342 register unsigned EMUSHORT
*p
, *q
;
3343 unsigned EMUSHORT i
;
3348 make_nan (b
, eiisneg (a
), XFmode
);
3352 /* Shift denormal long double Intel format significand down one bit. */
3353 if ((a
[E
] == 0) && ! REAL_WORDS_BIG_ENDIAN
)
3363 if (REAL_WORDS_BIG_ENDIAN
)
3367 q
= b
+ 4; /* point to output exponent */
3368 #if LONG_DOUBLE_TYPE_SIZE == 96
3369 /* Clear the last two bytes of 12-byte Intel format */
3375 /* combine sign and exponent */
3379 *q
++ = *p
++ | 0x8000;
3386 *q
-- = *p
++ | 0x8000;
3391 if (REAL_WORDS_BIG_ENDIAN
)
3394 *q
++ = *p
++ | 0x8000;
3402 *q
-- = *p
++ | 0x8000;
3407 /* skip over guard word */
3409 /* move the significand */
3411 for (i
= 0; i
< 4; i
++)
3415 for (i
= 0; i
< 4; i
++)
3419 if (REAL_WORDS_BIG_ENDIAN
)
3421 for (i
= 0; i
< 4; i
++)
3429 /* Intel long double infinity significand. */
3437 for (i
= 0; i
< 4; i
++)
3443 /* e type to double precision. */
3446 /* Convert e-type X to DEC-format double E. */
3450 unsigned EMUSHORT
*x
, *e
;
3452 etodec (x
, e
); /* see etodec.c */
3455 /* Convert exploded e-type X, that has already been rounded to
3456 56-bit double precision, to DEC double Y. */
3460 unsigned EMUSHORT
*x
, *y
;
3467 /* Convert e-type X to IBM 370-format double E. */
3471 unsigned EMUSHORT
*x
, *e
;
3473 etoibm (x
, e
, DFmode
);
3476 /* Convert exploded e-type X, that has already been rounded to
3477 56-bit precision, to IBM 370 double Y. */
3481 unsigned EMUSHORT
*x
, *y
;
3483 toibm (x
, y
, DFmode
);
3486 #else /* it's neither DEC nor IBM */
3488 /* Convert e-type X to IEEE double E. */
3492 unsigned EMUSHORT
*x
, *e
;
3494 unsigned EMUSHORT xi
[NI
];
3501 make_nan (e
, eisneg (x
), DFmode
);
3506 /* adjust exponent for offsets */
3507 exp
= (EMULONG
) xi
[E
] - (EXONE
- 0x3ff);
3512 /* round off to nearest or even */
3515 emdnorm (xi
, 0, 0, exp
, 64);
3521 /* Convert exploded e-type X, that has already been rounded to
3522 53-bit precision, to IEEE double Y. */
3526 unsigned EMUSHORT
*x
, *y
;
3528 unsigned EMUSHORT i
;
3529 unsigned EMUSHORT
*p
;
3534 make_nan (y
, eiisneg (x
), DFmode
);
3540 if (! REAL_WORDS_BIG_ENDIAN
)
3543 *y
= 0; /* output high order */
3545 *y
= 0x8000; /* output sign bit */
3548 if (i
>= (unsigned int) 2047)
3549 { /* Saturate at largest number less than infinity. */
3552 if (! REAL_WORDS_BIG_ENDIAN
)
3566 *y
|= (unsigned EMUSHORT
) 0x7fef;
3567 if (! REAL_WORDS_BIG_ENDIAN
)
3592 i
|= *p
++ & (unsigned EMUSHORT
) 0x0f; /* *p = xi[M] */
3593 *y
|= (unsigned EMUSHORT
) i
; /* high order output already has sign bit set */
3594 if (! REAL_WORDS_BIG_ENDIAN
)
3609 #endif /* not IBM */
3610 #endif /* not DEC */
3614 /* e type to single precision. */
3617 /* Convert e-type X to IBM 370 float E. */
3621 unsigned EMUSHORT
*x
, *e
;
3623 etoibm (x
, e
, SFmode
);
3626 /* Convert exploded e-type X, that has already been rounded to
3627 float precision, to IBM 370 float Y. */
3631 unsigned EMUSHORT
*x
, *y
;
3633 toibm (x
, y
, SFmode
);
3637 /* Convert e-type X to IEEE float E. DEC float is the same as IEEE float. */
3641 unsigned EMUSHORT
*x
, *e
;
3644 unsigned EMUSHORT xi
[NI
];
3650 make_nan (e
, eisneg (x
), SFmode
);
3655 /* adjust exponent for offsets */
3656 exp
= (EMULONG
) xi
[E
] - (EXONE
- 0177);
3661 /* round off to nearest or even */
3664 emdnorm (xi
, 0, 0, exp
, 64);
3670 /* Convert exploded e-type X, that has already been rounded to
3671 float precision, to IEEE float Y. */
3675 unsigned EMUSHORT
*x
, *y
;
3677 unsigned EMUSHORT i
;
3678 unsigned EMUSHORT
*p
;
3683 make_nan (y
, eiisneg (x
), SFmode
);
3689 if (! REAL_WORDS_BIG_ENDIAN
)
3695 *y
= 0; /* output high order */
3697 *y
= 0x8000; /* output sign bit */
3700 /* Handle overflow cases. */
3704 *y
|= (unsigned EMUSHORT
) 0x7f80;
3709 if (! REAL_WORDS_BIG_ENDIAN
)
3717 #else /* no INFINITY */
3718 *y
|= (unsigned EMUSHORT
) 0x7f7f;
3723 if (! REAL_WORDS_BIG_ENDIAN
)
3734 #endif /* no INFINITY */
3746 i
|= *p
++ & (unsigned EMUSHORT
) 0x7f; /* *p = xi[M] */
3747 /* High order output already has sign bit set. */
3753 if (! REAL_WORDS_BIG_ENDIAN
)
3762 #endif /* not IBM */
3764 /* Compare two e type numbers.
3768 -2 if either a or b is a NaN. */
3772 unsigned EMUSHORT
*a
, *b
;
3774 unsigned EMUSHORT ai
[NI
], bi
[NI
];
3775 register unsigned EMUSHORT
*p
, *q
;
3780 if (eisnan (a
) || eisnan (b
))
3789 { /* the signs are different */
3791 for (i
= 1; i
< NI
- 1; i
++)
3805 /* both are the same sign */
3820 return (0); /* equality */
3824 if (*(--p
) > *(--q
))
3825 return (msign
); /* p is bigger */
3827 return (-msign
); /* p is littler */
3830 /* Find e-type nearest integer to X, as floor (X + 0.5). */
3834 unsigned EMUSHORT
*x
, *y
;
3840 /* Convert HOST_WIDE_INT LP to e type Y. */
3845 unsigned EMUSHORT
*y
;
3847 unsigned EMUSHORT yi
[NI
];
3848 unsigned HOST_WIDE_INT ll
;
3854 /* make it positive */
3855 ll
= (unsigned HOST_WIDE_INT
) (-(*lp
));
3856 yi
[0] = 0xffff; /* put correct sign in the e type number */
3860 ll
= (unsigned HOST_WIDE_INT
) (*lp
);
3862 /* move the long integer to yi significand area */
3863 #if HOST_BITS_PER_WIDE_INT == 64
3864 yi
[M
] = (unsigned EMUSHORT
) (ll
>> 48);
3865 yi
[M
+ 1] = (unsigned EMUSHORT
) (ll
>> 32);
3866 yi
[M
+ 2] = (unsigned EMUSHORT
) (ll
>> 16);
3867 yi
[M
+ 3] = (unsigned EMUSHORT
) ll
;
3868 yi
[E
] = EXONE
+ 47; /* exponent if normalize shift count were 0 */
3870 yi
[M
] = (unsigned EMUSHORT
) (ll
>> 16);
3871 yi
[M
+ 1] = (unsigned EMUSHORT
) ll
;
3872 yi
[E
] = EXONE
+ 15; /* exponent if normalize shift count were 0 */
3875 if ((k
= enormlz (yi
)) > NBITS
)/* normalize the significand */
3876 ecleaz (yi
); /* it was zero */
3878 yi
[E
] -= (unsigned EMUSHORT
) k
;/* subtract shift count from exponent */
3879 emovo (yi
, y
); /* output the answer */
3882 /* Convert unsigned HOST_WIDE_INT LP to e type Y. */
3886 unsigned HOST_WIDE_INT
*lp
;
3887 unsigned EMUSHORT
*y
;
3889 unsigned EMUSHORT yi
[NI
];
3890 unsigned HOST_WIDE_INT ll
;
3896 /* move the long integer to ayi significand area */
3897 #if HOST_BITS_PER_WIDE_INT == 64
3898 yi
[M
] = (unsigned EMUSHORT
) (ll
>> 48);
3899 yi
[M
+ 1] = (unsigned EMUSHORT
) (ll
>> 32);
3900 yi
[M
+ 2] = (unsigned EMUSHORT
) (ll
>> 16);
3901 yi
[M
+ 3] = (unsigned EMUSHORT
) ll
;
3902 yi
[E
] = EXONE
+ 47; /* exponent if normalize shift count were 0 */
3904 yi
[M
] = (unsigned EMUSHORT
) (ll
>> 16);
3905 yi
[M
+ 1] = (unsigned EMUSHORT
) ll
;
3906 yi
[E
] = EXONE
+ 15; /* exponent if normalize shift count were 0 */
3909 if ((k
= enormlz (yi
)) > NBITS
)/* normalize the significand */
3910 ecleaz (yi
); /* it was zero */
3912 yi
[E
] -= (unsigned EMUSHORT
) k
; /* subtract shift count from exponent */
3913 emovo (yi
, y
); /* output the answer */
3917 /* Find signed HOST_WIDE_INT integer I and floating point fractional
3918 part FRAC of e-type (packed internal format) floating point input X.
3919 The integer output I has the sign of the input, except that
3920 positive overflow is permitted if FIXUNS_TRUNC_LIKE_FIX_TRUNC.
3921 The output e-type fraction FRAC is the positive fractional
3926 unsigned EMUSHORT
*x
;
3928 unsigned EMUSHORT
*frac
;
3930 unsigned EMUSHORT xi
[NI
];
3932 unsigned HOST_WIDE_INT ll
;
3935 k
= (int) xi
[E
] - (EXONE
- 1);
3938 /* if exponent <= 0, integer = 0 and real output is fraction */
3943 if (k
> (HOST_BITS_PER_WIDE_INT
- 1))
3945 /* long integer overflow: output large integer
3946 and correct fraction */
3948 *i
= ((unsigned HOST_WIDE_INT
) 1) << (HOST_BITS_PER_WIDE_INT
- 1);
3951 #ifdef FIXUNS_TRUNC_LIKE_FIX_TRUNC
3952 /* In this case, let it overflow and convert as if unsigned. */
3953 euifrac (x
, &ll
, frac
);
3954 *i
= (HOST_WIDE_INT
) ll
;
3957 /* In other cases, return the largest positive integer. */
3958 *i
= (((unsigned HOST_WIDE_INT
) 1) << (HOST_BITS_PER_WIDE_INT
- 1)) - 1;
3963 warning ("overflow on truncation to integer");
3967 /* Shift more than 16 bits: first shift up k-16 mod 16,
3968 then shift up by 16's. */
3969 j
= k
- ((k
>> 4) << 4);
3976 ll
= (ll
<< 16) | xi
[M
];
3978 while ((k
-= 16) > 0);
3985 /* shift not more than 16 bits */
3987 *i
= (HOST_WIDE_INT
) xi
[M
] & 0xffff;
3994 if ((k
= enormlz (xi
)) > NBITS
)
3997 xi
[E
] -= (unsigned EMUSHORT
) k
;
4003 /* Find unsigned HOST_WIDE_INT integer I and floating point fractional part
4004 FRAC of e-type X. A negative input yields integer output = 0 but
4005 correct fraction. */
4008 euifrac (x
, i
, frac
)
4009 unsigned EMUSHORT
*x
;
4010 unsigned HOST_WIDE_INT
*i
;
4011 unsigned EMUSHORT
*frac
;
4013 unsigned HOST_WIDE_INT ll
;
4014 unsigned EMUSHORT xi
[NI
];
4018 k
= (int) xi
[E
] - (EXONE
- 1);
4021 /* if exponent <= 0, integer = 0 and argument is fraction */
4026 if (k
> HOST_BITS_PER_WIDE_INT
)
4028 /* Long integer overflow: output large integer
4029 and correct fraction.
4030 Note, the BSD microvax compiler says that ~(0UL)
4031 is a syntax error. */
4035 warning ("overflow on truncation to unsigned integer");
4039 /* Shift more than 16 bits: first shift up k-16 mod 16,
4040 then shift up by 16's. */
4041 j
= k
- ((k
>> 4) << 4);
4048 ll
= (ll
<< 16) | xi
[M
];
4050 while ((k
-= 16) > 0);
4055 /* shift not more than 16 bits */
4057 *i
= (HOST_WIDE_INT
) xi
[M
] & 0xffff;
4060 if (xi
[0]) /* A negative value yields unsigned integer 0. */
4066 if ((k
= enormlz (xi
)) > NBITS
)
4069 xi
[E
] -= (unsigned EMUSHORT
) k
;
4074 /* Shift the significand of exploded e-type X up or down by SC bits. */
4078 unsigned EMUSHORT
*x
;
4081 unsigned EMUSHORT lost
;
4082 unsigned EMUSHORT
*p
;
4095 lost
|= *p
; /* remember lost bits */
4136 return ((int) lost
);
4139 /* Shift normalize the significand area of exploded e-type X.
4140 Return the shift count (up = positive). */
4144 unsigned EMUSHORT x
[];
4146 register unsigned EMUSHORT
*p
;
4155 return (0); /* already normalized */
4161 /* With guard word, there are NBITS+16 bits available.
4162 Return true if all are zero. */
4166 /* see if high byte is zero */
4167 while ((*p
& 0xff00) == 0)
4172 /* now shift 1 bit at a time */
4173 while ((*p
& 0x8000) == 0)
4179 mtherr ("enormlz", UNDERFLOW
);
4185 /* Normalize by shifting down out of the high guard word
4186 of the significand */
4201 mtherr ("enormlz", OVERFLOW
);
4208 /* Powers of ten used in decimal <-> binary conversions. */
4213 #if LONG_DOUBLE_TYPE_SIZE == 128
4214 static unsigned EMUSHORT etens
[NTEN
+ 1][NE
] =
4216 {0x6576, 0x4a92, 0x804a, 0x153f,
4217 0xc94c, 0x979a, 0x8a20, 0x5202, 0xc460, 0x7525,}, /* 10**4096 */
4218 {0x6a32, 0xce52, 0x329a, 0x28ce,
4219 0xa74d, 0x5de4, 0xc53d, 0x3b5d, 0x9e8b, 0x5a92,}, /* 10**2048 */
4220 {0x526c, 0x50ce, 0xf18b, 0x3d28,
4221 0x650d, 0x0c17, 0x8175, 0x7586, 0xc976, 0x4d48,},
4222 {0x9c66, 0x58f8, 0xbc50, 0x5c54,
4223 0xcc65, 0x91c6, 0xa60e, 0xa0ae, 0xe319, 0x46a3,},
4224 {0x851e, 0xeab7, 0x98fe, 0x901b,
4225 0xddbb, 0xde8d, 0x9df9, 0xebfb, 0xaa7e, 0x4351,},
4226 {0x0235, 0x0137, 0x36b1, 0x336c,
4227 0xc66f, 0x8cdf, 0x80e9, 0x47c9, 0x93ba, 0x41a8,},
4228 {0x50f8, 0x25fb, 0xc76b, 0x6b71,
4229 0x3cbf, 0xa6d5, 0xffcf, 0x1f49, 0xc278, 0x40d3,},
4230 {0x0000, 0x0000, 0x0000, 0x0000,
4231 0xf020, 0xb59d, 0x2b70, 0xada8, 0x9dc5, 0x4069,},
4232 {0x0000, 0x0000, 0x0000, 0x0000,
4233 0x0000, 0x0000, 0x0400, 0xc9bf, 0x8e1b, 0x4034,},
4234 {0x0000, 0x0000, 0x0000, 0x0000,
4235 0x0000, 0x0000, 0x0000, 0x2000, 0xbebc, 0x4019,},
4236 {0x0000, 0x0000, 0x0000, 0x0000,
4237 0x0000, 0x0000, 0x0000, 0x0000, 0x9c40, 0x400c,},
4238 {0x0000, 0x0000, 0x0000, 0x0000,
4239 0x0000, 0x0000, 0x0000, 0x0000, 0xc800, 0x4005,},
4240 {0x0000, 0x0000, 0x0000, 0x0000,
4241 0x0000, 0x0000, 0x0000, 0x0000, 0xa000, 0x4002,}, /* 10**1 */
4244 static unsigned EMUSHORT emtens
[NTEN
+ 1][NE
] =
4246 {0x2030, 0xcffc, 0xa1c3, 0x8123,
4247 0x2de3, 0x9fde, 0xd2ce, 0x04c8, 0xa6dd, 0x0ad8,}, /* 10**-4096 */
4248 {0x8264, 0xd2cb, 0xf2ea, 0x12d4,
4249 0x4925, 0x2de4, 0x3436, 0x534f, 0xceae, 0x256b,}, /* 10**-2048 */
4250 {0xf53f, 0xf698, 0x6bd3, 0x0158,
4251 0x87a6, 0xc0bd, 0xda57, 0x82a5, 0xa2a6, 0x32b5,},
4252 {0xe731, 0x04d4, 0xe3f2, 0xd332,
4253 0x7132, 0xd21c, 0xdb23, 0xee32, 0x9049, 0x395a,},
4254 {0xa23e, 0x5308, 0xfefb, 0x1155,
4255 0xfa91, 0x1939, 0x637a, 0x4325, 0xc031, 0x3cac,},
4256 {0xe26d, 0xdbde, 0xd05d, 0xb3f6,
4257 0xac7c, 0xe4a0, 0x64bc, 0x467c, 0xddd0, 0x3e55,},
4258 {0x2a20, 0x6224, 0x47b3, 0x98d7,
4259 0x3f23, 0xe9a5, 0xa539, 0xea27, 0xa87f, 0x3f2a,},
4260 {0x0b5b, 0x4af2, 0xa581, 0x18ed,
4261 0x67de, 0x94ba, 0x4539, 0x1ead, 0xcfb1, 0x3f94,},
4262 {0xbf71, 0xa9b3, 0x7989, 0xbe68,
4263 0x4c2e, 0xe15b, 0xc44d, 0x94be, 0xe695, 0x3fc9,},
4264 {0x3d4d, 0x7c3d, 0x36ba, 0x0d2b,
4265 0xfdc2, 0xcefc, 0x8461, 0x7711, 0xabcc, 0x3fe4,},
4266 {0xc155, 0xa4a8, 0x404e, 0x6113,
4267 0xd3c3, 0x652b, 0xe219, 0x1758, 0xd1b7, 0x3ff1,},
4268 {0xd70a, 0x70a3, 0x0a3d, 0xa3d7,
4269 0x3d70, 0xd70a, 0x70a3, 0x0a3d, 0xa3d7, 0x3ff8,},
4270 {0xcccd, 0xcccc, 0xcccc, 0xcccc,
4271 0xcccc, 0xcccc, 0xcccc, 0xcccc, 0xcccc, 0x3ffb,}, /* 10**-1 */
4274 /* LONG_DOUBLE_TYPE_SIZE is other than 128 */
4275 static unsigned EMUSHORT etens
[NTEN
+ 1][NE
] =
4277 {0xc94c, 0x979a, 0x8a20, 0x5202, 0xc460, 0x7525,}, /* 10**4096 */
4278 {0xa74d, 0x5de4, 0xc53d, 0x3b5d, 0x9e8b, 0x5a92,}, /* 10**2048 */
4279 {0x650d, 0x0c17, 0x8175, 0x7586, 0xc976, 0x4d48,},
4280 {0xcc65, 0x91c6, 0xa60e, 0xa0ae, 0xe319, 0x46a3,},
4281 {0xddbc, 0xde8d, 0x9df9, 0xebfb, 0xaa7e, 0x4351,},
4282 {0xc66f, 0x8cdf, 0x80e9, 0x47c9, 0x93ba, 0x41a8,},
4283 {0x3cbf, 0xa6d5, 0xffcf, 0x1f49, 0xc278, 0x40d3,},
4284 {0xf020, 0xb59d, 0x2b70, 0xada8, 0x9dc5, 0x4069,},
4285 {0x0000, 0x0000, 0x0400, 0xc9bf, 0x8e1b, 0x4034,},
4286 {0x0000, 0x0000, 0x0000, 0x2000, 0xbebc, 0x4019,},
4287 {0x0000, 0x0000, 0x0000, 0x0000, 0x9c40, 0x400c,},
4288 {0x0000, 0x0000, 0x0000, 0x0000, 0xc800, 0x4005,},
4289 {0x0000, 0x0000, 0x0000, 0x0000, 0xa000, 0x4002,}, /* 10**1 */
4292 static unsigned EMUSHORT emtens
[NTEN
+ 1][NE
] =
4294 {0x2de4, 0x9fde, 0xd2ce, 0x04c8, 0xa6dd, 0x0ad8,}, /* 10**-4096 */
4295 {0x4925, 0x2de4, 0x3436, 0x534f, 0xceae, 0x256b,}, /* 10**-2048 */
4296 {0x87a6, 0xc0bd, 0xda57, 0x82a5, 0xa2a6, 0x32b5,},
4297 {0x7133, 0xd21c, 0xdb23, 0xee32, 0x9049, 0x395a,},
4298 {0xfa91, 0x1939, 0x637a, 0x4325, 0xc031, 0x3cac,},
4299 {0xac7d, 0xe4a0, 0x64bc, 0x467c, 0xddd0, 0x3e55,},
4300 {0x3f24, 0xe9a5, 0xa539, 0xea27, 0xa87f, 0x3f2a,},
4301 {0x67de, 0x94ba, 0x4539, 0x1ead, 0xcfb1, 0x3f94,},
4302 {0x4c2f, 0xe15b, 0xc44d, 0x94be, 0xe695, 0x3fc9,},
4303 {0xfdc2, 0xcefc, 0x8461, 0x7711, 0xabcc, 0x3fe4,},
4304 {0xd3c3, 0x652b, 0xe219, 0x1758, 0xd1b7, 0x3ff1,},
4305 {0x3d71, 0xd70a, 0x70a3, 0x0a3d, 0xa3d7, 0x3ff8,},
4306 {0xcccd, 0xcccc, 0xcccc, 0xcccc, 0xcccc, 0x3ffb,}, /* 10**-1 */
4310 /* Convert float value X to ASCII string STRING with NDIG digits after
4311 the decimal point. */
4314 e24toasc (x
, string
, ndigs
)
4315 unsigned EMUSHORT x
[];
4319 unsigned EMUSHORT w
[NI
];
4322 etoasc (w
, string
, ndigs
);
4325 /* Convert double value X to ASCII string STRING with NDIG digits after
4326 the decimal point. */
4329 e53toasc (x
, string
, ndigs
)
4330 unsigned EMUSHORT x
[];
4334 unsigned EMUSHORT w
[NI
];
4337 etoasc (w
, string
, ndigs
);
4340 /* Convert double extended value X to ASCII string STRING with NDIG digits
4341 after the decimal point. */
4344 e64toasc (x
, string
, ndigs
)
4345 unsigned EMUSHORT x
[];
4349 unsigned EMUSHORT w
[NI
];
4352 etoasc (w
, string
, ndigs
);
4355 /* Convert 128-bit long double value X to ASCII string STRING with NDIG digits
4356 after the decimal point. */
4359 e113toasc (x
, string
, ndigs
)
4360 unsigned EMUSHORT x
[];
4364 unsigned EMUSHORT w
[NI
];
4367 etoasc (w
, string
, ndigs
);
4370 /* Convert e-type X to ASCII string STRING with NDIGS digits after
4371 the decimal point. */
4373 static char wstring
[80]; /* working storage for ASCII output */
4376 etoasc (x
, string
, ndigs
)
4377 unsigned EMUSHORT x
[];
4382 unsigned EMUSHORT y
[NI
], t
[NI
], u
[NI
], w
[NI
];
4383 unsigned EMUSHORT
*p
, *r
, *ten
;
4384 unsigned EMUSHORT sign
;
4385 int i
, j
, k
, expon
, rndsav
;
4387 unsigned EMUSHORT m
;
4398 sprintf (wstring
, " NaN ");
4402 rndprc
= NBITS
; /* set to full precision */
4403 emov (x
, y
); /* retain external format */
4404 if (y
[NE
- 1] & 0x8000)
4407 y
[NE
- 1] &= 0x7fff;
4414 ten
= &etens
[NTEN
][0];
4416 /* Test for zero exponent */
4419 for (k
= 0; k
< NE
- 1; k
++)
4422 goto tnzro
; /* denormalized number */
4424 goto isone
; /* valid all zeros */
4428 /* Test for infinity. */
4429 if (y
[NE
- 1] == 0x7fff)
4432 sprintf (wstring
, " -Infinity ");
4434 sprintf (wstring
, " Infinity ");
4438 /* Test for exponent nonzero but significand denormalized.
4439 * This is an error condition.
4441 if ((y
[NE
- 1] != 0) && ((y
[NE
- 2] & 0x8000) == 0))
4443 mtherr ("etoasc", DOMAIN
);
4444 sprintf (wstring
, "NaN");
4448 /* Compare to 1.0 */
4457 { /* Number is greater than 1 */
4458 /* Convert significand to an integer and strip trailing decimal zeros. */
4460 u
[NE
- 1] = EXONE
+ NBITS
- 1;
4462 p
= &etens
[NTEN
- 4][0];
4468 for (j
= 0; j
< NE
- 1; j
++)
4481 /* Rescale from integer significand */
4482 u
[NE
- 1] += y
[NE
- 1] - (unsigned int) (EXONE
+ NBITS
- 1);
4484 /* Find power of 10 */
4488 /* An unordered compare result shouldn't happen here. */
4489 while (ecmp (ten
, u
) <= 0)
4491 if (ecmp (p
, u
) <= 0)
4504 { /* Number is less than 1.0 */
4505 /* Pad significand with trailing decimal zeros. */
4508 while ((y
[NE
- 2] & 0x8000) == 0)
4517 for (i
= 0; i
< NDEC
+ 1; i
++)
4519 if ((w
[NI
- 1] & 0x7) != 0)
4521 /* multiply by 10 */
4534 if (eone
[NE
- 1] <= u
[1])
4546 while (ecmp (eone
, w
) > 0)
4548 if (ecmp (p
, w
) >= 0)
4563 /* Find the first (leading) digit. */
4569 digit
= equot
[NI
- 1];
4570 while ((digit
== 0) && (ecmp (y
, ezero
) != 0))
4578 digit
= equot
[NI
- 1];
4586 /* Examine number of digits requested by caller. */
4604 *s
++ = (char)digit
+ '0';
4607 /* Generate digits after the decimal point. */
4608 for (k
= 0; k
<= ndigs
; k
++)
4610 /* multiply current number by 10, without normalizing */
4617 *s
++ = (char) equot
[NI
- 1] + '0';
4619 digit
= equot
[NI
- 1];
4622 /* round off the ASCII string */
4625 /* Test for critical rounding case in ASCII output. */
4629 if (ecmp (t
, ezero
) != 0)
4630 goto roun
; /* round to nearest */
4631 if ((*(s
- 1) & 1) == 0)
4632 goto doexp
; /* round to even */
4634 /* Round up and propagate carry-outs */
4638 /* Carry out to most significant digit? */
4645 /* Most significant digit carries to 10? */
4653 /* Round up and carry out from less significant digits */
4665 sprintf (ss, "e+%d", expon);
4667 sprintf (ss, "e%d", expon);
4669 sprintf (ss
, "e%d", expon
);
4672 /* copy out the working string */
4675 while (*ss
== ' ') /* strip possible leading space */
4677 while ((*s
++ = *ss
++) != '\0')
4682 /* Convert ASCII string to floating point.
4684 Numeric input is a free format decimal number of any length, with
4685 or without decimal point. Entering E after the number followed by an
4686 integer number causes the second number to be interpreted as a power of
4687 10 to be multiplied by the first number (i.e., "scientific" notation). */
4689 /* Convert ASCII string S to single precision float value Y. */
4694 unsigned EMUSHORT
*y
;
4700 /* Convert ASCII string S to double precision value Y. */
4705 unsigned EMUSHORT
*y
;
4707 #if defined(DEC) || defined(IBM)
4715 /* Convert ASCII string S to double extended value Y. */
4720 unsigned EMUSHORT
*y
;
4725 /* Convert ASCII string S to 128-bit long double Y. */
4730 unsigned EMUSHORT
*y
;
4732 asctoeg (s
, y
, 113);
4735 /* Convert ASCII string S to e type Y. */
4740 unsigned EMUSHORT
*y
;
4742 asctoeg (s
, y
, NBITS
);
4745 /* Convert ASCII string SS to e type Y, with a specified rounding precision
4749 asctoeg (ss
, y
, oprec
)
4751 unsigned EMUSHORT
*y
;
4754 unsigned EMUSHORT yy
[NI
], xt
[NI
], tt
[NI
];
4755 int esign
, decflg
, sgnflg
, nexp
, exp
, prec
, lost
;
4756 int k
, trail
, c
, rndsav
;
4758 unsigned EMUSHORT nsign
, *p
;
4759 char *sp
, *s
, *lstr
;
4761 /* Copy the input string. */
4762 lstr
= (char *) alloca (strlen (ss
) + 1);
4764 while (*s
== ' ') /* skip leading spaces */
4767 while ((*sp
++ = *s
++) != '\0')
4772 rndprc
= NBITS
; /* Set to full precision */
4785 if ((k
>= 0) && (k
<= 9))
4787 /* Ignore leading zeros */
4788 if ((prec
== 0) && (decflg
== 0) && (k
== 0))
4790 /* Identify and strip trailing zeros after the decimal point. */
4791 if ((trail
== 0) && (decflg
!= 0))
4794 while ((*sp
>= '0') && (*sp
<= '9'))
4796 /* Check for syntax error */
4798 if ((c
!= 'e') && (c
!= 'E') && (c
!= '\0')
4799 && (c
!= '\n') && (c
!= '\r') && (c
!= ' ')
4810 /* If enough digits were given to more than fill up the yy register,
4811 continuing until overflow into the high guard word yy[2]
4812 guarantees that there will be a roundoff bit at the top
4813 of the low guard word after normalization. */
4818 nexp
+= 1; /* count digits after decimal point */
4819 eshup1 (yy
); /* multiply current number by 10 */
4825 xt
[NI
- 2] = (unsigned EMUSHORT
) k
;
4830 /* Mark any lost non-zero digit. */
4832 /* Count lost digits before the decimal point. */
4847 case '.': /* decimal point */
4877 mtherr ("asctoe", DOMAIN
);
4886 /* Exponent interpretation */
4892 /* check for + or - */
4900 while ((*s
>= '0') && (*s
<= '9'))
4904 if (exp
> -(MINDECEXP
))
4914 if (exp
> MAXDECEXP
)
4918 yy
[E
] = 0x7fff; /* infinity */
4921 if (exp
< MINDECEXP
)
4930 /* Pad trailing zeros to minimize power of 10, per IEEE spec. */
4931 while ((nexp
> 0) && (yy
[2] == 0))
4943 if ((k
= enormlz (yy
)) > NBITS
)
4948 lexp
= (EXONE
- 1 + NBITS
) - k
;
4949 emdnorm (yy
, lost
, 0, lexp
, 64);
4951 /* Convert to external format:
4953 Multiply by 10**nexp. If precision is 64 bits,
4954 the maximum relative error incurred in forming 10**n
4955 for 0 <= n <= 324 is 8.2e-20, at 10**180.
4956 For 0 <= n <= 999, the peak relative error is 1.4e-19 at 10**947.
4957 For 0 >= n >= -999, it is -1.55e-19 at 10**-435. */
4972 /* Punt. Can't handle this without 2 divides. */
4973 emovi (etens
[0], tt
);
4980 p
= &etens
[NTEN
][0];
4990 while (exp
<= MAXP
);
5008 /* Round and convert directly to the destination type */
5010 lexp
-= EXONE
- 0x3ff;
5012 else if (oprec
== 24 || oprec
== 56)
5013 lexp
-= EXONE
- (0x41 << 2);
5015 else if (oprec
== 24)
5016 lexp
-= EXONE
- 0177;
5019 else if (oprec
== 56)
5020 lexp
-= EXONE
- 0201;
5023 emdnorm (yy
, k
, 0, lexp
, 64);
5033 todec (yy
, y
); /* see etodec.c */
5038 toibm (yy
, y
, DFmode
);
5061 /* Return Y = largest integer not greater than X (truncated toward minus
5064 static unsigned EMUSHORT bmask
[] =
5087 unsigned EMUSHORT x
[], y
[];
5089 register unsigned EMUSHORT
*p
;
5091 unsigned EMUSHORT f
[NE
];
5093 emov (x
, f
); /* leave in external format */
5094 expon
= (int) f
[NE
- 1];
5095 e
= (expon
& 0x7fff) - (EXONE
- 1);
5101 /* number of bits to clear out */
5113 /* clear the remaining bits */
5115 /* truncate negatives toward minus infinity */
5118 if ((unsigned EMUSHORT
) expon
& (unsigned EMUSHORT
) 0x8000)
5120 for (i
= 0; i
< NE
- 1; i
++)
5132 /* Return S and EXP such that S * 2^EXP = X and .5 <= S < 1.
5133 For example, 1.1 = 0.55 * 2^1. */
5137 unsigned EMUSHORT x
[];
5139 unsigned EMUSHORT s
[];
5141 unsigned EMUSHORT xi
[NI
];
5145 /* Handle denormalized numbers properly using long integer exponent. */
5146 li
= (EMULONG
) ((EMUSHORT
) xi
[1]);
5154 *exp
= (int) (li
- 0x3ffe);
5157 /* Return e type Y = X * 2^PWR2. */
5161 unsigned EMUSHORT x
[];
5163 unsigned EMUSHORT y
[];
5165 unsigned EMUSHORT xi
[NI
];
5173 emdnorm (xi
, i
, i
, li
, 64);
5178 /* C = remainder after dividing B by A, all e type values.
5179 Least significant integer quotient bits left in EQUOT. */
5183 unsigned EMUSHORT a
[], b
[], c
[];
5185 unsigned EMUSHORT den
[NI
], num
[NI
];
5189 || (ecmp (a
, ezero
) == 0)
5197 if (ecmp (a
, ezero
) == 0)
5199 mtherr ("eremain", SING
);
5205 eiremain (den
, num
);
5206 /* Sign of remainder = sign of quotient */
5214 /* Return quotient of exploded e-types NUM / DEN in EQUOT,
5215 remainder in NUM. */
5219 unsigned EMUSHORT den
[], num
[];
5222 unsigned EMUSHORT j
;
5225 ld
-= enormlz (den
);
5227 ln
-= enormlz (num
);
5231 if (ecmpm (den
, num
) <= 0)
5243 emdnorm (num
, 0, 0, ln
, 0);
5246 /* Report an error condition CODE encountered in function NAME.
5247 CODE is one of the following:
5249 Mnemonic Value Significance
5251 DOMAIN 1 argument domain error
5252 SING 2 function singularity
5253 OVERFLOW 3 overflow range error
5254 UNDERFLOW 4 underflow range error
5255 TLOSS 5 total loss of precision
5256 PLOSS 6 partial loss of precision
5257 INVALID 7 NaN - producing operation
5258 EDOM 33 Unix domain error code
5259 ERANGE 34 Unix range error code
5261 The order of appearance of the following messages is bound to the
5262 error codes defined above. */
5265 static char *ermsg
[NMSGS
] =
5267 "unknown", /* error code 0 */
5268 "domain", /* error code 1 */
5269 "singularity", /* et seq. */
5272 "total loss of precision",
5273 "partial loss of precision",
5287 /* The string passed by the calling program is supposed to be the
5288 name of the function in which the error occurred.
5289 The code argument selects which error message string will be printed. */
5291 if ((code
<= 0) || (code
>= NMSGS
))
5293 sprintf (errstr
, " %s %s error", name
, ermsg
[code
]);
5296 /* Set global error message word */
5301 /* Convert DEC double precision D to e type E. */
5305 unsigned EMUSHORT
*d
;
5306 unsigned EMUSHORT
*e
;
5308 unsigned EMUSHORT y
[NI
];
5309 register unsigned EMUSHORT r
, *p
;
5311 ecleaz (y
); /* start with a zero */
5312 p
= y
; /* point to our number */
5313 r
= *d
; /* get DEC exponent word */
5314 if (*d
& (unsigned int) 0x8000)
5315 *p
= 0xffff; /* fill in our sign */
5316 ++p
; /* bump pointer to our exponent word */
5317 r
&= 0x7fff; /* strip the sign bit */
5318 if (r
== 0) /* answer = 0 if high order DEC word = 0 */
5322 r
>>= 7; /* shift exponent word down 7 bits */
5323 r
+= EXONE
- 0201; /* subtract DEC exponent offset */
5324 /* add our e type exponent offset */
5325 *p
++ = r
; /* to form our exponent */
5327 r
= *d
++; /* now do the high order mantissa */
5328 r
&= 0177; /* strip off the DEC exponent and sign bits */
5329 r
|= 0200; /* the DEC understood high order mantissa bit */
5330 *p
++ = r
; /* put result in our high guard word */
5332 *p
++ = *d
++; /* fill in the rest of our mantissa */
5336 eshdn8 (y
); /* shift our mantissa down 8 bits */
5341 /* Convert e type X to DEC double precision D. */
5345 unsigned EMUSHORT
*x
, *d
;
5347 unsigned EMUSHORT xi
[NI
];
5352 /* Adjust exponent for offsets. */
5353 exp
= (EMULONG
) xi
[E
] - (EXONE
- 0201);
5354 /* Round off to nearest or even. */
5357 emdnorm (xi
, 0, 0, exp
, 64);
5362 /* Convert exploded e-type X, that has already been rounded to
5363 56-bit precision, to DEC format double Y. */
5367 unsigned EMUSHORT
*x
, *y
;
5369 unsigned EMUSHORT i
;
5370 unsigned EMUSHORT
*p
;
5409 /* Convert IBM single/double precision to e type. */
5413 unsigned EMUSHORT
*d
;
5414 unsigned EMUSHORT
*e
;
5415 enum machine_mode mode
;
5417 unsigned EMUSHORT y
[NI
];
5418 register unsigned EMUSHORT r
, *p
;
5421 ecleaz (y
); /* start with a zero */
5422 p
= y
; /* point to our number */
5423 r
= *d
; /* get IBM exponent word */
5424 if (*d
& (unsigned int) 0x8000)
5425 *p
= 0xffff; /* fill in our sign */
5426 ++p
; /* bump pointer to our exponent word */
5427 r
&= 0x7f00; /* strip the sign bit */
5428 r
>>= 6; /* shift exponent word down 6 bits */
5429 /* in fact shift by 8 right and 2 left */
5430 r
+= EXONE
- (0x41 << 2); /* subtract IBM exponent offset */
5431 /* add our e type exponent offset */
5432 *p
++ = r
; /* to form our exponent */
5434 *p
++ = *d
++ & 0xff; /* now do the high order mantissa */
5435 /* strip off the IBM exponent and sign bits */
5436 if (mode
!= SFmode
) /* there are only 2 words in SFmode */
5438 *p
++ = *d
++; /* fill in the rest of our mantissa */
5443 if (y
[M
] == 0 && y
[M
+1] == 0 && y
[M
+2] == 0 && y
[M
+3] == 0)
5446 y
[E
] -= 5 + enormlz (y
); /* now normalise the mantissa */
5447 /* handle change in RADIX */
5453 /* Convert e type to IBM single/double precision. */
5457 unsigned EMUSHORT
*x
, *d
;
5458 enum machine_mode mode
;
5460 unsigned EMUSHORT xi
[NI
];
5465 exp
= (EMULONG
) xi
[E
] - (EXONE
- (0x41 << 2)); /* adjust exponent for offsets */
5466 /* round off to nearest or even */
5469 emdnorm (xi
, 0, 0, exp
, 64);
5471 toibm (xi
, d
, mode
);
5476 unsigned EMUSHORT
*x
, *y
;
5477 enum machine_mode mode
;
5479 unsigned EMUSHORT i
;
5480 unsigned EMUSHORT
*p
;
5528 /* Output a binary NaN bit pattern in the target machine's format. */
5530 /* If special NaN bit patterns are required, define them in tm.h
5531 as arrays of unsigned 16-bit shorts. Otherwise, use the default
5537 unsigned EMUSHORT TFbignan
[8] =
5538 {0x7fff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff};
5539 unsigned EMUSHORT TFlittlenan
[8] = {0, 0, 0, 0, 0, 0, 0x8000, 0xffff};
5547 unsigned EMUSHORT XFbignan
[6] =
5548 {0x7fff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff};
5549 unsigned EMUSHORT XFlittlenan
[6] = {0, 0, 0, 0xc000, 0xffff, 0};
5557 unsigned EMUSHORT DFbignan
[4] = {0x7fff, 0xffff, 0xffff, 0xffff};
5558 unsigned EMUSHORT DFlittlenan
[4] = {0, 0, 0, 0xfff8};
5566 unsigned EMUSHORT SFbignan
[2] = {0x7fff, 0xffff};
5567 unsigned EMUSHORT SFlittlenan
[2] = {0, 0xffc0};
5573 make_nan (nan
, sign
, mode
)
5574 unsigned EMUSHORT
*nan
;
5576 enum machine_mode mode
;
5579 unsigned EMUSHORT
*p
;
5583 /* Possibly the `reserved operand' patterns on a VAX can be
5584 used like NaN's, but probably not in the same way as IEEE. */
5585 #if !defined(DEC) && !defined(IBM)
5588 if (REAL_WORDS_BIG_ENDIAN
)
5595 if (REAL_WORDS_BIG_ENDIAN
)
5602 if (REAL_WORDS_BIG_ENDIAN
)
5610 if (REAL_WORDS_BIG_ENDIAN
)
5619 if (REAL_WORDS_BIG_ENDIAN
)
5620 *nan
++ = (sign
<< 15) | *p
++;
5623 if (! REAL_WORDS_BIG_ENDIAN
)
5624 *nan
= (sign
<< 15) | *p
;
5627 /* Convert an SFmode target `float' value to a REAL_VALUE_TYPE.
5628 This is the inverse of the function `etarsingle' invoked by
5629 REAL_VALUE_TO_TARGET_SINGLE. */
5632 ereal_from_float (f
)
5636 unsigned EMUSHORT s
[2];
5637 unsigned EMUSHORT e
[NE
];
5639 /* Convert 32 bit integer to array of 16 bit pieces in target machine order.
5640 This is the inverse operation to what the function `endian' does. */
5641 if (REAL_WORDS_BIG_ENDIAN
)
5643 s
[0] = (unsigned EMUSHORT
) (f
>> 16);
5644 s
[1] = (unsigned EMUSHORT
) f
;
5648 s
[0] = (unsigned EMUSHORT
) f
;
5649 s
[1] = (unsigned EMUSHORT
) (f
>> 16);
5651 /* Convert and promote the target float to E-type. */
5653 /* Output E-type to REAL_VALUE_TYPE. */
5659 /* Convert a DFmode target `double' value to a REAL_VALUE_TYPE.
5660 This is the inverse of the function `etardouble' invoked by
5661 REAL_VALUE_TO_TARGET_DOUBLE.
5663 The DFmode is stored as an array of HOST_WIDE_INT in the target's
5664 data format, with no holes in the bit packing. The first element
5665 of the input array holds the bits that would come first in the
5666 target computer's memory. */
5669 ereal_from_double (d
)
5673 unsigned EMUSHORT s
[4];
5674 unsigned EMUSHORT e
[NE
];
5676 /* Convert array of HOST_WIDE_INT to equivalent array of 16-bit pieces. */
5677 if (REAL_WORDS_BIG_ENDIAN
)
5679 s
[0] = (unsigned EMUSHORT
) (d
[0] >> 16);
5680 s
[1] = (unsigned EMUSHORT
) d
[0];
5681 #if HOST_BITS_PER_WIDE_INT == 32
5682 s
[2] = (unsigned EMUSHORT
) (d
[1] >> 16);
5683 s
[3] = (unsigned EMUSHORT
) d
[1];
5685 /* In this case the entire target double is contained in the
5686 first array element. The second element of the input is
5688 s
[2] = (unsigned EMUSHORT
) (d
[0] >> 48);
5689 s
[3] = (unsigned EMUSHORT
) (d
[0] >> 32);
5694 /* Target float words are little-endian. */
5695 s
[0] = (unsigned EMUSHORT
) d
[0];
5696 s
[1] = (unsigned EMUSHORT
) (d
[0] >> 16);
5697 #if HOST_BITS_PER_WIDE_INT == 32
5698 s
[2] = (unsigned EMUSHORT
) d
[1];
5699 s
[3] = (unsigned EMUSHORT
) (d
[1] >> 16);
5701 s
[2] = (unsigned EMUSHORT
) (d
[0] >> 32);
5702 s
[3] = (unsigned EMUSHORT
) (d
[0] >> 48);
5705 /* Convert target double to E-type. */
5707 /* Output E-type to REAL_VALUE_TYPE. */
5713 /* Convert target computer unsigned 64-bit integer to e-type.
5714 The endian-ness of DImode follows the convention for integers,
5715 so we use WORDS_BIG_ENDIAN here, not REAL_WORDS_BIG_ENDIAN. */
5719 unsigned EMUSHORT
*di
; /* Address of the 64-bit int. */
5720 unsigned EMUSHORT
*e
;
5722 unsigned EMUSHORT yi
[NI
];
5726 if (WORDS_BIG_ENDIAN
)
5728 for (k
= M
; k
< M
+ 4; k
++)
5733 for (k
= M
+ 3; k
>= M
; k
--)
5736 yi
[E
] = EXONE
+ 47; /* exponent if normalize shift count were 0 */
5737 if ((k
= enormlz (yi
)) > NBITS
)/* normalize the significand */
5738 ecleaz (yi
); /* it was zero */
5740 yi
[E
] -= (unsigned EMUSHORT
) k
;/* subtract shift count from exponent */
5744 /* Convert target computer signed 64-bit integer to e-type. */
5748 unsigned EMUSHORT
*di
; /* Address of the 64-bit int. */
5749 unsigned EMUSHORT
*e
;
5751 unsigned EMULONG acc
;
5752 unsigned EMUSHORT yi
[NI
];
5753 unsigned EMUSHORT carry
;
5757 if (WORDS_BIG_ENDIAN
)
5759 for (k
= M
; k
< M
+ 4; k
++)
5764 for (k
= M
+ 3; k
>= M
; k
--)
5767 /* Take absolute value */
5773 for (k
= M
+ 3; k
>= M
; k
--)
5775 acc
= (unsigned EMULONG
) (~yi
[k
] & 0xffff) + carry
;
5782 yi
[E
] = EXONE
+ 47; /* exponent if normalize shift count were 0 */
5783 if ((k
= enormlz (yi
)) > NBITS
)/* normalize the significand */
5784 ecleaz (yi
); /* it was zero */
5786 yi
[E
] -= (unsigned EMUSHORT
) k
;/* subtract shift count from exponent */
5793 /* Convert e-type to unsigned 64-bit int. */
5797 unsigned EMUSHORT
*x
;
5798 unsigned EMUSHORT
*i
;
5800 unsigned EMUSHORT xi
[NI
];
5809 k
= (int) xi
[E
] - (EXONE
- 1);
5812 for (j
= 0; j
< 4; j
++)
5818 for (j
= 0; j
< 4; j
++)
5821 warning ("overflow on truncation to integer");
5826 /* Shift more than 16 bits: first shift up k-16 mod 16,
5827 then shift up by 16's. */
5828 j
= k
- ((k
>> 4) << 4);
5832 if (WORDS_BIG_ENDIAN
)
5843 if (WORDS_BIG_ENDIAN
)
5848 while ((k
-= 16) > 0);
5852 /* shift not more than 16 bits */
5857 if (WORDS_BIG_ENDIAN
)
5876 /* Convert e-type to signed 64-bit int. */
5880 unsigned EMUSHORT
*x
;
5881 unsigned EMUSHORT
*i
;
5883 unsigned EMULONG acc
;
5884 unsigned EMUSHORT xi
[NI
];
5885 unsigned EMUSHORT carry
;
5886 unsigned EMUSHORT
*isave
;
5890 k
= (int) xi
[E
] - (EXONE
- 1);
5893 for (j
= 0; j
< 4; j
++)
5899 for (j
= 0; j
< 4; j
++)
5902 warning ("overflow on truncation to integer");
5908 /* Shift more than 16 bits: first shift up k-16 mod 16,
5909 then shift up by 16's. */
5910 j
= k
- ((k
>> 4) << 4);
5914 if (WORDS_BIG_ENDIAN
)
5925 if (WORDS_BIG_ENDIAN
)
5930 while ((k
-= 16) > 0);
5934 /* shift not more than 16 bits */
5937 if (WORDS_BIG_ENDIAN
)
5953 /* Negate if negative */
5957 if (WORDS_BIG_ENDIAN
)
5959 for (k
= 0; k
< 4; k
++)
5961 acc
= (unsigned EMULONG
) (~(*isave
) & 0xffff) + carry
;
5962 if (WORDS_BIG_ENDIAN
)
5974 /* Longhand square root routine. */
5977 static int esqinited
= 0;
5978 static unsigned short sqrndbit
[NI
];
5982 unsigned EMUSHORT
*x
, *y
;
5984 unsigned EMUSHORT temp
[NI
], num
[NI
], sq
[NI
], xx
[NI
];
5986 int i
, j
, k
, n
, nlups
;
5991 sqrndbit
[NI
- 2] = 1;
5994 /* Check for arg <= 0 */
5995 i
= ecmp (x
, ezero
);
6000 mtherr ("esqrt", DOMAIN
);
6016 /* Bring in the arg and renormalize if it is denormal. */
6018 m
= (EMULONG
) xx
[1]; /* local long word exponent */
6022 /* Divide exponent by 2 */
6024 exp
= (unsigned short) ((m
/ 2) + 0x3ffe);
6026 /* Adjust if exponent odd */
6036 n
= 8; /* get 8 bits of result per inner loop */
6042 /* bring in next word of arg */
6044 num
[NI
- 1] = xx
[j
+ 3];
6045 /* Do additional bit on last outer loop, for roundoff. */
6048 for (i
= 0; i
< n
; i
++)
6050 /* Next 2 bits of arg */
6053 /* Shift up answer */
6055 /* Make trial divisor */
6056 for (k
= 0; k
< NI
; k
++)
6059 eaddm (sqrndbit
, temp
);
6060 /* Subtract and insert answer bit if it goes in */
6061 if (ecmpm (temp
, num
) <= 0)
6071 /* Adjust for extra, roundoff loop done. */
6072 exp
+= (NBITS
- 1) - rndprc
;
6074 /* Sticky bit = 1 if the remainder is nonzero. */
6076 for (i
= 3; i
< NI
; i
++)
6079 /* Renormalize and round off. */
6080 emdnorm (sq
, k
, 0, exp
, 64);
6083 #endif /* EMU_NON_COMPILE not defined */
6085 /* Return the binary precision of the significand for a given
6086 floating point mode. The mode can hold an integer value
6087 that many bits wide, without losing any bits. */
6090 significand_size (mode
)
6091 enum machine_mode mode
;
6100 #if TARGET_FLOAT_FORMAT == IEEE_FLOAT_FORMAT
6103 #if TARGET_FLOAT_FORMAT == IBM_FLOAT_FORMAT
6106 #if TARGET_FLOAT_FORMAT == VAX_FLOAT_FORMAT
This page took 0.370396 seconds and 5 git commands to generate.