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, 94-98, 1999 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. */
28 /* To enable support of XFmode extended real floating point, define
29 LONG_DOUBLE_TYPE_SIZE 96 in the tm.h file (m68k.h or i386.h).
31 To support cross compilation between IEEE, VAX and IBM floating
32 point formats, define REAL_ARITHMETIC in the tm.h file.
34 In either case the machine files (tm.h) must not contain any code
35 that tries to use host floating point arithmetic to convert
36 REAL_VALUE_TYPEs from `double' to `float', pass them to fprintf,
37 etc. In cross-compile situations a REAL_VALUE_TYPE may not
38 be intelligible to the host computer's native arithmetic.
40 The emulator defaults to the host's floating point format so that
41 its decimal conversion functions can be used if desired (see
44 The first part of this file interfaces gcc to a floating point
45 arithmetic suite that was not written with gcc in mind. Avoid
46 changing the low-level arithmetic routines unless you have suitable
47 test programs available. A special version of the PARANOIA floating
48 point arithmetic tester, modified for this purpose, can be found on
49 usc.edu: /pub/C-numanal/ieeetest.zoo. Other tests, and libraries of
50 XFmode and TFmode transcendental functions, can be obtained by ftp from
51 netlib.att.com: netlib/cephes. */
53 /* Type of computer arithmetic.
54 Only one of DEC, IBM, IEEE, C4X, or UNK should get defined.
56 `IEEE', when REAL_WORDS_BIG_ENDIAN is non-zero, refers generically
57 to big-endian IEEE floating-point data structure. This definition
58 should work in SFmode `float' type and DFmode `double' type on
59 virtually all big-endian IEEE machines. If LONG_DOUBLE_TYPE_SIZE
60 has been defined to be 96, then IEEE also invokes the particular
61 XFmode (`long double' type) data structure used by the Motorola
62 680x0 series processors.
64 `IEEE', when REAL_WORDS_BIG_ENDIAN is zero, refers generally to
65 little-endian IEEE machines. In this case, if LONG_DOUBLE_TYPE_SIZE
66 has been defined to be 96, then IEEE also invokes the particular
67 XFmode `long double' data structure used by the Intel 80x86 series
70 `DEC' refers specifically to the Digital Equipment Corp PDP-11
71 and VAX floating point data structure. This model currently
72 supports no type wider than DFmode.
74 `IBM' refers specifically to the IBM System/370 and compatible
75 floating point data structure. This model currently supports
76 no type wider than DFmode. The IBM conversions were contributed by
77 frank@atom.ansto.gov.au (Frank Crawford).
79 `C4X' refers specifically to the floating point format used on
80 Texas Instruments TMS320C3x and TMS320C4x digital signal
81 processors. This supports QFmode (32-bit float, double) and HFmode
82 (40-bit long double) where BITS_PER_BYTE is 32. Unlike IEEE
83 floats, C4x floats are not rounded to be even. The C4x conversions
84 were contributed by m.hayes@elec.canterbury.ac.nz (Michael Hayes) and
85 Haj.Ten.Brugge@net.HCC.nl (Herman ten Brugge).
87 If LONG_DOUBLE_TYPE_SIZE = 64 (the default, unless tm.h defines it)
88 then `long double' and `double' are both implemented, but they
89 both mean DFmode. In this case, the software floating-point
90 support available here is activated by writing
91 #define REAL_ARITHMETIC
94 The case LONG_DOUBLE_TYPE_SIZE = 128 activates TFmode support
95 and may deactivate XFmode since `long double' is used to refer
98 The macros FLOAT_WORDS_BIG_ENDIAN, HOST_FLOAT_WORDS_BIG_ENDIAN,
99 contributed by Richard Earnshaw <Richard.Earnshaw@cl.cam.ac.uk>,
100 separate the floating point unit's endian-ness from that of
101 the integer addressing. This permits one to define a big-endian
102 FPU on a little-endian machine (e.g., ARM). An extension to
103 BYTES_BIG_ENDIAN may be required for some machines in the future.
104 These optional macros may be defined in tm.h. In real.h, they
105 default to WORDS_BIG_ENDIAN, etc., so there is no need to define
106 them for any normal host or target machine on which the floats
107 and the integers have the same endian-ness. */
110 /* The following converts gcc macros into the ones used by this file. */
112 /* REAL_ARITHMETIC defined means that macros in real.h are
113 defined to call emulator functions. */
114 #ifdef REAL_ARITHMETIC
116 #if TARGET_FLOAT_FORMAT == VAX_FLOAT_FORMAT
117 /* PDP-11, Pro350, VAX: */
119 #else /* it's not VAX */
120 #if TARGET_FLOAT_FORMAT == IBM_FLOAT_FORMAT
121 /* IBM System/370 style */
123 #else /* it's also not an IBM */
124 #if TARGET_FLOAT_FORMAT == C4X_FLOAT_FORMAT
125 /* TMS320C3x/C4x style */
127 #else /* it's also not a C4X */
128 #if TARGET_FLOAT_FORMAT == IEEE_FLOAT_FORMAT
130 #else /* it's not IEEE either */
131 /* UNKnown arithmetic. We don't support this and can't go on. */
132 unknown arithmetic type
134 #endif /* not IEEE */
139 #define REAL_WORDS_BIG_ENDIAN FLOAT_WORDS_BIG_ENDIAN
142 /* REAL_ARITHMETIC not defined means that the *host's* data
143 structure will be used. It may differ by endian-ness from the
144 target machine's structure and will get its ends swapped
145 accordingly (but not here). Probably only the decimal <-> binary
146 functions in this file will actually be used in this case. */
148 #if HOST_FLOAT_FORMAT == VAX_FLOAT_FORMAT
150 #else /* it's not VAX */
151 #if HOST_FLOAT_FORMAT == IBM_FLOAT_FORMAT
152 /* IBM System/370 style */
154 #else /* it's also not an IBM */
155 #if HOST_FLOAT_FORMAT == IEEE_FLOAT_FORMAT
157 #else /* it's not IEEE either */
158 unknown arithmetic type
160 #endif /* not IEEE */
164 #define REAL_WORDS_BIG_ENDIAN HOST_FLOAT_WORDS_BIG_ENDIAN
166 #endif /* REAL_ARITHMETIC not defined */
168 /* Define INFINITY for support of infinity.
169 Define NANS for support of Not-a-Number's (NaN's). */
170 #if !defined(DEC) && !defined(IBM) && !defined(C4X)
175 /* Support of NaNs requires support of infinity. */
182 /* Find a host integer type that is at least 16 bits wide,
183 and another type at least twice whatever that size is. */
185 #if HOST_BITS_PER_CHAR >= 16
186 #define EMUSHORT char
187 #define EMUSHORT_SIZE HOST_BITS_PER_CHAR
188 #define EMULONG_SIZE (2 * HOST_BITS_PER_CHAR)
190 #if HOST_BITS_PER_SHORT >= 16
191 #define EMUSHORT short
192 #define EMUSHORT_SIZE HOST_BITS_PER_SHORT
193 #define EMULONG_SIZE (2 * HOST_BITS_PER_SHORT)
195 #if HOST_BITS_PER_INT >= 16
197 #define EMUSHORT_SIZE HOST_BITS_PER_INT
198 #define EMULONG_SIZE (2 * HOST_BITS_PER_INT)
200 #if HOST_BITS_PER_LONG >= 16
201 #define EMUSHORT long
202 #define EMUSHORT_SIZE HOST_BITS_PER_LONG
203 #define EMULONG_SIZE (2 * HOST_BITS_PER_LONG)
205 /* You will have to modify this program to have a smaller unit size. */
206 #define EMU_NON_COMPILE
212 #if HOST_BITS_PER_SHORT >= EMULONG_SIZE
213 #define EMULONG short
215 #if HOST_BITS_PER_INT >= EMULONG_SIZE
218 #if HOST_BITS_PER_LONG >= EMULONG_SIZE
221 #if HOST_BITS_PER_LONGLONG >= EMULONG_SIZE
222 #define EMULONG long long int
224 /* You will have to modify this program to have a smaller unit size. */
225 #define EMU_NON_COMPILE
232 /* The host interface doesn't work if no 16-bit size exists. */
233 #if EMUSHORT_SIZE != 16
234 #define EMU_NON_COMPILE
237 /* OK to continue compilation. */
238 #ifndef EMU_NON_COMPILE
240 /* Construct macros to translate between REAL_VALUE_TYPE and e type.
241 In GET_REAL and PUT_REAL, r and e are pointers.
242 A REAL_VALUE_TYPE is guaranteed to occupy contiguous locations
243 in memory, with no holes. */
245 #if LONG_DOUBLE_TYPE_SIZE == 96
246 /* Number of 16 bit words in external e type format */
248 #define MAXDECEXP 4932
249 #define MINDECEXP -4956
250 #define GET_REAL(r,e) bcopy ((char *) r, (char *) e, 2*NE)
251 #define PUT_REAL(e,r) \
253 if (2*NE < sizeof(*r)) \
254 bzero((char *)r, sizeof(*r)); \
255 bcopy ((char *) e, (char *) r, 2*NE); \
257 #else /* no XFmode */
258 #if LONG_DOUBLE_TYPE_SIZE == 128
260 #define MAXDECEXP 4932
261 #define MINDECEXP -4977
262 #define GET_REAL(r,e) bcopy ((char *) r, (char *) e, 2*NE)
263 #define PUT_REAL(e,r) bcopy ((char *) e, (char *) r, 2*NE)
266 #define MAXDECEXP 4932
267 #define MINDECEXP -4956
268 #ifdef REAL_ARITHMETIC
269 /* Emulator uses target format internally
270 but host stores it in host endian-ness. */
272 #define GET_REAL(r,e) \
274 if (HOST_FLOAT_WORDS_BIG_ENDIAN == REAL_WORDS_BIG_ENDIAN) \
275 e53toe ((unsigned EMUSHORT *) (r), (e)); \
278 unsigned EMUSHORT w[4]; \
279 w[3] = ((EMUSHORT *) r)[0]; \
280 w[2] = ((EMUSHORT *) r)[1]; \
281 w[1] = ((EMUSHORT *) r)[2]; \
282 w[0] = ((EMUSHORT *) r)[3]; \
287 #define PUT_REAL(e,r) \
289 if (HOST_FLOAT_WORDS_BIG_ENDIAN == REAL_WORDS_BIG_ENDIAN) \
290 etoe53 ((e), (unsigned EMUSHORT *) (r)); \
293 unsigned EMUSHORT w[4]; \
295 *((EMUSHORT *) r) = w[3]; \
296 *((EMUSHORT *) r + 1) = w[2]; \
297 *((EMUSHORT *) r + 2) = w[1]; \
298 *((EMUSHORT *) r + 3) = w[0]; \
302 #else /* not REAL_ARITHMETIC */
304 /* emulator uses host format */
305 #define GET_REAL(r,e) e53toe ((unsigned EMUSHORT *) (r), (e))
306 #define PUT_REAL(e,r) etoe53 ((e), (unsigned EMUSHORT *) (r))
308 #endif /* not REAL_ARITHMETIC */
309 #endif /* not TFmode */
310 #endif /* not XFmode */
313 /* Number of 16 bit words in internal format */
316 /* Array offset to exponent */
319 /* Array offset to high guard word */
322 /* Number of bits of precision */
323 #define NBITS ((NI-4)*16)
325 /* Maximum number of decimal digits in ASCII conversion
328 #define NDEC (NBITS*8/27)
330 /* The exponent of 1.0 */
331 #define EXONE (0x3fff)
333 extern int extra_warnings;
334 extern unsigned EMUSHORT ezero[], ehalf[], eone[], etwo[];
335 extern unsigned EMUSHORT elog2[], esqrt2[];
337 static void endian PROTO((unsigned EMUSHORT *, long *,
339 static void eclear PROTO((unsigned EMUSHORT *));
340 static void emov PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
342 static void eabs PROTO((unsigned EMUSHORT *));
344 static void eneg PROTO((unsigned EMUSHORT *));
345 static int eisneg PROTO((unsigned EMUSHORT *));
346 static int eisinf PROTO((unsigned EMUSHORT *));
347 static int eisnan PROTO((unsigned EMUSHORT *));
348 static void einfin PROTO((unsigned EMUSHORT *));
349 static void enan PROTO((unsigned EMUSHORT *, int));
350 static void emovi PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
351 static void emovo PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
352 static void ecleaz PROTO((unsigned EMUSHORT *));
353 static void ecleazs PROTO((unsigned EMUSHORT *));
354 static void emovz PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
355 static void einan PROTO((unsigned EMUSHORT *));
356 static int eiisnan PROTO((unsigned EMUSHORT *));
357 static int eiisneg PROTO((unsigned EMUSHORT *));
359 static void eiinfin PROTO((unsigned EMUSHORT *));
361 static int eiisinf PROTO((unsigned EMUSHORT *));
362 static int ecmpm PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
363 static void eshdn1 PROTO((unsigned EMUSHORT *));
364 static void eshup1 PROTO((unsigned EMUSHORT *));
365 static void eshdn8 PROTO((unsigned EMUSHORT *));
366 static void eshup8 PROTO((unsigned EMUSHORT *));
367 static void eshup6 PROTO((unsigned EMUSHORT *));
368 static void eshdn6 PROTO((unsigned EMUSHORT *));
369 static void eaddm PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
\f
370 static void esubm PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
371 static void m16m PROTO((unsigned int, unsigned short *,
373 static int edivm PROTO((unsigned short *, unsigned short *));
374 static int emulm PROTO((unsigned short *, unsigned short *));
375 static void emdnorm PROTO((unsigned EMUSHORT *, int, int, EMULONG, int));
376 static void esub PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *,
377 unsigned EMUSHORT *));
378 static void eadd PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *,
379 unsigned EMUSHORT *));
380 static void eadd1 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *,
381 unsigned EMUSHORT *));
382 static void ediv PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *,
383 unsigned EMUSHORT *));
384 static void emul PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *,
385 unsigned EMUSHORT *));
386 static void e53toe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
387 static void e64toe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
388 static void e113toe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
389 static void e24toe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
390 static void etoe113 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
391 static void toe113 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
392 static void etoe64 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
393 static void toe64 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
394 static void etoe53 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
395 static void toe53 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
396 static void etoe24 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
397 static void toe24 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
398 static int ecmp PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
400 static void eround PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
402 static void ltoe PROTO((HOST_WIDE_INT *, unsigned EMUSHORT *));
403 static void ultoe PROTO((unsigned HOST_WIDE_INT *, unsigned EMUSHORT *));
404 static void eifrac PROTO((unsigned EMUSHORT *, HOST_WIDE_INT *,
405 unsigned EMUSHORT *));
406 static void euifrac PROTO((unsigned EMUSHORT *, unsigned HOST_WIDE_INT *,
407 unsigned EMUSHORT *));
408 static int eshift PROTO((unsigned EMUSHORT *, int));
409 static int enormlz PROTO((unsigned EMUSHORT *));
411 static void e24toasc PROTO((unsigned EMUSHORT *, char *, int));
412 static void e53toasc PROTO((unsigned EMUSHORT *, char *, int));
413 static void e64toasc PROTO((unsigned EMUSHORT *, char *, int));
414 static void e113toasc PROTO((unsigned EMUSHORT *, char *, int));
416 static void etoasc PROTO((unsigned EMUSHORT *, char *, int));
417 static void asctoe24 PROTO((const char *, unsigned EMUSHORT *));
418 static void asctoe53 PROTO((const char *, unsigned EMUSHORT *));
419 static void asctoe64 PROTO((const char *, unsigned EMUSHORT *));
420 static void asctoe113 PROTO((const char *, unsigned EMUSHORT *));
421 static void asctoe PROTO((const char *, unsigned EMUSHORT *));
422 static void asctoeg PROTO((const char *, unsigned EMUSHORT *, int));
423 static void efloor PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
425 static void efrexp PROTO((unsigned EMUSHORT *, int *,
426 unsigned EMUSHORT *));
428 static void eldexp PROTO((unsigned EMUSHORT *, int, unsigned EMUSHORT *));
430 static void eremain PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *,
431 unsigned EMUSHORT *));
433 static void eiremain PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
434 static void mtherr PROTO((const char *, int));
436 static void dectoe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
437 static void etodec PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
438 static void todec PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
441 static void ibmtoe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *,
443 static void etoibm PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *,
445 static void toibm PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *,
449 static void c4xtoe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *,
451 static void etoc4x PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *,
453 static void toc4x PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *,
456 static void make_nan PROTO((unsigned EMUSHORT *, int, enum machine_mode));
458 static void uditoe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
459 static void ditoe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
460 static void etoudi PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
461 static void etodi PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
462 static void esqrt PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
465 /* Copy 32-bit numbers obtained from array containing 16-bit numbers,
466 swapping ends if required, into output array of longs. The
467 result is normally passed to fprintf by the ASM_OUTPUT_ macros. */
471 unsigned EMUSHORT e[];
473 enum machine_mode mode;
477 if (REAL_WORDS_BIG_ENDIAN)
482 /* Swap halfwords in the fourth long. */
483 th = (unsigned long) e[6] & 0xffff;
484 t = (unsigned long) e[7] & 0xffff;
489 /* Swap halfwords in the third long. */
490 th = (unsigned long) e[4] & 0xffff;
491 t = (unsigned long) e[5] & 0xffff;
494 /* fall into the double case */
497 /* Swap halfwords in the second word. */
498 th = (unsigned long) e[2] & 0xffff;
499 t = (unsigned long) e[3] & 0xffff;
502 /* fall into the float case */
506 /* Swap halfwords in the first word. */
507 th = (unsigned long) e[0] & 0xffff;
508 t = (unsigned long) e[1] & 0xffff;
519 /* Pack the output array without swapping. */
524 /* Pack the fourth long. */
525 th = (unsigned long) e[7] & 0xffff;
526 t = (unsigned long) e[6] & 0xffff;
531 /* Pack the third long.
532 Each element of the input REAL_VALUE_TYPE array has 16 useful bits
534 th = (unsigned long) e[5] & 0xffff;
535 t = (unsigned long) e[4] & 0xffff;
538 /* fall into the double case */
541 /* Pack the second long */
542 th = (unsigned long) e[3] & 0xffff;
543 t = (unsigned long) e[2] & 0xffff;
546 /* fall into the float case */
550 /* Pack the first long */
551 th = (unsigned long) e[1] & 0xffff;
552 t = (unsigned long) e[0] & 0xffff;
564 /* This is the implementation of the REAL_ARITHMETIC macro. */
567 earith (value, icode, r1, r2)
568 REAL_VALUE_TYPE *value;
573 unsigned EMUSHORT d1[NE], d2[NE], v[NE];
579 /* Return NaN input back to the caller. */
582 PUT_REAL (d1, value);
587 PUT_REAL (d2, value);
591 code = (enum tree_code) icode;
599 esub (d2, d1, v); /* d1 - d2 */
607 #ifndef REAL_INFINITY
608 if (ecmp (d2, ezero) == 0)
611 enan (v, eisneg (d1) ^ eisneg (d2));
618 ediv (d2, d1, v); /* d1/d2 */
621 case MIN_EXPR: /* min (d1,d2) */
622 if (ecmp (d1, d2) < 0)
628 case MAX_EXPR: /* max (d1,d2) */
629 if (ecmp (d1, d2) > 0)
642 /* Truncate REAL_VALUE_TYPE toward zero to signed HOST_WIDE_INT.
643 implements REAL_VALUE_RNDZINT (x) (etrunci (x)). */
649 unsigned EMUSHORT f[NE], g[NE];
665 /* Truncate REAL_VALUE_TYPE toward zero to unsigned HOST_WIDE_INT;
666 implements REAL_VALUE_UNSIGNED_RNDZINT (x) (etruncui (x)). */
672 unsigned EMUSHORT f[NE], g[NE];
674 unsigned HOST_WIDE_INT l;
688 /* This is the REAL_VALUE_ATOF function. It converts a decimal or hexadecimal
689 string to binary, rounding off as indicated by the machine_mode argument.
690 Then it promotes the rounded value to REAL_VALUE_TYPE. */
697 unsigned EMUSHORT tem[NE], e[NE];
740 /* Expansion of REAL_NEGATE. */
746 unsigned EMUSHORT e[NE];
756 /* Round real toward zero to HOST_WIDE_INT;
757 implements REAL_VALUE_FIX (x). */
763 unsigned EMUSHORT f[NE], g[NE];
770 warning ("conversion from NaN to int");
778 /* Round real toward zero to unsigned HOST_WIDE_INT
779 implements REAL_VALUE_UNSIGNED_FIX (x).
780 Negative input returns zero. */
782 unsigned HOST_WIDE_INT
786 unsigned EMUSHORT f[NE], g[NE];
787 unsigned HOST_WIDE_INT l;
793 warning ("conversion from NaN to unsigned int");
802 /* REAL_VALUE_FROM_INT macro. */
805 ereal_from_int (d, i, j, mode)
808 enum machine_mode mode;
810 unsigned EMUSHORT df[NE], dg[NE];
811 HOST_WIDE_INT low, high;
814 if (GET_MODE_CLASS (mode) != MODE_FLOAT)
821 /* complement and add 1 */
828 eldexp (eone, HOST_BITS_PER_WIDE_INT, df);
829 ultoe ((unsigned HOST_WIDE_INT *) &high, dg);
831 ultoe ((unsigned HOST_WIDE_INT *) &low, df);
836 /* A REAL_VALUE_TYPE may not be wide enough to hold the two HOST_WIDE_INTS.
837 Avoid double-rounding errors later by rounding off now from the
838 extra-wide internal format to the requested precision. */
839 switch (GET_MODE_BITSIZE (mode))
869 /* REAL_VALUE_FROM_UNSIGNED_INT macro. */
872 ereal_from_uint (d, i, j, mode)
874 unsigned HOST_WIDE_INT i, j;
875 enum machine_mode mode;
877 unsigned EMUSHORT df[NE], dg[NE];
878 unsigned HOST_WIDE_INT low, high;
880 if (GET_MODE_CLASS (mode) != MODE_FLOAT)
884 eldexp (eone, HOST_BITS_PER_WIDE_INT, df);
890 /* A REAL_VALUE_TYPE may not be wide enough to hold the two HOST_WIDE_INTS.
891 Avoid double-rounding errors later by rounding off now from the
892 extra-wide internal format to the requested precision. */
893 switch (GET_MODE_BITSIZE (mode))
923 /* REAL_VALUE_TO_INT macro. */
926 ereal_to_int (low, high, rr)
927 HOST_WIDE_INT *low, *high;
930 unsigned EMUSHORT d[NE], df[NE], dg[NE], dh[NE];
937 warning ("conversion from NaN to int");
943 /* convert positive value */
950 eldexp (eone, HOST_BITS_PER_WIDE_INT, df);
951 ediv (df, d, dg); /* dg = d / 2^32 is the high word */
952 euifrac (dg, (unsigned HOST_WIDE_INT *) high, dh);
953 emul (df, dh, dg); /* fractional part is the low word */
954 euifrac (dg, (unsigned HOST_WIDE_INT *)low, dh);
957 /* complement and add 1 */
967 /* REAL_VALUE_LDEXP macro. */
974 unsigned EMUSHORT e[NE], y[NE];
987 /* These routines are conditionally compiled because functions
988 of the same names may be defined in fold-const.c. */
990 #ifdef REAL_ARITHMETIC
992 /* Check for infinity in a REAL_VALUE_TYPE. */
998 unsigned EMUSHORT e[NE];
1002 return (eisinf (e));
1008 /* Check whether a REAL_VALUE_TYPE item is a NaN. */
1014 unsigned EMUSHORT e[NE];
1018 return (eisnan (e));
1025 /* Check for a negative REAL_VALUE_TYPE number.
1026 This just checks the sign bit, so that -0 counts as negative. */
1032 return ereal_isneg (x);
1035 /* Expansion of REAL_VALUE_TRUNCATE.
1036 The result is in floating point, rounded to nearest or even. */
1039 real_value_truncate (mode, arg)
1040 enum machine_mode mode;
1041 REAL_VALUE_TYPE arg;
1043 unsigned EMUSHORT e[NE], t[NE];
1089 /* If an unsupported type was requested, presume that
1090 the machine files know something useful to do with
1091 the unmodified value. */
1100 /* Try to change R into its exact multiplicative inverse in machine mode
1101 MODE. Return nonzero function value if successful. */
1104 exact_real_inverse (mode, r)
1105 enum machine_mode mode;
1108 unsigned EMUSHORT e[NE], einv[NE];
1109 REAL_VALUE_TYPE rinv;
1114 /* Test for input in range. Don't transform IEEE special values. */
1115 if (eisinf (e) || eisnan (e) || (ecmp (e, ezero) == 0))
1118 /* Test for a power of 2: all significand bits zero except the MSB.
1119 We are assuming the target has binary (or hex) arithmetic. */
1120 if (e[NE - 2] != 0x8000)
1123 for (i = 0; i < NE - 2; i++)
1129 /* Compute the inverse and truncate it to the required mode. */
1130 ediv (e, eone, einv);
1131 PUT_REAL (einv, &rinv);
1132 rinv = real_value_truncate (mode, rinv);
1134 #ifdef CHECK_FLOAT_VALUE
1135 /* This check is not redundant. It may, for example, flush
1136 a supposedly IEEE denormal value to zero. */
1138 if (CHECK_FLOAT_VALUE (mode, rinv, i))
1141 GET_REAL (&rinv, einv);
1143 /* Check the bits again, because the truncation might have
1144 generated an arbitrary saturation value on overflow. */
1145 if (einv[NE - 2] != 0x8000)
1148 for (i = 0; i < NE - 2; i++)
1154 /* Fail if the computed inverse is out of range. */
1155 if (eisinf (einv) || eisnan (einv) || (ecmp (einv, ezero) == 0))
1158 /* Output the reciprocal and return success flag. */
1162 #endif /* REAL_ARITHMETIC defined */
1164 /* Used for debugging--print the value of R in human-readable format
1173 REAL_VALUE_TO_DECIMAL (r, "%.20g", dstr);
1174 fprintf (stderr, "%s", dstr);
1178 /* The following routines convert REAL_VALUE_TYPE to the various floating
1179 point formats that are meaningful to supported computers.
1181 The results are returned in 32-bit pieces, each piece stored in a `long'.
1182 This is so they can be printed by statements like
1184 fprintf (file, "%lx, %lx", L[0], L[1]);
1186 that will work on both narrow- and wide-word host computers. */
1188 /* Convert R to a 128-bit long double precision value. The output array L
1189 contains four 32-bit pieces of the result, in the order they would appear
1197 unsigned EMUSHORT e[NE];
1201 endian (e, l, TFmode);
1204 /* Convert R to a double extended precision value. The output array L
1205 contains three 32-bit pieces of the result, in the order they would
1206 appear in memory. */
1213 unsigned EMUSHORT e[NE];
1217 endian (e, l, XFmode);
1220 /* Convert R to a double precision value. The output array L contains two
1221 32-bit pieces of the result, in the order they would appear in memory. */
1228 unsigned EMUSHORT e[NE];
1232 endian (e, l, DFmode);
1235 /* Convert R to a single precision float value stored in the least-significant
1236 bits of a `long'. */
1242 unsigned EMUSHORT e[NE];
1247 endian (e, &l, SFmode);
1251 /* Convert X to a decimal ASCII string S for output to an assembly
1252 language file. Note, there is no standard way to spell infinity or
1253 a NaN, so these values may require special treatment in the tm.h
1257 ereal_to_decimal (x, s)
1261 unsigned EMUSHORT e[NE];
1267 /* Compare X and Y. Return 1 if X > Y, 0 if X == Y, -1 if X < Y,
1268 or -2 if either is a NaN. */
1272 REAL_VALUE_TYPE x, y;
1274 unsigned EMUSHORT ex[NE], ey[NE];
1278 return (ecmp (ex, ey));
1281 /* Return 1 if the sign bit of X is set, else return 0. */
1287 unsigned EMUSHORT ex[NE];
1290 return (eisneg (ex));
1293 /* End of REAL_ARITHMETIC interface */
1296 Extended precision IEEE binary floating point arithmetic routines
1298 Numbers are stored in C language as arrays of 16-bit unsigned
1299 short integers. The arguments of the routines are pointers to
1302 External e type data structure, similar to Intel 8087 chip
1303 temporary real format but possibly with a larger significand:
1305 NE-1 significand words (least significant word first,
1306 most significant bit is normally set)
1307 exponent (value = EXONE for 1.0,
1308 top bit is the sign)
1311 Internal exploded e-type data structure of a number (a "word" is 16 bits):
1313 ei[0] sign word (0 for positive, 0xffff for negative)
1314 ei[1] biased exponent (value = EXONE for the number 1.0)
1315 ei[2] high guard word (always zero after normalization)
1317 to ei[NI-2] significand (NI-4 significand words,
1318 most significant word first,
1319 most significant bit is set)
1320 ei[NI-1] low guard word (0x8000 bit is rounding place)
1324 Routines for external format e-type numbers
1326 asctoe (string, e) ASCII string to extended double e type
1327 asctoe64 (string, &d) ASCII string to long double
1328 asctoe53 (string, &d) ASCII string to double
1329 asctoe24 (string, &f) ASCII string to single
1330 asctoeg (string, e, prec) ASCII string to specified precision
1331 e24toe (&f, e) IEEE single precision to e type
1332 e53toe (&d, e) IEEE double precision to e type
1333 e64toe (&d, e) IEEE long double precision to e type
1334 e113toe (&d, e) 128-bit long double precision to e type
1336 eabs (e) absolute value
1338 eadd (a, b, c) c = b + a
1340 ecmp (a, b) Returns 1 if a > b, 0 if a == b,
1341 -1 if a < b, -2 if either a or b is a NaN.
1342 ediv (a, b, c) c = b / a
1343 efloor (a, b) truncate to integer, toward -infinity
1344 efrexp (a, exp, s) extract exponent and significand
1345 eifrac (e, &l, frac) e to HOST_WIDE_INT and e type fraction
1346 euifrac (e, &l, frac) e to unsigned HOST_WIDE_INT and e type fraction
1347 einfin (e) set e to infinity, leaving its sign alone
1348 eldexp (a, n, b) multiply by 2**n
1350 emul (a, b, c) c = b * a
1353 eround (a, b) b = nearest integer value to a
1355 esub (a, b, c) c = b - a
1357 e24toasc (&f, str, n) single to ASCII string, n digits after decimal
1358 e53toasc (&d, str, n) double to ASCII string, n digits after decimal
1359 e64toasc (&d, str, n) 80-bit long double to ASCII string
1360 e113toasc (&d, str, n) 128-bit long double to ASCII string
1362 etoasc (e, str, n) e to ASCII string, n digits after decimal
1363 etoe24 (e, &f) convert e type to IEEE single precision
1364 etoe53 (e, &d) convert e type to IEEE double precision
1365 etoe64 (e, &d) convert e type to IEEE long double precision
1366 ltoe (&l, e) HOST_WIDE_INT to e type
1367 ultoe (&l, e) unsigned HOST_WIDE_INT to e type
1368 eisneg (e) 1 if sign bit of e != 0, else 0
1369 eisinf (e) 1 if e has maximum exponent (non-IEEE)
1370 or is infinite (IEEE)
1371 eisnan (e) 1 if e is a NaN
1374 Routines for internal format exploded e-type numbers
1376 eaddm (ai, bi) add significands, bi = bi + ai
1378 ecleazs (ei) set ei = 0 but leave its sign alone
1379 ecmpm (ai, bi) compare significands, return 1, 0, or -1
1380 edivm (ai, bi) divide significands, bi = bi / ai
1381 emdnorm (ai,l,s,exp) normalize and round off
1382 emovi (a, ai) convert external a to internal ai
1383 emovo (ai, a) convert internal ai to external a
1384 emovz (ai, bi) bi = ai, low guard word of bi = 0
1385 emulm (ai, bi) multiply significands, bi = bi * ai
1386 enormlz (ei) left-justify the significand
1387 eshdn1 (ai) shift significand and guards down 1 bit
1388 eshdn8 (ai) shift down 8 bits
1389 eshdn6 (ai) shift down 16 bits
1390 eshift (ai, n) shift ai n bits up (or down if n < 0)
1391 eshup1 (ai) shift significand and guards up 1 bit
1392 eshup8 (ai) shift up 8 bits
1393 eshup6 (ai) shift up 16 bits
1394 esubm (ai, bi) subtract significands, bi = bi - ai
1395 eiisinf (ai) 1 if infinite
1396 eiisnan (ai) 1 if a NaN
1397 eiisneg (ai) 1 if sign bit of ai != 0, else 0
1398 einan (ai) set ai = NaN
1400 eiinfin (ai) set ai = infinity
1403 The result is always normalized and rounded to NI-4 word precision
1404 after each arithmetic operation.
1406 Exception flags are NOT fully supported.
1408 Signaling NaN's are NOT supported; they are treated the same
1411 Define INFINITY for support of infinity; otherwise a
1412 saturation arithmetic is implemented.
1414 Define NANS for support of Not-a-Number items; otherwise the
1415 arithmetic will never produce a NaN output, and might be confused
1417 If NaN's are supported, the output of `ecmp (a,b)' is -2 if
1418 either a or b is a NaN. This means asking `if (ecmp (a,b) < 0)'
1419 may not be legitimate. Use `if (ecmp (a,b) == -1)' for `less than'
1422 Denormals are always supported here where appropriate (e.g., not
1423 for conversion to DEC numbers). */
1425 /* Definitions for error codes that are passed to the common error handling
1428 For Digital Equipment PDP-11 and VAX computers, certain
1429 IBM systems, and others that use numbers with a 56-bit
1430 significand, the symbol DEC should be defined. In this
1431 mode, most floating point constants are given as arrays
1432 of octal integers to eliminate decimal to binary conversion
1433 errors that might be introduced by the compiler.
1435 For computers, such as IBM PC, that follow the IEEE
1436 Standard for Binary Floating Point Arithmetic (ANSI/IEEE
1437 Std 754-1985), the symbol IEEE should be defined.
1438 These numbers have 53-bit significands. In this mode, constants
1439 are provided as arrays of hexadecimal 16 bit integers.
1440 The endian-ness of generated values is controlled by
1441 REAL_WORDS_BIG_ENDIAN.
1443 To accommodate other types of computer arithmetic, all
1444 constants are also provided in a normal decimal radix
1445 which one can hope are correctly converted to a suitable
1446 format by the available C language compiler. To invoke
1447 this mode, the symbol UNK is defined.
1449 An important difference among these modes is a predefined
1450 set of machine arithmetic constants for each. The numbers
1451 MACHEP (the machine roundoff error), MAXNUM (largest number
1452 represented), and several other parameters are preset by
1453 the configuration symbol. Check the file const.c to
1454 ensure that these values are correct for your computer.
1456 For ANSI C compatibility, define ANSIC equal to 1. Currently
1457 this affects only the atan2 function and others that use it. */
1459 /* Constant definitions for math error conditions. */
1461 #define DOMAIN 1 /* argument domain error */
1462 #define SING 2 /* argument singularity */
1463 #define OVERFLOW 3 /* overflow range error */
1464 #define UNDERFLOW 4 /* underflow range error */
1465 #define TLOSS 5 /* total loss of precision */
1466 #define PLOSS 6 /* partial loss of precision */
1467 #define INVALID 7 /* NaN-producing operation */
1469 /* e type constants used by high precision check routines */
1471 #if LONG_DOUBLE_TYPE_SIZE == 128
1473 unsigned EMUSHORT ezero[NE] =
1474 {0x0000, 0x0000, 0x0000, 0x0000,
1475 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000,};
1476 extern unsigned EMUSHORT ezero[];
1479 unsigned EMUSHORT ehalf[NE] =
1480 {0x0000, 0x0000, 0x0000, 0x0000,
1481 0x0000, 0x0000, 0x0000, 0x0000, 0x8000, 0x3ffe,};
1482 extern unsigned EMUSHORT ehalf[];
1485 unsigned EMUSHORT eone[NE] =
1486 {0x0000, 0x0000, 0x0000, 0x0000,
1487 0x0000, 0x0000, 0x0000, 0x0000, 0x8000, 0x3fff,};
1488 extern unsigned EMUSHORT eone[];
1491 unsigned EMUSHORT etwo[NE] =
1492 {0x0000, 0x0000, 0x0000, 0x0000,
1493 0x0000, 0x0000, 0x0000, 0x0000, 0x8000, 0x4000,};
1494 extern unsigned EMUSHORT etwo[];
1497 unsigned EMUSHORT e32[NE] =
1498 {0x0000, 0x0000, 0x0000, 0x0000,
1499 0x0000, 0x0000, 0x0000, 0x0000, 0x8000, 0x4004,};
1500 extern unsigned EMUSHORT e32[];
1502 /* 6.93147180559945309417232121458176568075500134360255E-1 */
1503 unsigned EMUSHORT elog2[NE] =
1504 {0x40f3, 0xf6af, 0x03f2, 0xb398,
1505 0xc9e3, 0x79ab, 0150717, 0013767, 0130562, 0x3ffe,};
1506 extern unsigned EMUSHORT elog2[];
1508 /* 1.41421356237309504880168872420969807856967187537695E0 */
1509 unsigned EMUSHORT esqrt2[NE] =
1510 {0x1d6f, 0xbe9f, 0x754a, 0x89b3,
1511 0x597d, 0x6484, 0174736, 0171463, 0132404, 0x3fff,};
1512 extern unsigned EMUSHORT esqrt2[];
1514 /* 3.14159265358979323846264338327950288419716939937511E0 */
1515 unsigned EMUSHORT epi[NE] =
1516 {0x2902, 0x1cd1, 0x80dc, 0x628b,
1517 0xc4c6, 0xc234, 0020550, 0155242, 0144417, 0040000,};
1518 extern unsigned EMUSHORT epi[];
1521 /* LONG_DOUBLE_TYPE_SIZE is other than 128 */
1522 unsigned EMUSHORT ezero[NE] =
1523 {0, 0000000, 0000000, 0000000, 0000000, 0000000,};
1524 unsigned EMUSHORT ehalf[NE] =
1525 {0, 0000000, 0000000, 0000000, 0100000, 0x3ffe,};
1526 unsigned EMUSHORT eone[NE] =
1527 {0, 0000000, 0000000, 0000000, 0100000, 0x3fff,};
1528 unsigned EMUSHORT etwo[NE] =
1529 {0, 0000000, 0000000, 0000000, 0100000, 0040000,};
1530 unsigned EMUSHORT e32[NE] =
1531 {0, 0000000, 0000000, 0000000, 0100000, 0040004,};
1532 unsigned EMUSHORT elog2[NE] =
1533 {0xc9e4, 0x79ab, 0150717, 0013767, 0130562, 0x3ffe,};
1534 unsigned EMUSHORT esqrt2[NE] =
1535 {0x597e, 0x6484, 0174736, 0171463, 0132404, 0x3fff,};
1536 unsigned EMUSHORT epi[NE] =
1537 {0xc4c6, 0xc234, 0020550, 0155242, 0144417, 0040000,};
1540 /* Control register for rounding precision.
1541 This can be set to 113 (if NE=10), 80 (if NE=6), 64, 56, 53, or 24 bits. */
1546 /* Clear out entire e-type number X. */
1550 register unsigned EMUSHORT *x;
1554 for (i = 0; i < NE; i++)
1558 /* Move e-type number from A to B. */
1562 register unsigned EMUSHORT *a, *b;
1566 for (i = 0; i < NE; i++)
1572 /* Absolute value of e-type X. */
1576 unsigned EMUSHORT x[];
1578 /* sign is top bit of last word of external format */
1579 x[NE - 1] &= 0x7fff;
1583 /* Negate the e-type number X. */
1587 unsigned EMUSHORT x[];
1590 x[NE - 1] ^= 0x8000; /* Toggle the sign bit */
1593 /* Return 1 if sign bit of e-type number X is nonzero, else zero. */
1597 unsigned EMUSHORT x[];
1600 if (x[NE - 1] & 0x8000)
1606 /* Return 1 if e-type number X is infinity, else return zero. */
1610 unsigned EMUSHORT x[];
1617 if ((x[NE - 1] & 0x7fff) == 0x7fff)
1623 /* Check if e-type number is not a number. The bit pattern is one that we
1624 defined, so we know for sure how to detect it. */
1628 unsigned EMUSHORT x[];
1633 /* NaN has maximum exponent */
1634 if ((x[NE - 1] & 0x7fff) != 0x7fff)
1636 /* ... and non-zero significand field. */
1637 for (i = 0; i < NE - 1; i++)
1647 /* Fill e-type number X with infinity pattern (IEEE)
1648 or largest possible number (non-IEEE). */
1652 register unsigned EMUSHORT *x;
1657 for (i = 0; i < NE - 1; i++)
1661 for (i = 0; i < NE - 1; i++)
1689 /* Output an e-type NaN.
1690 This generates Intel's quiet NaN pattern for extended real.
1691 The exponent is 7fff, the leading mantissa word is c000. */
1695 register unsigned EMUSHORT *x;
1700 for (i = 0; i < NE - 2; i++)
1703 *x = (sign << 15) | 0x7fff;
1706 /* Move in an e-type number A, converting it to exploded e-type B. */
1710 unsigned EMUSHORT *a, *b;
1712 register unsigned EMUSHORT *p, *q;
1716 p = a + (NE - 1); /* point to last word of external number */
1717 /* get the sign bit */
1722 /* get the exponent */
1724 *q++ &= 0x7fff; /* delete the sign bit */
1726 if ((*(q - 1) & 0x7fff) == 0x7fff)
1732 for (i = 3; i < NI; i++)
1738 for (i = 2; i < NI; i++)
1744 /* clear high guard word */
1746 /* move in the significand */
1747 for (i = 0; i < NE - 1; i++)
1749 /* clear low guard word */
1753 /* Move out exploded e-type number A, converting it to e type B. */
1757 unsigned EMUSHORT *a, *b;
1759 register unsigned EMUSHORT *p, *q;
1760 unsigned EMUSHORT i;
1764 q = b + (NE - 1); /* point to output exponent */
1765 /* combine sign and exponent */
1768 *q-- = *p++ | 0x8000;
1772 if (*(p - 1) == 0x7fff)
1777 enan (b, eiisneg (a));
1785 /* skip over guard word */
1787 /* move the significand */
1788 for (j = 0; j < NE - 1; j++)
1792 /* Clear out exploded e-type number XI. */
1796 register unsigned EMUSHORT *xi;
1800 for (i = 0; i < NI; i++)
1804 /* Clear out exploded e-type XI, but don't touch the sign. */
1808 register unsigned EMUSHORT *xi;
1813 for (i = 0; i < NI - 1; i++)
1817 /* Move exploded e-type number from A to B. */
1821 register unsigned EMUSHORT *a, *b;
1825 for (i = 0; i < NI - 1; i++)
1827 /* clear low guard word */
1831 /* Generate exploded e-type NaN.
1832 The explicit pattern for this is maximum exponent and
1833 top two significant bits set. */
1837 unsigned EMUSHORT x[];
1845 /* Return nonzero if exploded e-type X is a NaN. */
1849 unsigned EMUSHORT x[];
1853 if ((x[E] & 0x7fff) == 0x7fff)
1855 for (i = M + 1; i < NI; i++)
1864 /* Return nonzero if sign of exploded e-type X is nonzero. */
1868 unsigned EMUSHORT x[];
1875 /* Fill exploded e-type X with infinity pattern.
1876 This has maximum exponent and significand all zeros. */
1880 unsigned EMUSHORT x[];
1888 /* Return nonzero if exploded e-type X is infinite. */
1892 unsigned EMUSHORT x[];
1899 if ((x[E] & 0x7fff) == 0x7fff)
1905 /* Compare significands of numbers in internal exploded e-type format.
1906 Guard words are included in the comparison.
1914 register unsigned EMUSHORT *a, *b;
1918 a += M; /* skip up to significand area */
1920 for (i = M; i < NI; i++)
1928 if (*(--a) > *(--b))
1934 /* Shift significand of exploded e-type X down by 1 bit. */
1938 register unsigned EMUSHORT *x;
1940 register unsigned EMUSHORT bits;
1943 x += M; /* point to significand area */
1946 for (i = M; i < NI; i++)
1958 /* Shift significand of exploded e-type X up by 1 bit. */
1962 register unsigned EMUSHORT *x;
1964 register unsigned EMUSHORT bits;
1970 for (i = M; i < NI; i++)
1983 /* Shift significand of exploded e-type X down by 8 bits. */
1987 register unsigned EMUSHORT *x;
1989 register unsigned EMUSHORT newbyt, oldbyt;
1994 for (i = M; i < NI; i++)
2004 /* Shift significand of exploded e-type X up by 8 bits. */
2008 register unsigned EMUSHORT *x;
2011 register unsigned EMUSHORT newbyt, oldbyt;
2016 for (i = M; i < NI; i++)
2026 /* Shift significand of exploded e-type X up by 16 bits. */
2030 register unsigned EMUSHORT *x;
2033 register unsigned EMUSHORT *p;
2038 for (i = M; i < NI - 1; i++)
2044 /* Shift significand of exploded e-type X down by 16 bits. */
2048 register unsigned EMUSHORT *x;
2051 register unsigned EMUSHORT *p;
2056 for (i = M; i < NI - 1; i++)
2062 /* Add significands of exploded e-type X and Y. X + Y replaces Y. */
2066 unsigned EMUSHORT *x, *y;
2068 register unsigned EMULONG a;
2075 for (i = M; i < NI; i++)
2077 a = (unsigned EMULONG) (*x) + (unsigned EMULONG) (*y) + carry;
2082 *y = (unsigned EMUSHORT) a;
2088 /* Subtract significands of exploded e-type X and Y. Y - X replaces Y. */
2092 unsigned EMUSHORT *x, *y;
2101 for (i = M; i < NI; i++)
2103 a = (unsigned EMULONG) (*y) - (unsigned EMULONG) (*x) - carry;
2108 *y = (unsigned EMUSHORT) a;
2115 static unsigned EMUSHORT equot[NI];
2119 /* Radix 2 shift-and-add versions of multiply and divide */
2122 /* Divide significands */
2126 unsigned EMUSHORT den[], num[];
2129 register unsigned EMUSHORT *p, *q;
2130 unsigned EMUSHORT j;
2136 for (i = M; i < NI; i++)
2141 /* Use faster compare and subtraction if denominator has only 15 bits of
2147 for (i = M + 3; i < NI; i++)
2152 if ((den[M + 1] & 1) != 0)
2160 for (i = 0; i < NBITS + 2; i++)
2178 /* The number of quotient bits to calculate is NBITS + 1 scaling guard
2179 bit + 1 roundoff bit. */
2184 for (i = 0; i < NBITS + 2; i++)
2186 if (ecmpm (den, num) <= 0)
2189 j = 1; /* quotient bit = 1 */
2203 /* test for nonzero remainder after roundoff bit */
2206 for (i = M; i < NI; i++)
2214 for (i = 0; i < NI; i++)
2220 /* Multiply significands */
2224 unsigned EMUSHORT a[], b[];
2226 unsigned EMUSHORT *p, *q;
2231 for (i = M; i < NI; i++)
2236 while (*p == 0) /* significand is not supposed to be zero */
2241 if ((*p & 0xff) == 0)
2249 for (i = 0; i < k; i++)
2253 /* remember if there were any nonzero bits shifted out */
2260 for (i = 0; i < NI; i++)
2263 /* return flag for lost nonzero bits */
2269 /* Radix 65536 versions of multiply and divide. */
2271 /* Multiply significand of e-type number B
2272 by 16-bit quantity A, return e-type result to C. */
2277 unsigned EMUSHORT b[], c[];
2279 register unsigned EMUSHORT *pp;
2280 register unsigned EMULONG carry;
2281 unsigned EMUSHORT *ps;
2282 unsigned EMUSHORT p[NI];
2283 unsigned EMULONG aa, m;
2292 for (i=M+1; i<NI; i++)
2302 m = (unsigned EMULONG) aa * *ps--;
2303 carry = (m & 0xffff) + *pp;
2304 *pp-- = (unsigned EMUSHORT)carry;
2305 carry = (carry >> 16) + (m >> 16) + *pp;
2306 *pp = (unsigned EMUSHORT)carry;
2307 *(pp-1) = carry >> 16;
2310 for (i=M; i<NI; i++)
2314 /* Divide significands of exploded e-types NUM / DEN. Neither the
2315 numerator NUM nor the denominator DEN is permitted to have its high guard
2320 unsigned EMUSHORT den[], num[];
2323 register unsigned EMUSHORT *p;
2324 unsigned EMULONG tnum;
2325 unsigned EMUSHORT j, tdenm, tquot;
2326 unsigned EMUSHORT tprod[NI+1];
2332 for (i=M; i<NI; i++)
2338 for (i=M; i<NI; i++)
2340 /* Find trial quotient digit (the radix is 65536). */
2341 tnum = (((unsigned EMULONG) num[M]) << 16) + num[M+1];
2343 /* Do not execute the divide instruction if it will overflow. */
2344 if ((tdenm * (unsigned long)0xffff) < tnum)
2347 tquot = tnum / tdenm;
2348 /* Multiply denominator by trial quotient digit. */
2349 m16m ((unsigned int)tquot, den, tprod);
2350 /* The quotient digit may have been overestimated. */
2351 if (ecmpm (tprod, num) > 0)
2355 if (ecmpm (tprod, num) > 0)
2365 /* test for nonzero remainder after roundoff bit */
2368 for (i=M; i<NI; i++)
2375 for (i=0; i<NI; i++)
2381 /* Multiply significands of exploded e-type A and B, result in B. */
2385 unsigned EMUSHORT a[], b[];
2387 unsigned EMUSHORT *p, *q;
2388 unsigned EMUSHORT pprod[NI];
2389 unsigned EMUSHORT j;
2394 for (i=M; i<NI; i++)
2400 for (i=M+1; i<NI; i++)
2408 m16m ((unsigned int) *p--, b, pprod);
2409 eaddm(pprod, equot);
2415 for (i=0; i<NI; i++)
2418 /* return flag for lost nonzero bits */
2424 /* Normalize and round off.
2426 The internal format number to be rounded is S.
2427 Input LOST is 0 if the value is exact. This is the so-called sticky bit.
2429 Input SUBFLG indicates whether the number was obtained
2430 by a subtraction operation. In that case if LOST is nonzero
2431 then the number is slightly smaller than indicated.
2433 Input EXP is the biased exponent, which may be negative.
2434 the exponent field of S is ignored but is replaced by
2435 EXP as adjusted by normalization and rounding.
2437 Input RCNTRL is the rounding control. If it is nonzero, the
2438 returned value will be rounded to RNDPRC bits.
2440 For future reference: In order for emdnorm to round off denormal
2441 significands at the right point, the input exponent must be
2442 adjusted to be the actual value it would have after conversion to
2443 the final floating point type. This adjustment has been
2444 implemented for all type conversions (etoe53, etc.) and decimal
2445 conversions, but not for the arithmetic functions (eadd, etc.).
2446 Data types having standard 15-bit exponents are not affected by
2447 this, but SFmode and DFmode are affected. For example, ediv with
2448 rndprc = 24 will not round correctly to 24-bit precision if the
2449 result is denormal. */
2451 static int rlast = -1;
2453 static unsigned EMUSHORT rmsk = 0;
2454 static unsigned EMUSHORT rmbit = 0;
2455 static unsigned EMUSHORT rebit = 0;
2457 static unsigned EMUSHORT rbit[NI];
2460 emdnorm (s, lost, subflg, exp, rcntrl)
2461 unsigned EMUSHORT s[];
2468 unsigned EMUSHORT r;
2473 /* a blank significand could mean either zero or infinity. */
2486 if ((j > NBITS) && (exp < 32767))
2494 if (exp > (EMULONG) (-NBITS - 1))
2507 /* Round off, unless told not to by rcntrl. */
2510 /* Set up rounding parameters if the control register changed. */
2511 if (rndprc != rlast)
2518 rw = NI - 1; /* low guard word */
2541 /* For DEC or IBM arithmetic */
2558 /* For C4x arithmetic */
2579 /* Shift down 1 temporarily if the data structure has an implied
2580 most significant bit and the number is denormal.
2581 Intel long double denormals also lose one bit of precision. */
2582 if ((exp <= 0) && (rndprc != NBITS)
2583 && ((rndprc != 64) || ((rndprc == 64) && ! REAL_WORDS_BIG_ENDIAN)))
2585 lost |= s[NI - 1] & 1;
2588 /* Clear out all bits below the rounding bit,
2589 remembering in r if any were nonzero. */
2603 if ((r & rmbit) != 0)
2609 { /* round to even */
2610 if ((s[re] & rebit) == 0)
2623 /* Undo the temporary shift for denormal values. */
2624 if ((exp <= 0) && (rndprc != NBITS)
2625 && ((rndprc != 64) || ((rndprc == 64) && ! REAL_WORDS_BIG_ENDIAN)))
2630 { /* overflow on roundoff */
2643 for (i = 2; i < NI - 1; i++)
2646 warning ("floating point overflow");
2650 for (i = M + 1; i < NI - 1; i++)
2653 if ((rndprc < 64) || (rndprc == 113))
2668 s[1] = (unsigned EMUSHORT) exp;
2671 /* Subtract. C = B - A, all e type numbers. */
2673 static int subflg = 0;
2677 unsigned EMUSHORT *a, *b, *c;
2691 /* Infinity minus infinity is a NaN.
2692 Test for subtracting infinities of the same sign. */
2693 if (eisinf (a) && eisinf (b)
2694 && ((eisneg (a) ^ eisneg (b)) == 0))
2696 mtherr ("esub", INVALID);
2705 /* Add. C = A + B, all e type. */
2709 unsigned EMUSHORT *a, *b, *c;
2713 /* NaN plus anything is a NaN. */
2724 /* Infinity minus infinity is a NaN.
2725 Test for adding infinities of opposite signs. */
2726 if (eisinf (a) && eisinf (b)
2727 && ((eisneg (a) ^ eisneg (b)) != 0))
2729 mtherr ("esub", INVALID);
2738 /* Arithmetic common to both addition and subtraction. */
2742 unsigned EMUSHORT *a, *b, *c;
2744 unsigned EMUSHORT ai[NI], bi[NI], ci[NI];
2746 EMULONG lt, lta, ltb;
2767 /* compare exponents */
2772 { /* put the larger number in bi */
2782 if (lt < (EMULONG) (-NBITS - 1))
2783 goto done; /* answer same as larger addend */
2785 lost = eshift (ai, k); /* shift the smaller number down */
2789 /* exponents were the same, so must compare significands */
2792 { /* the numbers are identical in magnitude */
2793 /* if different signs, result is zero */
2799 /* if same sign, result is double */
2800 /* double denormalized tiny number */
2801 if ((bi[E] == 0) && ((bi[3] & 0x8000) == 0))
2806 /* add 1 to exponent unless both are zero! */
2807 for (j = 1; j < NI - 1; j++)
2823 bi[E] = (unsigned EMUSHORT) ltb;
2827 { /* put the larger number in bi */
2843 emdnorm (bi, lost, subflg, ltb, 64);
2849 /* Divide: C = B/A, all e type. */
2853 unsigned EMUSHORT *a, *b, *c;
2855 unsigned EMUSHORT ai[NI], bi[NI];
2857 EMULONG lt, lta, ltb;
2859 /* IEEE says if result is not a NaN, the sign is "-" if and only if
2860 operands have opposite signs -- but flush -0 to 0 later if not IEEE. */
2861 sign = eisneg(a) ^ eisneg(b);
2864 /* Return any NaN input. */
2875 /* Zero over zero, or infinity over infinity, is a NaN. */
2876 if (((ecmp (a, ezero) == 0) && (ecmp (b, ezero) == 0))
2877 || (eisinf (a) && eisinf (b)))
2879 mtherr ("ediv", INVALID);
2884 /* Infinity over anything else is infinity. */
2891 /* Anything else over infinity is zero. */
2903 { /* See if numerator is zero. */
2904 for (i = 1; i < NI - 1; i++)
2908 ltb -= enormlz (bi);
2918 { /* possible divide by zero */
2919 for (i = 1; i < NI - 1; i++)
2923 lta -= enormlz (ai);
2927 /* Divide by zero is not an invalid operation.
2928 It is a divide-by-zero operation! */
2930 mtherr ("ediv", SING);
2936 /* calculate exponent */
2937 lt = ltb - lta + EXONE;
2938 emdnorm (bi, i, 0, lt, 64);
2945 && (ecmp (c, ezero) != 0)
2948 *(c+(NE-1)) |= 0x8000;
2950 *(c+(NE-1)) &= ~0x8000;
2953 /* Multiply e-types A and B, return e-type product C. */
2957 unsigned EMUSHORT *a, *b, *c;
2959 unsigned EMUSHORT ai[NI], bi[NI];
2961 EMULONG lt, lta, ltb;
2963 /* IEEE says if result is not a NaN, the sign is "-" if and only if
2964 operands have opposite signs -- but flush -0 to 0 later if not IEEE. */
2965 sign = eisneg(a) ^ eisneg(b);
2968 /* NaN times anything is the same NaN. */
2979 /* Zero times infinity is a NaN. */
2980 if ((eisinf (a) && (ecmp (b, ezero) == 0))
2981 || (eisinf (b) && (ecmp (a, ezero) == 0)))
2983 mtherr ("emul", INVALID);
2988 /* Infinity times anything else is infinity. */
2990 if (eisinf (a) || eisinf (b))
3002 for (i = 1; i < NI - 1; i++)
3006 lta -= enormlz (ai);
3017 for (i = 1; i < NI - 1; i++)
3021 ltb -= enormlz (bi);
3030 /* Multiply significands */
3032 /* calculate exponent */
3033 lt = lta + ltb - (EXONE - 1);
3034 emdnorm (bi, j, 0, lt, 64);
3041 && (ecmp (c, ezero) != 0)
3044 *(c+(NE-1)) |= 0x8000;
3046 *(c+(NE-1)) &= ~0x8000;
3049 /* Convert double precision PE to e-type Y. */
3053 unsigned EMUSHORT *pe, *y;
3062 ibmtoe (pe, y, DFmode);
3067 c4xtoe (pe, y, HFmode);
3070 register unsigned EMUSHORT r;
3071 register unsigned EMUSHORT *e, *p;
3072 unsigned EMUSHORT yy[NI];
3076 denorm = 0; /* flag if denormalized number */
3078 if (! REAL_WORDS_BIG_ENDIAN)
3084 yy[M] = (r & 0x0f) | 0x10;
3085 r &= ~0x800f; /* strip sign and 4 significand bits */
3090 if (! REAL_WORDS_BIG_ENDIAN)
3092 if (((pe[3] & 0xf) != 0) || (pe[2] != 0)
3093 || (pe[1] != 0) || (pe[0] != 0))
3095 enan (y, yy[0] != 0);
3101 if (((pe[0] & 0xf) != 0) || (pe[1] != 0)
3102 || (pe[2] != 0) || (pe[3] != 0))
3104 enan (y, yy[0] != 0);
3115 #endif /* INFINITY */
3117 /* If zero exponent, then the significand is denormalized.
3118 So take back the understood high significand bit. */
3129 if (! REAL_WORDS_BIG_ENDIAN)
3146 /* If zero exponent, then normalize the significand. */
3147 if ((k = enormlz (yy)) > NBITS)
3150 yy[E] -= (unsigned EMUSHORT) (k - 1);
3153 #endif /* not C4X */
3154 #endif /* not IBM */
3155 #endif /* not DEC */
3158 /* Convert double extended precision float PE to e type Y. */
3162 unsigned EMUSHORT *pe, *y;
3164 unsigned EMUSHORT yy[NI];
3165 unsigned EMUSHORT *e, *p, *q;
3170 for (i = 0; i < NE - 5; i++)
3172 /* This precision is not ordinarily supported on DEC or IBM. */
3174 for (i = 0; i < 5; i++)
3178 p = &yy[0] + (NE - 1);
3181 for (i = 0; i < 5; i++)
3185 if (! REAL_WORDS_BIG_ENDIAN)
3187 for (i = 0; i < 5; i++)
3190 /* For denormal long double Intel format, shift significand up one
3191 -- but only if the top significand bit is zero. A top bit of 1
3192 is "pseudodenormal" when the exponent is zero. */
3193 if((yy[NE-1] & 0x7fff) == 0 && (yy[NE-2] & 0x8000) == 0)
3195 unsigned EMUSHORT temp[NI];
3205 p = &yy[0] + (NE - 1);
3206 #ifdef ARM_EXTENDED_IEEE_FORMAT
3207 /* For ARMs, the exponent is in the lowest 15 bits of the word. */
3208 *p-- = (e[0] & 0x8000) | (e[1] & 0x7ffff);
3214 for (i = 0; i < 4; i++)
3219 /* Point to the exponent field and check max exponent cases. */
3221 if ((*p & 0x7fff) == 0x7fff)
3224 if (! REAL_WORDS_BIG_ENDIAN)
3226 for (i = 0; i < 4; i++)
3228 if ((i != 3 && pe[i] != 0)
3229 /* Anything but 0x8000 here, including 0, is a NaN. */
3230 || (i == 3 && pe[i] != 0x8000))
3232 enan (y, (*p & 0x8000) != 0);
3239 #ifdef ARM_EXTENDED_IEEE_FORMAT
3240 for (i = 2; i <= 5; i++)
3244 enan (y, (*p & 0x8000) != 0);
3249 /* In Motorola extended precision format, the most significant
3250 bit of an infinity mantissa could be either 1 or 0. It is
3251 the lower order bits that tell whether the value is a NaN. */
3252 if ((pe[2] & 0x7fff) != 0)
3255 for (i = 3; i <= 5; i++)
3260 enan (y, (*p & 0x8000) != 0);
3264 #endif /* not ARM */
3273 #endif /* INFINITY */
3276 for (i = 0; i < NE; i++)
3280 /* Convert 128-bit long double precision float PE to e type Y. */
3284 unsigned EMUSHORT *pe, *y;
3286 register unsigned EMUSHORT r;
3287 unsigned EMUSHORT *e, *p;
3288 unsigned EMUSHORT yy[NI];
3295 if (! REAL_WORDS_BIG_ENDIAN)
3307 if (! REAL_WORDS_BIG_ENDIAN)
3309 for (i = 0; i < 7; i++)
3313 enan (y, yy[0] != 0);
3320 for (i = 1; i < 8; i++)
3324 enan (y, yy[0] != 0);
3336 #endif /* INFINITY */
3340 if (! REAL_WORDS_BIG_ENDIAN)
3342 for (i = 0; i < 7; i++)
3348 for (i = 0; i < 7; i++)
3352 /* If denormal, remove the implied bit; else shift down 1. */
3365 /* Convert single precision float PE to e type Y. */
3369 unsigned EMUSHORT *pe, *y;
3373 ibmtoe (pe, y, SFmode);
3379 c4xtoe (pe, y, QFmode);
3383 register unsigned EMUSHORT r;
3384 register unsigned EMUSHORT *e, *p;
3385 unsigned EMUSHORT yy[NI];
3389 denorm = 0; /* flag if denormalized number */
3392 if (! REAL_WORDS_BIG_ENDIAN)
3402 yy[M] = (r & 0x7f) | 0200;
3403 r &= ~0x807f; /* strip sign and 7 significand bits */
3408 if (REAL_WORDS_BIG_ENDIAN)
3410 if (((pe[0] & 0x7f) != 0) || (pe[1] != 0))
3412 enan (y, yy[0] != 0);
3418 if (((pe[1] & 0x7f) != 0) || (pe[0] != 0))
3420 enan (y, yy[0] != 0);
3431 #endif /* INFINITY */
3433 /* If zero exponent, then the significand is denormalized.
3434 So take back the understood high significand bit. */
3447 if (! REAL_WORDS_BIG_ENDIAN)
3457 { /* if zero exponent, then normalize the significand */
3458 if ((k = enormlz (yy)) > NBITS)
3461 yy[E] -= (unsigned EMUSHORT) (k - 1);
3464 #endif /* not C4X */
3465 #endif /* not IBM */
3468 /* Convert e-type X to IEEE 128-bit long double format E. */
3472 unsigned EMUSHORT *x, *e;
3474 unsigned EMUSHORT xi[NI];
3481 make_nan (e, eisneg (x), TFmode);
3486 exp = (EMULONG) xi[E];
3491 /* round off to nearest or even */
3494 emdnorm (xi, 0, 0, exp, 64);
3500 /* Convert exploded e-type X, that has already been rounded to
3501 113-bit precision, to IEEE 128-bit long double format Y. */
3505 unsigned EMUSHORT *a, *b;
3507 register unsigned EMUSHORT *p, *q;
3508 unsigned EMUSHORT i;
3513 make_nan (b, eiisneg (a), TFmode);
3518 if (REAL_WORDS_BIG_ENDIAN)
3521 q = b + 7; /* point to output exponent */
3523 /* If not denormal, delete the implied bit. */
3528 /* combine sign and exponent */
3530 if (REAL_WORDS_BIG_ENDIAN)
3533 *q++ = *p++ | 0x8000;
3540 *q-- = *p++ | 0x8000;
3544 /* skip over guard word */
3546 /* move the significand */
3547 if (REAL_WORDS_BIG_ENDIAN)
3549 for (i = 0; i < 7; i++)
3554 for (i = 0; i < 7; i++)
3559 /* Convert e-type X to IEEE double extended format E. */
3563 unsigned EMUSHORT *x, *e;
3565 unsigned EMUSHORT xi[NI];
3572 make_nan (e, eisneg (x), XFmode);
3577 /* adjust exponent for offset */
3578 exp = (EMULONG) xi[E];
3583 /* round off to nearest or even */
3586 emdnorm (xi, 0, 0, exp, 64);
3592 /* Convert exploded e-type X, that has already been rounded to
3593 64-bit precision, to IEEE double extended format Y. */
3597 unsigned EMUSHORT *a, *b;
3599 register unsigned EMUSHORT *p, *q;
3600 unsigned EMUSHORT i;
3605 make_nan (b, eiisneg (a), XFmode);
3609 /* Shift denormal long double Intel format significand down one bit. */
3610 if ((a[E] == 0) && ! REAL_WORDS_BIG_ENDIAN)
3620 if (REAL_WORDS_BIG_ENDIAN)
3624 q = b + 4; /* point to output exponent */
3625 #if LONG_DOUBLE_TYPE_SIZE == 96
3626 /* Clear the last two bytes of 12-byte Intel format */
3632 /* combine sign and exponent */
3636 *q++ = *p++ | 0x8000;
3643 *q-- = *p++ | 0x8000;
3648 if (REAL_WORDS_BIG_ENDIAN)
3650 #ifdef ARM_EXTENDED_IEEE_FORMAT
3651 /* The exponent is in the lowest 15 bits of the first word. */
3652 *q++ = i ? 0x8000 : 0;
3656 *q++ = *p++ | 0x8000;
3665 *q-- = *p++ | 0x8000;
3670 /* skip over guard word */
3672 /* move the significand */
3674 for (i = 0; i < 4; i++)
3678 for (i = 0; i < 4; i++)
3682 if (REAL_WORDS_BIG_ENDIAN)
3684 for (i = 0; i < 4; i++)
3692 /* Intel long double infinity significand. */
3700 for (i = 0; i < 4; i++)
3706 /* e type to double precision. */
3709 /* Convert e-type X to DEC-format double E. */
3713 unsigned EMUSHORT *x, *e;
3715 etodec (x, e); /* see etodec.c */
3718 /* Convert exploded e-type X, that has already been rounded to
3719 56-bit double precision, to DEC double Y. */
3723 unsigned EMUSHORT *x, *y;
3730 /* Convert e-type X to IBM 370-format double E. */
3734 unsigned EMUSHORT *x, *e;
3736 etoibm (x, e, DFmode);
3739 /* Convert exploded e-type X, that has already been rounded to
3740 56-bit precision, to IBM 370 double Y. */
3744 unsigned EMUSHORT *x, *y;
3746 toibm (x, y, DFmode);
3749 #else /* it's neither DEC nor IBM */
3751 /* Convert e-type X to C4X-format long double E. */
3755 unsigned EMUSHORT *x, *e;
3757 etoc4x (x, e, HFmode);
3760 /* Convert exploded e-type X, that has already been rounded to
3761 56-bit precision, to IBM 370 double Y. */
3765 unsigned EMUSHORT *x, *y;
3767 toc4x (x, y, HFmode);
3770 #else /* it's neither DEC nor IBM nor C4X */
3772 /* Convert e-type X to IEEE double E. */
3776 unsigned EMUSHORT *x, *e;
3778 unsigned EMUSHORT xi[NI];
3785 make_nan (e, eisneg (x), DFmode);
3790 /* adjust exponent for offsets */
3791 exp = (EMULONG) xi[E] - (EXONE - 0x3ff);
3796 /* round off to nearest or even */
3799 emdnorm (xi, 0, 0, exp, 64);
3805 /* Convert exploded e-type X, that has already been rounded to
3806 53-bit precision, to IEEE double Y. */
3810 unsigned EMUSHORT *x, *y;
3812 unsigned EMUSHORT i;
3813 unsigned EMUSHORT *p;
3818 make_nan (y, eiisneg (x), DFmode);
3824 if (! REAL_WORDS_BIG_ENDIAN)
3827 *y = 0; /* output high order */
3829 *y = 0x8000; /* output sign bit */
3832 if (i >= (unsigned int) 2047)
3834 /* Saturate at largest number less than infinity. */
3837 if (! REAL_WORDS_BIG_ENDIAN)
3851 *y |= (unsigned EMUSHORT) 0x7fef;
3852 if (! REAL_WORDS_BIG_ENDIAN)
3877 i |= *p++ & (unsigned EMUSHORT) 0x0f; /* *p = xi[M] */
3878 *y |= (unsigned EMUSHORT) i; /* high order output already has sign bit set */
3879 if (! REAL_WORDS_BIG_ENDIAN)
3894 #endif /* not C4X */
3895 #endif /* not IBM */
3896 #endif /* not DEC */
3900 /* e type to single precision. */
3903 /* Convert e-type X to IBM 370 float E. */
3907 unsigned EMUSHORT *x, *e;
3909 etoibm (x, e, SFmode);
3912 /* Convert exploded e-type X, that has already been rounded to
3913 float precision, to IBM 370 float Y. */
3917 unsigned EMUSHORT *x, *y;
3919 toibm (x, y, SFmode);
3925 /* Convert e-type X to C4X float E. */
3929 unsigned EMUSHORT *x, *e;
3931 etoc4x (x, e, QFmode);
3934 /* Convert exploded e-type X, that has already been rounded to
3935 float precision, to IBM 370 float Y. */
3939 unsigned EMUSHORT *x, *y;
3941 toc4x (x, y, QFmode);
3946 /* Convert e-type X to IEEE float E. DEC float is the same as IEEE float. */
3950 unsigned EMUSHORT *x, *e;
3953 unsigned EMUSHORT xi[NI];
3959 make_nan (e, eisneg (x), SFmode);
3964 /* adjust exponent for offsets */
3965 exp = (EMULONG) xi[E] - (EXONE - 0177);
3970 /* round off to nearest or even */
3973 emdnorm (xi, 0, 0, exp, 64);
3979 /* Convert exploded e-type X, that has already been rounded to
3980 float precision, to IEEE float Y. */
3984 unsigned EMUSHORT *x, *y;
3986 unsigned EMUSHORT i;
3987 unsigned EMUSHORT *p;
3992 make_nan (y, eiisneg (x), SFmode);
3998 if (! REAL_WORDS_BIG_ENDIAN)
4004 *y = 0; /* output high order */
4006 *y = 0x8000; /* output sign bit */
4009 /* Handle overflow cases. */
4013 *y |= (unsigned EMUSHORT) 0x7f80;
4018 if (! REAL_WORDS_BIG_ENDIAN)
4026 #else /* no INFINITY */
4027 *y |= (unsigned EMUSHORT) 0x7f7f;
4032 if (! REAL_WORDS_BIG_ENDIAN)
4043 #endif /* no INFINITY */
4055 i |= *p++ & (unsigned EMUSHORT) 0x7f; /* *p = xi[M] */
4056 /* High order output already has sign bit set. */
4062 if (! REAL_WORDS_BIG_ENDIAN)
4071 #endif /* not C4X */
4072 #endif /* not IBM */
4074 /* Compare two e type numbers.
4078 -2 if either a or b is a NaN. */
4082 unsigned EMUSHORT *a, *b;
4084 unsigned EMUSHORT ai[NI], bi[NI];
4085 register unsigned EMUSHORT *p, *q;
4090 if (eisnan (a) || eisnan (b))
4099 { /* the signs are different */
4101 for (i = 1; i < NI - 1; i++)
4115 /* both are the same sign */
4130 return (0); /* equality */
4134 if (*(--p) > *(--q))
4135 return (msign); /* p is bigger */
4137 return (-msign); /* p is littler */
4141 /* Find e-type nearest integer to X, as floor (X + 0.5). */
4145 unsigned EMUSHORT *x, *y;
4152 /* Convert HOST_WIDE_INT LP to e type Y. */
4157 unsigned EMUSHORT *y;
4159 unsigned EMUSHORT yi[NI];
4160 unsigned HOST_WIDE_INT ll;
4166 /* make it positive */
4167 ll = (unsigned HOST_WIDE_INT) (-(*lp));
4168 yi[0] = 0xffff; /* put correct sign in the e type number */
4172 ll = (unsigned HOST_WIDE_INT) (*lp);
4174 /* move the long integer to yi significand area */
4175 #if HOST_BITS_PER_WIDE_INT == 64
4176 yi[M] = (unsigned EMUSHORT) (ll >> 48);
4177 yi[M + 1] = (unsigned EMUSHORT) (ll >> 32);
4178 yi[M + 2] = (unsigned EMUSHORT) (ll >> 16);
4179 yi[M + 3] = (unsigned EMUSHORT) ll;
4180 yi[E] = EXONE + 47; /* exponent if normalize shift count were 0 */
4182 yi[M] = (unsigned EMUSHORT) (ll >> 16);
4183 yi[M + 1] = (unsigned EMUSHORT) ll;
4184 yi[E] = EXONE + 15; /* exponent if normalize shift count were 0 */
4187 if ((k = enormlz (yi)) > NBITS)/* normalize the significand */
4188 ecleaz (yi); /* it was zero */
4190 yi[E] -= (unsigned EMUSHORT) k;/* subtract shift count from exponent */
4191 emovo (yi, y); /* output the answer */
4194 /* Convert unsigned HOST_WIDE_INT LP to e type Y. */
4198 unsigned HOST_WIDE_INT *lp;
4199 unsigned EMUSHORT *y;
4201 unsigned EMUSHORT yi[NI];
4202 unsigned HOST_WIDE_INT ll;
4208 /* move the long integer to ayi significand area */
4209 #if HOST_BITS_PER_WIDE_INT == 64
4210 yi[M] = (unsigned EMUSHORT) (ll >> 48);
4211 yi[M + 1] = (unsigned EMUSHORT) (ll >> 32);
4212 yi[M + 2] = (unsigned EMUSHORT) (ll >> 16);
4213 yi[M + 3] = (unsigned EMUSHORT) ll;
4214 yi[E] = EXONE + 47; /* exponent if normalize shift count were 0 */
4216 yi[M] = (unsigned EMUSHORT) (ll >> 16);
4217 yi[M + 1] = (unsigned EMUSHORT) ll;
4218 yi[E] = EXONE + 15; /* exponent if normalize shift count were 0 */
4221 if ((k = enormlz (yi)) > NBITS)/* normalize the significand */
4222 ecleaz (yi); /* it was zero */
4224 yi[E] -= (unsigned EMUSHORT) k; /* subtract shift count from exponent */
4225 emovo (yi, y); /* output the answer */
4229 /* Find signed HOST_WIDE_INT integer I and floating point fractional
4230 part FRAC of e-type (packed internal format) floating point input X.
4231 The integer output I has the sign of the input, except that
4232 positive overflow is permitted if FIXUNS_TRUNC_LIKE_FIX_TRUNC.
4233 The output e-type fraction FRAC is the positive fractional
4238 unsigned EMUSHORT *x;
4240 unsigned EMUSHORT *frac;
4242 unsigned EMUSHORT xi[NI];
4244 unsigned HOST_WIDE_INT ll;
4247 k = (int) xi[E] - (EXONE - 1);
4250 /* if exponent <= 0, integer = 0 and real output is fraction */
4255 if (k > (HOST_BITS_PER_WIDE_INT - 1))
4257 /* long integer overflow: output large integer
4258 and correct fraction */
4260 *i = ((unsigned HOST_WIDE_INT) 1) << (HOST_BITS_PER_WIDE_INT - 1);
4263 #ifdef FIXUNS_TRUNC_LIKE_FIX_TRUNC
4264 /* In this case, let it overflow and convert as if unsigned. */
4265 euifrac (x, &ll, frac);
4266 *i = (HOST_WIDE_INT) ll;
4269 /* In other cases, return the largest positive integer. */
4270 *i = (((unsigned HOST_WIDE_INT) 1) << (HOST_BITS_PER_WIDE_INT - 1)) - 1;
4275 warning ("overflow on truncation to integer");
4279 /* Shift more than 16 bits: first shift up k-16 mod 16,
4280 then shift up by 16's. */
4281 j = k - ((k >> 4) << 4);
4288 ll = (ll << 16) | xi[M];
4290 while ((k -= 16) > 0);
4297 /* shift not more than 16 bits */
4299 *i = (HOST_WIDE_INT) xi[M] & 0xffff;
4306 if ((k = enormlz (xi)) > NBITS)
4309 xi[E] -= (unsigned EMUSHORT) k;
4315 /* Find unsigned HOST_WIDE_INT integer I and floating point fractional part
4316 FRAC of e-type X. A negative input yields integer output = 0 but
4317 correct fraction. */
4320 euifrac (x, i, frac)
4321 unsigned EMUSHORT *x;
4322 unsigned HOST_WIDE_INT *i;
4323 unsigned EMUSHORT *frac;
4325 unsigned HOST_WIDE_INT ll;
4326 unsigned EMUSHORT xi[NI];
4330 k = (int) xi[E] - (EXONE - 1);
4333 /* if exponent <= 0, integer = 0 and argument is fraction */
4338 if (k > HOST_BITS_PER_WIDE_INT)
4340 /* Long integer overflow: output large integer
4341 and correct fraction.
4342 Note, the BSD microvax compiler says that ~(0UL)
4343 is a syntax error. */
4347 warning ("overflow on truncation to unsigned integer");
4351 /* Shift more than 16 bits: first shift up k-16 mod 16,
4352 then shift up by 16's. */
4353 j = k - ((k >> 4) << 4);
4360 ll = (ll << 16) | xi[M];
4362 while ((k -= 16) > 0);
4367 /* shift not more than 16 bits */
4369 *i = (HOST_WIDE_INT) xi[M] & 0xffff;
4372 if (xi[0]) /* A negative value yields unsigned integer 0. */
4378 if ((k = enormlz (xi)) > NBITS)
4381 xi[E] -= (unsigned EMUSHORT) k;
4386 /* Shift the significand of exploded e-type X up or down by SC bits. */
4390 unsigned EMUSHORT *x;
4393 unsigned EMUSHORT lost;
4394 unsigned EMUSHORT *p;
4407 lost |= *p; /* remember lost bits */
4448 return ((int) lost);
4451 /* Shift normalize the significand area of exploded e-type X.
4452 Return the shift count (up = positive). */
4456 unsigned EMUSHORT x[];
4458 register unsigned EMUSHORT *p;
4467 return (0); /* already normalized */
4473 /* With guard word, there are NBITS+16 bits available.
4474 Return true if all are zero. */
4478 /* see if high byte is zero */
4479 while ((*p & 0xff00) == 0)
4484 /* now shift 1 bit at a time */
4485 while ((*p & 0x8000) == 0)
4491 mtherr ("enormlz", UNDERFLOW);
4497 /* Normalize by shifting down out of the high guard word
4498 of the significand */
4513 mtherr ("enormlz", OVERFLOW);
4520 /* Powers of ten used in decimal <-> binary conversions. */
4525 #if LONG_DOUBLE_TYPE_SIZE == 128
4526 static unsigned EMUSHORT etens[NTEN + 1][NE] =
4528 {0x6576, 0x4a92, 0x804a, 0x153f,
4529 0xc94c, 0x979a, 0x8a20, 0x5202, 0xc460, 0x7525,}, /* 10**4096 */
4530 {0x6a32, 0xce52, 0x329a, 0x28ce,
4531 0xa74d, 0x5de4, 0xc53d, 0x3b5d, 0x9e8b, 0x5a92,}, /* 10**2048 */
4532 {0x526c, 0x50ce, 0xf18b, 0x3d28,
4533 0x650d, 0x0c17, 0x8175, 0x7586, 0xc976, 0x4d48,},
4534 {0x9c66, 0x58f8, 0xbc50, 0x5c54,
4535 0xcc65, 0x91c6, 0xa60e, 0xa0ae, 0xe319, 0x46a3,},
4536 {0x851e, 0xeab7, 0x98fe, 0x901b,
4537 0xddbb, 0xde8d, 0x9df9, 0xebfb, 0xaa7e, 0x4351,},
4538 {0x0235, 0x0137, 0x36b1, 0x336c,
4539 0xc66f, 0x8cdf, 0x80e9, 0x47c9, 0x93ba, 0x41a8,},
4540 {0x50f8, 0x25fb, 0xc76b, 0x6b71,
4541 0x3cbf, 0xa6d5, 0xffcf, 0x1f49, 0xc278, 0x40d3,},
4542 {0x0000, 0x0000, 0x0000, 0x0000,
4543 0xf020, 0xb59d, 0x2b70, 0xada8, 0x9dc5, 0x4069,},
4544 {0x0000, 0x0000, 0x0000, 0x0000,
4545 0x0000, 0x0000, 0x0400, 0xc9bf, 0x8e1b, 0x4034,},
4546 {0x0000, 0x0000, 0x0000, 0x0000,
4547 0x0000, 0x0000, 0x0000, 0x2000, 0xbebc, 0x4019,},
4548 {0x0000, 0x0000, 0x0000, 0x0000,
4549 0x0000, 0x0000, 0x0000, 0x0000, 0x9c40, 0x400c,},
4550 {0x0000, 0x0000, 0x0000, 0x0000,
4551 0x0000, 0x0000, 0x0000, 0x0000, 0xc800, 0x4005,},
4552 {0x0000, 0x0000, 0x0000, 0x0000,
4553 0x0000, 0x0000, 0x0000, 0x0000, 0xa000, 0x4002,}, /* 10**1 */
4556 static unsigned EMUSHORT emtens[NTEN + 1][NE] =
4558 {0x2030, 0xcffc, 0xa1c3, 0x8123,
4559 0x2de3, 0x9fde, 0xd2ce, 0x04c8, 0xa6dd, 0x0ad8,}, /* 10**-4096 */
4560 {0x8264, 0xd2cb, 0xf2ea, 0x12d4,
4561 0x4925, 0x2de4, 0x3436, 0x534f, 0xceae, 0x256b,}, /* 10**-2048 */
4562 {0xf53f, 0xf698, 0x6bd3, 0x0158,
4563 0x87a6, 0xc0bd, 0xda57, 0x82a5, 0xa2a6, 0x32b5,},
4564 {0xe731, 0x04d4, 0xe3f2, 0xd332,
4565 0x7132, 0xd21c, 0xdb23, 0xee32, 0x9049, 0x395a,},
4566 {0xa23e, 0x5308, 0xfefb, 0x1155,
4567 0xfa91, 0x1939, 0x637a, 0x4325, 0xc031, 0x3cac,},
4568 {0xe26d, 0xdbde, 0xd05d, 0xb3f6,
4569 0xac7c, 0xe4a0, 0x64bc, 0x467c, 0xddd0, 0x3e55,},
4570 {0x2a20, 0x6224, 0x47b3, 0x98d7,
4571 0x3f23, 0xe9a5, 0xa539, 0xea27, 0xa87f, 0x3f2a,},
4572 {0x0b5b, 0x4af2, 0xa581, 0x18ed,
4573 0x67de, 0x94ba, 0x4539, 0x1ead, 0xcfb1, 0x3f94,},
4574 {0xbf71, 0xa9b3, 0x7989, 0xbe68,
4575 0x4c2e, 0xe15b, 0xc44d, 0x94be, 0xe695, 0x3fc9,},
4576 {0x3d4d, 0x7c3d, 0x36ba, 0x0d2b,
4577 0xfdc2, 0xcefc, 0x8461, 0x7711, 0xabcc, 0x3fe4,},
4578 {0xc155, 0xa4a8, 0x404e, 0x6113,
4579 0xd3c3, 0x652b, 0xe219, 0x1758, 0xd1b7, 0x3ff1,},
4580 {0xd70a, 0x70a3, 0x0a3d, 0xa3d7,
4581 0x3d70, 0xd70a, 0x70a3, 0x0a3d, 0xa3d7, 0x3ff8,},
4582 {0xcccd, 0xcccc, 0xcccc, 0xcccc,
4583 0xcccc, 0xcccc, 0xcccc, 0xcccc, 0xcccc, 0x3ffb,}, /* 10**-1 */
4586 /* LONG_DOUBLE_TYPE_SIZE is other than 128 */
4587 static unsigned EMUSHORT etens[NTEN + 1][NE] =
4589 {0xc94c, 0x979a, 0x8a20, 0x5202, 0xc460, 0x7525,}, /* 10**4096 */
4590 {0xa74d, 0x5de4, 0xc53d, 0x3b5d, 0x9e8b, 0x5a92,}, /* 10**2048 */
4591 {0x650d, 0x0c17, 0x8175, 0x7586, 0xc976, 0x4d48,},
4592 {0xcc65, 0x91c6, 0xa60e, 0xa0ae, 0xe319, 0x46a3,},
4593 {0xddbc, 0xde8d, 0x9df9, 0xebfb, 0xaa7e, 0x4351,},
4594 {0xc66f, 0x8cdf, 0x80e9, 0x47c9, 0x93ba, 0x41a8,},
4595 {0x3cbf, 0xa6d5, 0xffcf, 0x1f49, 0xc278, 0x40d3,},
4596 {0xf020, 0xb59d, 0x2b70, 0xada8, 0x9dc5, 0x4069,},
4597 {0x0000, 0x0000, 0x0400, 0xc9bf, 0x8e1b, 0x4034,},
4598 {0x0000, 0x0000, 0x0000, 0x2000, 0xbebc, 0x4019,},
4599 {0x0000, 0x0000, 0x0000, 0x0000, 0x9c40, 0x400c,},
4600 {0x0000, 0x0000, 0x0000, 0x0000, 0xc800, 0x4005,},
4601 {0x0000, 0x0000, 0x0000, 0x0000, 0xa000, 0x4002,}, /* 10**1 */
4604 static unsigned EMUSHORT emtens[NTEN + 1][NE] =
4606 {0x2de4, 0x9fde, 0xd2ce, 0x04c8, 0xa6dd, 0x0ad8,}, /* 10**-4096 */
4607 {0x4925, 0x2de4, 0x3436, 0x534f, 0xceae, 0x256b,}, /* 10**-2048 */
4608 {0x87a6, 0xc0bd, 0xda57, 0x82a5, 0xa2a6, 0x32b5,},
4609 {0x7133, 0xd21c, 0xdb23, 0xee32, 0x9049, 0x395a,},
4610 {0xfa91, 0x1939, 0x637a, 0x4325, 0xc031, 0x3cac,},
4611 {0xac7d, 0xe4a0, 0x64bc, 0x467c, 0xddd0, 0x3e55,},
4612 {0x3f24, 0xe9a5, 0xa539, 0xea27, 0xa87f, 0x3f2a,},
4613 {0x67de, 0x94ba, 0x4539, 0x1ead, 0xcfb1, 0x3f94,},
4614 {0x4c2f, 0xe15b, 0xc44d, 0x94be, 0xe695, 0x3fc9,},
4615 {0xfdc2, 0xcefc, 0x8461, 0x7711, 0xabcc, 0x3fe4,},
4616 {0xd3c3, 0x652b, 0xe219, 0x1758, 0xd1b7, 0x3ff1,},
4617 {0x3d71, 0xd70a, 0x70a3, 0x0a3d, 0xa3d7, 0x3ff8,},
4618 {0xcccd, 0xcccc, 0xcccc, 0xcccc, 0xcccc, 0x3ffb,}, /* 10**-1 */
4623 /* Convert float value X to ASCII string STRING with NDIG digits after
4624 the decimal point. */
4627 e24toasc (x, string, ndigs)
4628 unsigned EMUSHORT x[];
4632 unsigned EMUSHORT w[NI];
4635 etoasc (w, string, ndigs);
4638 /* Convert double value X to ASCII string STRING with NDIG digits after
4639 the decimal point. */
4642 e53toasc (x, string, ndigs)
4643 unsigned EMUSHORT x[];
4647 unsigned EMUSHORT w[NI];
4650 etoasc (w, string, ndigs);
4653 /* Convert double extended value X to ASCII string STRING with NDIG digits
4654 after the decimal point. */
4657 e64toasc (x, string, ndigs)
4658 unsigned EMUSHORT x[];
4662 unsigned EMUSHORT w[NI];
4665 etoasc (w, string, ndigs);
4668 /* Convert 128-bit long double value X to ASCII string STRING with NDIG digits
4669 after the decimal point. */
4672 e113toasc (x, string, ndigs)
4673 unsigned EMUSHORT x[];
4677 unsigned EMUSHORT w[NI];
4680 etoasc (w, string, ndigs);
4684 /* Convert e-type X to ASCII string STRING with NDIGS digits after
4685 the decimal point. */
4687 static char wstring[80]; /* working storage for ASCII output */
4690 etoasc (x, string, ndigs)
4691 unsigned EMUSHORT x[];
4696 unsigned EMUSHORT y[NI], t[NI], u[NI], w[NI];
4697 unsigned EMUSHORT *p, *r, *ten;
4698 unsigned EMUSHORT sign;
4699 int i, j, k, expon, rndsav;
4701 unsigned EMUSHORT m;
4712 sprintf (wstring, " NaN ");
4716 rndprc = NBITS; /* set to full precision */
4717 emov (x, y); /* retain external format */
4718 if (y[NE - 1] & 0x8000)
4721 y[NE - 1] &= 0x7fff;
4728 ten = &etens[NTEN][0];
4730 /* Test for zero exponent */
4733 for (k = 0; k < NE - 1; k++)
4736 goto tnzro; /* denormalized number */
4738 goto isone; /* valid all zeros */
4742 /* Test for infinity. */
4743 if (y[NE - 1] == 0x7fff)
4746 sprintf (wstring, " -Infinity ");
4748 sprintf (wstring, " Infinity ");
4752 /* Test for exponent nonzero but significand denormalized.
4753 * This is an error condition.
4755 if ((y[NE - 1] != 0) && ((y[NE - 2] & 0x8000) == 0))
4757 mtherr ("etoasc", DOMAIN);
4758 sprintf (wstring, "NaN");
4762 /* Compare to 1.0 */
4771 { /* Number is greater than 1 */
4772 /* Convert significand to an integer and strip trailing decimal zeros. */
4774 u[NE - 1] = EXONE + NBITS - 1;
4776 p = &etens[NTEN - 4][0];
4782 for (j = 0; j < NE - 1; j++)
4795 /* Rescale from integer significand */
4796 u[NE - 1] += y[NE - 1] - (unsigned int) (EXONE + NBITS - 1);
4798 /* Find power of 10 */
4802 /* An unordered compare result shouldn't happen here. */
4803 while (ecmp (ten, u) <= 0)
4805 if (ecmp (p, u) <= 0)
4818 { /* Number is less than 1.0 */
4819 /* Pad significand with trailing decimal zeros. */
4822 while ((y[NE - 2] & 0x8000) == 0)
4831 for (i = 0; i < NDEC + 1; i++)
4833 if ((w[NI - 1] & 0x7) != 0)
4835 /* multiply by 10 */
4848 if (eone[NE - 1] <= u[1])
4860 while (ecmp (eone, w) > 0)
4862 if (ecmp (p, w) >= 0)
4877 /* Find the first (leading) digit. */
4883 digit = equot[NI - 1];
4884 while ((digit == 0) && (ecmp (y, ezero) != 0))
4892 digit = equot[NI - 1];
4900 /* Examine number of digits requested by caller. */
4918 *s++ = (char)digit + '0';
4921 /* Generate digits after the decimal point. */
4922 for (k = 0; k <= ndigs; k++)
4924 /* multiply current number by 10, without normalizing */
4931 *s++ = (char) equot[NI - 1] + '0';
4933 digit = equot[NI - 1];
4936 /* round off the ASCII string */
4939 /* Test for critical rounding case in ASCII output. */
4943 if (ecmp (t, ezero) != 0)
4944 goto roun; /* round to nearest */
4946 if ((*(s - 1) & 1) == 0)
4947 goto doexp; /* round to even */
4950 /* Round up and propagate carry-outs */
4954 /* Carry out to most significant digit? */
4961 /* Most significant digit carries to 10? */
4969 /* Round up and carry out from less significant digits */
4981 sprintf (ss, "e+%d", expon);
4983 sprintf (ss, "e%d", expon);
4985 sprintf (ss, "e%d", expon);
4988 /* copy out the working string */
4991 while (*ss == ' ') /* strip possible leading space */
4993 while ((*s++ = *ss++) != '\0')
4998 /* Convert ASCII string to floating point.
5000 Numeric input is a free format decimal number of any length, with
5001 or without decimal point. Entering E after the number followed by an
5002 integer number causes the second number to be interpreted as a power of
5003 10 to be multiplied by the first number (i.e., "scientific" notation). */
5005 /* Convert ASCII string S to single precision float value Y. */
5010 unsigned EMUSHORT *y;
5016 /* Convert ASCII string S to double precision value Y. */
5021 unsigned EMUSHORT *y;
5023 #if defined(DEC) || defined(IBM)
5035 /* Convert ASCII string S to double extended value Y. */
5040 unsigned EMUSHORT *y;
5045 /* Convert ASCII string S to 128-bit long double Y. */
5050 unsigned EMUSHORT *y;
5052 asctoeg (s, y, 113);
5055 /* Convert ASCII string S to e type Y. */
5060 unsigned EMUSHORT *y;
5062 asctoeg (s, y, NBITS);
5065 /* Convert ASCII string SS to e type Y, with a specified rounding precision
5066 of OPREC bits. BASE is 16 for C9X hexadecimal floating constants. */
5069 asctoeg (ss, y, oprec)
5071 unsigned EMUSHORT *y;
5074 unsigned EMUSHORT yy[NI], xt[NI], tt[NI];
5075 int esign, decflg, sgnflg, nexp, exp, prec, lost;
5076 int k, trail, c, rndsav;
5078 unsigned EMUSHORT nsign, *p;
5079 char *sp, *s, *lstr;
5082 /* Copy the input string. */
5083 lstr = (char *) alloca (strlen (ss) + 1);
5085 while (*ss == ' ') /* skip leading spaces */
5089 while ((*sp++ = *ss++) != '\0')
5093 if (s[0] == '0' && (s[1] == 'x' || s[1] == 'X'))
5100 rndprc = NBITS; /* Set to full precision */
5112 if (*s >= '0' && *s <= '9')
5118 if ((k >= 0) && (k < base))
5120 /* Ignore leading zeros */
5121 if ((prec == 0) && (decflg == 0) && (k == 0))
5123 /* Identify and strip trailing zeros after the decimal point. */
5124 if ((trail == 0) && (decflg != 0))
5127 while ((*sp >= '0' && *sp <= '9')
5128 || (base == 16 && ((*sp >= 'a' && *sp <= 'f')
5129 || (*sp >= 'A' && *sp <= 'F'))))
5131 /* Check for syntax error */
5133 if ((base != 10 || ((c != 'e') && (c != 'E')))
5134 && (base != 16 || ((c != 'p') && (c != 'P')))
5136 && (c != '\n') && (c != '\r') && (c != ' ')
5147 /* If enough digits were given to more than fill up the yy register,
5148 continuing until overflow into the high guard word yy[2]
5149 guarantees that there will be a roundoff bit at the top
5150 of the low guard word after normalization. */
5157 nexp += 4; /* count digits after decimal point */
5159 eshup1 (yy); /* multiply current number by 16 */
5167 nexp += 1; /* count digits after decimal point */
5169 eshup1 (yy); /* multiply current number by 10 */
5175 /* Insert the current digit. */
5177 xt[NI - 2] = (unsigned EMUSHORT) k;
5182 /* Mark any lost non-zero digit. */
5184 /* Count lost digits before the decimal point. */
5206 case '.': /* decimal point */
5236 mtherr ("asctoe", DOMAIN);
5245 /* Exponent interpretation */
5247 /* 0.0eXXX is zero, regardless of XXX. Check for the 0.0. */
5248 for (k = 0; k < NI; k++)
5259 /* check for + or - */
5267 while ((*s >= '0') && (*s <= '9'))
5276 if ((exp > MAXDECEXP) && (base == 10))
5280 yy[E] = 0x7fff; /* infinity */
5283 if ((exp < MINDECEXP) && (base == 10))
5293 /* Base 16 hexadecimal floating constant. */
5294 if ((k = enormlz (yy)) > NBITS)
5299 /* Adjust the exponent. NEXP is the number of hex digits,
5300 EXP is a power of 2. */
5301 lexp = (EXONE - 1 + NBITS) - k + yy[E] + exp - nexp;
5311 /* Pad trailing zeros to minimize power of 10, per IEEE spec. */
5312 while ((nexp > 0) && (yy[2] == 0))
5324 if ((k = enormlz (yy)) > NBITS)
5329 lexp = (EXONE - 1 + NBITS) - k;
5330 emdnorm (yy, lost, 0, lexp, 64);
5333 /* Convert to external format:
5335 Multiply by 10**nexp. If precision is 64 bits,
5336 the maximum relative error incurred in forming 10**n
5337 for 0 <= n <= 324 is 8.2e-20, at 10**180.
5338 For 0 <= n <= 999, the peak relative error is 1.4e-19 at 10**947.
5339 For 0 >= n >= -999, it is -1.55e-19 at 10**-435. */
5354 /* Punt. Can't handle this without 2 divides. */
5355 emovi (etens[0], tt);
5362 p = &etens[NTEN][0];
5372 while (exp <= MAXP);
5391 /* Round and convert directly to the destination type */
5393 lexp -= EXONE - 0x3ff;
5395 else if (oprec == 24 || oprec == 32)
5396 lexp -= (EXONE - 0x7f);
5399 else if (oprec == 24 || oprec == 56)
5400 lexp -= EXONE - (0x41 << 2);
5402 else if (oprec == 24)
5403 lexp -= EXONE - 0177;
5407 else if (oprec == 56)
5408 lexp -= EXONE - 0201;
5411 emdnorm (yy, lost, 0, lexp, 64);
5421 todec (yy, y); /* see etodec.c */
5426 toibm (yy, y, DFmode);
5431 toc4x (yy, y, HFmode);
5455 /* Return Y = largest integer not greater than X (truncated toward minus
5458 static unsigned EMUSHORT bmask[] =
5481 unsigned EMUSHORT x[], y[];
5483 register unsigned EMUSHORT *p;
5485 unsigned EMUSHORT f[NE];
5487 emov (x, f); /* leave in external format */
5488 expon = (int) f[NE - 1];
5489 e = (expon & 0x7fff) - (EXONE - 1);
5495 /* number of bits to clear out */
5507 /* clear the remaining bits */
5509 /* truncate negatives toward minus infinity */
5512 if ((unsigned EMUSHORT) expon & (unsigned EMUSHORT) 0x8000)
5514 for (i = 0; i < NE - 1; i++)
5527 /* Return S and EXP such that S * 2^EXP = X and .5 <= S < 1.
5528 For example, 1.1 = 0.55 * 2^1. */
5532 unsigned EMUSHORT x[];
5534 unsigned EMUSHORT s[];
5536 unsigned EMUSHORT xi[NI];
5540 /* Handle denormalized numbers properly using long integer exponent. */
5541 li = (EMULONG) ((EMUSHORT) xi[1]);
5549 *exp = (int) (li - 0x3ffe);
5553 /* Return e type Y = X * 2^PWR2. */
5557 unsigned EMUSHORT x[];
5559 unsigned EMUSHORT y[];
5561 unsigned EMUSHORT xi[NI];
5569 emdnorm (xi, i, i, li, 64);
5575 /* C = remainder after dividing B by A, all e type values.
5576 Least significant integer quotient bits left in EQUOT. */
5580 unsigned EMUSHORT a[], b[], c[];
5582 unsigned EMUSHORT den[NI], num[NI];
5586 || (ecmp (a, ezero) == 0)
5594 if (ecmp (a, ezero) == 0)
5596 mtherr ("eremain", SING);
5602 eiremain (den, num);
5603 /* Sign of remainder = sign of quotient */
5612 /* Return quotient of exploded e-types NUM / DEN in EQUOT,
5613 remainder in NUM. */
5617 unsigned EMUSHORT den[], num[];
5620 unsigned EMUSHORT j;
5623 ld -= enormlz (den);
5625 ln -= enormlz (num);
5629 if (ecmpm (den, num) <= 0)
5641 emdnorm (num, 0, 0, ln, 0);
5644 /* Report an error condition CODE encountered in function NAME.
5646 Mnemonic Value Significance
5648 DOMAIN 1 argument domain error
5649 SING 2 function singularity
5650 OVERFLOW 3 overflow range error
5651 UNDERFLOW 4 underflow range error
5652 TLOSS 5 total loss of precision
5653 PLOSS 6 partial loss of precision
5654 INVALID 7 NaN - producing operation
5655 EDOM 33 Unix domain error code
5656 ERANGE 34 Unix range error code
5658 The order of appearance of the following messages is bound to the
5659 error codes defined above. */
5669 /* The string passed by the calling program is supposed to be the
5670 name of the function in which the error occurred.
5671 The code argument selects which error message string will be printed. */
5673 if (strcmp (name, "esub") == 0)
5674 name = "subtraction";
5675 else if (strcmp (name, "ediv") == 0)
5677 else if (strcmp (name, "emul") == 0)
5678 name = "multiplication";
5679 else if (strcmp (name, "enormlz") == 0)
5680 name = "normalization";
5681 else if (strcmp (name, "etoasc") == 0)
5682 name = "conversion to text";
5683 else if (strcmp (name, "asctoe") == 0)
5685 else if (strcmp (name, "eremain") == 0)
5687 else if (strcmp (name, "esqrt") == 0)
5688 name = "square root";
5693 case DOMAIN: warning ("%s: argument domain error" , name); break;
5694 case SING: warning ("%s: function singularity" , name); break;
5695 case OVERFLOW: warning ("%s: overflow range error" , name); break;
5696 case UNDERFLOW: warning ("%s: underflow range error" , name); break;
5697 case TLOSS: warning ("%s: total loss of precision" , name); break;
5698 case PLOSS: warning ("%s: partial loss of precision", name); break;
5699 case INVALID: warning ("%s: NaN - producing operation", name); break;
5704 /* Set global error message word */
5709 /* Convert DEC double precision D to e type E. */
5713 unsigned EMUSHORT *d;
5714 unsigned EMUSHORT *e;
5716 unsigned EMUSHORT y[NI];
5717 register unsigned EMUSHORT r, *p;
5719 ecleaz (y); /* start with a zero */
5720 p = y; /* point to our number */
5721 r = *d; /* get DEC exponent word */
5722 if (*d & (unsigned int) 0x8000)
5723 *p = 0xffff; /* fill in our sign */
5724 ++p; /* bump pointer to our exponent word */
5725 r &= 0x7fff; /* strip the sign bit */
5726 if (r == 0) /* answer = 0 if high order DEC word = 0 */
5730 r >>= 7; /* shift exponent word down 7 bits */
5731 r += EXONE - 0201; /* subtract DEC exponent offset */
5732 /* add our e type exponent offset */
5733 *p++ = r; /* to form our exponent */
5735 r = *d++; /* now do the high order mantissa */
5736 r &= 0177; /* strip off the DEC exponent and sign bits */
5737 r |= 0200; /* the DEC understood high order mantissa bit */
5738 *p++ = r; /* put result in our high guard word */
5740 *p++ = *d++; /* fill in the rest of our mantissa */
5744 eshdn8 (y); /* shift our mantissa down 8 bits */
5749 /* Convert e type X to DEC double precision D. */
5753 unsigned EMUSHORT *x, *d;
5755 unsigned EMUSHORT xi[NI];
5760 /* Adjust exponent for offsets. */
5761 exp = (EMULONG) xi[E] - (EXONE - 0201);
5762 /* Round off to nearest or even. */
5765 emdnorm (xi, 0, 0, exp, 64);
5770 /* Convert exploded e-type X, that has already been rounded to
5771 56-bit precision, to DEC format double Y. */
5775 unsigned EMUSHORT *x, *y;
5777 unsigned EMUSHORT i;
5778 unsigned EMUSHORT *p;
5817 /* Convert IBM single/double precision to e type. */
5821 unsigned EMUSHORT *d;
5822 unsigned EMUSHORT *e;
5823 enum machine_mode mode;
5825 unsigned EMUSHORT y[NI];
5826 register unsigned EMUSHORT r, *p;
5829 ecleaz (y); /* start with a zero */
5830 p = y; /* point to our number */
5831 r = *d; /* get IBM exponent word */
5832 if (*d & (unsigned int) 0x8000)
5833 *p = 0xffff; /* fill in our sign */
5834 ++p; /* bump pointer to our exponent word */
5835 r &= 0x7f00; /* strip the sign bit */
5836 r >>= 6; /* shift exponent word down 6 bits */
5837 /* in fact shift by 8 right and 2 left */
5838 r += EXONE - (0x41 << 2); /* subtract IBM exponent offset */
5839 /* add our e type exponent offset */
5840 *p++ = r; /* to form our exponent */
5842 *p++ = *d++ & 0xff; /* now do the high order mantissa */
5843 /* strip off the IBM exponent and sign bits */
5844 if (mode != SFmode) /* there are only 2 words in SFmode */
5846 *p++ = *d++; /* fill in the rest of our mantissa */
5851 if (y[M] == 0 && y[M+1] == 0 && y[M+2] == 0 && y[M+3] == 0)
5854 y[E] -= 5 + enormlz (y); /* now normalise the mantissa */
5855 /* handle change in RADIX */
5861 /* Convert e type to IBM single/double precision. */
5865 unsigned EMUSHORT *x, *d;
5866 enum machine_mode mode;
5868 unsigned EMUSHORT xi[NI];
5873 exp = (EMULONG) xi[E] - (EXONE - (0x41 << 2)); /* adjust exponent for offsets */
5874 /* round off to nearest or even */
5877 emdnorm (xi, 0, 0, exp, 64);
5879 toibm (xi, d, mode);
5884 unsigned EMUSHORT *x, *y;
5885 enum machine_mode mode;
5887 unsigned EMUSHORT i;
5888 unsigned EMUSHORT *p;
5938 /* Convert C4X single/double precision to e type. */
5942 unsigned EMUSHORT *d;
5943 unsigned EMUSHORT *e;
5944 enum machine_mode mode;
5946 unsigned EMUSHORT y[NI];
5953 /* Short-circuit the zero case. */
5954 if ((d[0] == 0x8000)
5956 && ((mode == QFmode) || ((d[2] == 0x0000) && (d[3] == 0x0000))))
5967 ecleaz (y); /* start with a zero */
5968 r = d[0]; /* get sign/exponent part */
5969 if (r & (unsigned int) 0x0080)
5971 y[0] = 0xffff; /* fill in our sign */
5979 r >>= 8; /* Shift exponent word down 8 bits. */
5980 if (r & 0x80) /* Make the exponent negative if it is. */
5982 r = r | (~0 & ~0xff);
5987 /* Now do the high order mantissa. We don't "or" on the high bit
5988 because it is 2 (not 1) and is handled a little differently
5993 if (mode != QFmode) /* There are only 2 words in QFmode. */
5995 y[M+2] = d[2]; /* Fill in the rest of our mantissa. */
6005 /* Now do the two's complement on the data. */
6007 carry = 1; /* Initially add 1 for the two's complement. */
6008 for (i=size + M; i > M; i--)
6010 if (carry && (y[i] == 0x0000))
6012 /* We overflowed into the next word, carry is the same. */
6013 y[i] = carry ? 0x0000 : 0xffff;
6017 /* No overflow, just invert and add carry. */
6018 y[i] = ((~y[i]) + carry) & 0xffff;
6033 /* Add our e type exponent offset to form our exponent. */
6037 /* Now do the high order mantissa strip off the exponent and sign
6038 bits and add the high 1 bit. */
6039 y[M] = (d[0] & 0x7f) | 0x80;
6042 if (mode != QFmode) /* There are only 2 words in QFmode. */
6044 y[M+2] = d[2]; /* Fill in the rest of our mantissa. */
6054 /* Convert e type to C4X single/double precision. */
6058 unsigned EMUSHORT *x, *d;
6059 enum machine_mode mode;
6061 unsigned EMUSHORT xi[NI];
6067 /* Adjust exponent for offsets. */
6068 exp = (EMULONG) xi[E] - (EXONE - 0x7f);
6070 /* Round off to nearest or even. */
6072 rndprc = mode == QFmode ? 24 : 32;
6073 emdnorm (xi, 0, 0, exp, 64);
6075 toc4x (xi, d, mode);
6080 unsigned EMUSHORT *x, *y;
6081 enum machine_mode mode;
6087 /* Short-circuit the zero case */
6088 if ((x[0] == 0) /* Zero exponent and sign */
6090 && (x[M] == 0) /* The rest is for zero mantissa */
6092 /* Only check for double if necessary */
6093 && ((mode == QFmode) || ((x[M+2] == 0) && (x[M+3] == 0))))
6095 /* We have a zero. Put it into the output and return. */
6108 /* Negative number require a two's complement conversion of the
6114 i = ((int) x[1]) - 0x7f;
6116 /* Now add 1 to the inverted data to do the two's complement. */
6126 x[v] = carry ? 0x0000 : 0xffff;
6130 x[v] = ((~x[v]) + carry) & 0xffff;
6136 /* The following is a special case. The C4X negative float requires
6137 a zero in the high bit (because the format is (2 - x) x 2^m), so
6138 if a one is in that bit, we have to shift left one to get rid
6139 of it. This only occurs if the number is -1 x 2^m. */
6140 if (x[M+1] & 0x8000)
6142 /* This is the case of -1 x 2^m, we have to rid ourselves of the
6143 high sign bit and shift the exponent. */
6150 i = ((int) x[1]) - 0x7f;
6153 if ((i < -128) || (i > 127))
6168 y[0] |= ((i & 0xff) << 8);
6172 y[0] |= x[M] & 0x7f;
6182 /* Output a binary NaN bit pattern in the target machine's format. */
6184 /* If special NaN bit patterns are required, define them in tm.h
6185 as arrays of unsigned 16-bit shorts. Otherwise, use the default
6191 unsigned EMUSHORT TFbignan[8] =
6192 {0x7fff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff};
6193 unsigned EMUSHORT TFlittlenan[8] = {0, 0, 0, 0, 0, 0, 0x8000, 0xffff};
6201 unsigned EMUSHORT XFbignan[6] =
6202 {0x7fff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff};
6203 unsigned EMUSHORT XFlittlenan[6] = {0, 0, 0, 0xc000, 0xffff, 0};
6211 unsigned EMUSHORT DFbignan[4] = {0x7fff, 0xffff, 0xffff, 0xffff};
6212 unsigned EMUSHORT DFlittlenan[4] = {0, 0, 0, 0xfff8};
6220 unsigned EMUSHORT SFbignan[2] = {0x7fff, 0xffff};
6221 unsigned EMUSHORT SFlittlenan[2] = {0, 0xffc0};
6227 make_nan (nan, sign, mode)
6228 unsigned EMUSHORT *nan;
6230 enum machine_mode mode;
6233 unsigned EMUSHORT *p;
6237 /* Possibly the `reserved operand' patterns on a VAX can be
6238 used like NaN's, but probably not in the same way as IEEE. */
6239 #if !defined(DEC) && !defined(IBM) && !defined(C4X)
6242 if (REAL_WORDS_BIG_ENDIAN)
6250 if (REAL_WORDS_BIG_ENDIAN)
6258 if (REAL_WORDS_BIG_ENDIAN)
6267 if (REAL_WORDS_BIG_ENDIAN)
6277 if (REAL_WORDS_BIG_ENDIAN)
6278 *nan++ = (sign << 15) | (*p++ & 0x7fff);
6281 if (! REAL_WORDS_BIG_ENDIAN)
6282 *nan = (sign << 15) | (*p & 0x7fff);
6285 /* This is the inverse of the function `etarsingle' invoked by
6286 REAL_VALUE_TO_TARGET_SINGLE. */
6289 ereal_unto_float (f)
6293 unsigned EMUSHORT s[2];
6294 unsigned EMUSHORT e[NE];
6296 /* Convert 32 bit integer to array of 16 bit pieces in target machine order.
6297 This is the inverse operation to what the function `endian' does. */
6298 if (REAL_WORDS_BIG_ENDIAN)
6300 s[0] = (unsigned EMUSHORT) (f >> 16);
6301 s[1] = (unsigned EMUSHORT) f;
6305 s[0] = (unsigned EMUSHORT) f;
6306 s[1] = (unsigned EMUSHORT) (f >> 16);
6308 /* Convert and promote the target float to E-type. */
6310 /* Output E-type to REAL_VALUE_TYPE. */
6316 /* This is the inverse of the function `etardouble' invoked by
6317 REAL_VALUE_TO_TARGET_DOUBLE. */
6320 ereal_unto_double (d)
6324 unsigned EMUSHORT s[4];
6325 unsigned EMUSHORT e[NE];
6327 /* Convert array of HOST_WIDE_INT to equivalent array of 16-bit pieces. */
6328 if (REAL_WORDS_BIG_ENDIAN)
6330 s[0] = (unsigned EMUSHORT) (d[0] >> 16);
6331 s[1] = (unsigned EMUSHORT) d[0];
6332 s[2] = (unsigned EMUSHORT) (d[1] >> 16);
6333 s[3] = (unsigned EMUSHORT) d[1];
6337 /* Target float words are little-endian. */
6338 s[0] = (unsigned EMUSHORT) d[0];
6339 s[1] = (unsigned EMUSHORT) (d[0] >> 16);
6340 s[2] = (unsigned EMUSHORT) d[1];
6341 s[3] = (unsigned EMUSHORT) (d[1] >> 16);
6343 /* Convert target double to E-type. */
6345 /* Output E-type to REAL_VALUE_TYPE. */
6351 /* Convert an SFmode target `float' value to a REAL_VALUE_TYPE.
6352 This is somewhat like ereal_unto_float, but the input types
6353 for these are different. */
6356 ereal_from_float (f)
6360 unsigned EMUSHORT s[2];
6361 unsigned EMUSHORT e[NE];
6363 /* Convert 32 bit integer to array of 16 bit pieces in target machine order.
6364 This is the inverse operation to what the function `endian' does. */
6365 if (REAL_WORDS_BIG_ENDIAN)
6367 s[0] = (unsigned EMUSHORT) (f >> 16);
6368 s[1] = (unsigned EMUSHORT) f;
6372 s[0] = (unsigned EMUSHORT) f;
6373 s[1] = (unsigned EMUSHORT) (f >> 16);
6375 /* Convert and promote the target float to E-type. */
6377 /* Output E-type to REAL_VALUE_TYPE. */
6383 /* Convert a DFmode target `double' value to a REAL_VALUE_TYPE.
6384 This is somewhat like ereal_unto_double, but the input types
6385 for these are different.
6387 The DFmode is stored as an array of HOST_WIDE_INT in the target's
6388 data format, with no holes in the bit packing. The first element
6389 of the input array holds the bits that would come first in the
6390 target computer's memory. */
6393 ereal_from_double (d)
6397 unsigned EMUSHORT s[4];
6398 unsigned EMUSHORT e[NE];
6400 /* Convert array of HOST_WIDE_INT to equivalent array of 16-bit pieces. */
6401 if (REAL_WORDS_BIG_ENDIAN)
6403 #if HOST_BITS_PER_WIDE_INT == 32
6404 s[0] = (unsigned EMUSHORT) (d[0] >> 16);
6405 s[1] = (unsigned EMUSHORT) d[0];
6406 s[2] = (unsigned EMUSHORT) (d[1] >> 16);
6407 s[3] = (unsigned EMUSHORT) d[1];
6409 /* In this case the entire target double is contained in the
6410 first array element. The second element of the input is
6412 s[0] = (unsigned EMUSHORT) (d[0] >> 48);
6413 s[1] = (unsigned EMUSHORT) (d[0] >> 32);
6414 s[2] = (unsigned EMUSHORT) (d[0] >> 16);
6415 s[3] = (unsigned EMUSHORT) d[0];
6420 /* Target float words are little-endian. */
6421 s[0] = (unsigned EMUSHORT) d[0];
6422 s[1] = (unsigned EMUSHORT) (d[0] >> 16);
6423 #if HOST_BITS_PER_WIDE_INT == 32
6424 s[2] = (unsigned EMUSHORT) d[1];
6425 s[3] = (unsigned EMUSHORT) (d[1] >> 16);
6427 s[2] = (unsigned EMUSHORT) (d[0] >> 32);
6428 s[3] = (unsigned EMUSHORT) (d[0] >> 48);
6431 /* Convert target double to E-type. */
6433 /* Output E-type to REAL_VALUE_TYPE. */
6440 /* Convert target computer unsigned 64-bit integer to e-type.
6441 The endian-ness of DImode follows the convention for integers,
6442 so we use WORDS_BIG_ENDIAN here, not REAL_WORDS_BIG_ENDIAN. */
6446 unsigned EMUSHORT *di; /* Address of the 64-bit int. */
6447 unsigned EMUSHORT *e;
6449 unsigned EMUSHORT yi[NI];
6453 if (WORDS_BIG_ENDIAN)
6455 for (k = M; k < M + 4; k++)
6460 for (k = M + 3; k >= M; k--)
6463 yi[E] = EXONE + 47; /* exponent if normalize shift count were 0 */
6464 if ((k = enormlz (yi)) > NBITS)/* normalize the significand */
6465 ecleaz (yi); /* it was zero */
6467 yi[E] -= (unsigned EMUSHORT) k;/* subtract shift count from exponent */
6471 /* Convert target computer signed 64-bit integer to e-type. */
6475 unsigned EMUSHORT *di; /* Address of the 64-bit int. */
6476 unsigned EMUSHORT *e;
6478 unsigned EMULONG acc;
6479 unsigned EMUSHORT yi[NI];
6480 unsigned EMUSHORT carry;
6484 if (WORDS_BIG_ENDIAN)
6486 for (k = M; k < M + 4; k++)
6491 for (k = M + 3; k >= M; k--)
6494 /* Take absolute value */
6500 for (k = M + 3; k >= M; k--)
6502 acc = (unsigned EMULONG) (~yi[k] & 0xffff) + carry;
6509 yi[E] = EXONE + 47; /* exponent if normalize shift count were 0 */
6510 if ((k = enormlz (yi)) > NBITS)/* normalize the significand */
6511 ecleaz (yi); /* it was zero */
6513 yi[E] -= (unsigned EMUSHORT) k;/* subtract shift count from exponent */
6520 /* Convert e-type to unsigned 64-bit int. */
6524 unsigned EMUSHORT *x;
6525 unsigned EMUSHORT *i;
6527 unsigned EMUSHORT xi[NI];
6536 k = (int) xi[E] - (EXONE - 1);
6539 for (j = 0; j < 4; j++)
6545 for (j = 0; j < 4; j++)
6548 warning ("overflow on truncation to integer");
6553 /* Shift more than 16 bits: first shift up k-16 mod 16,
6554 then shift up by 16's. */
6555 j = k - ((k >> 4) << 4);
6559 if (WORDS_BIG_ENDIAN)
6570 if (WORDS_BIG_ENDIAN)
6575 while ((k -= 16) > 0);
6579 /* shift not more than 16 bits */
6584 if (WORDS_BIG_ENDIAN)
6603 /* Convert e-type to signed 64-bit int. */
6607 unsigned EMUSHORT *x;
6608 unsigned EMUSHORT *i;
6610 unsigned EMULONG acc;
6611 unsigned EMUSHORT xi[NI];
6612 unsigned EMUSHORT carry;
6613 unsigned EMUSHORT *isave;
6617 k = (int) xi[E] - (EXONE - 1);
6620 for (j = 0; j < 4; j++)
6626 for (j = 0; j < 4; j++)
6629 warning ("overflow on truncation to integer");
6635 /* Shift more than 16 bits: first shift up k-16 mod 16,
6636 then shift up by 16's. */
6637 j = k - ((k >> 4) << 4);
6641 if (WORDS_BIG_ENDIAN)
6652 if (WORDS_BIG_ENDIAN)
6657 while ((k -= 16) > 0);
6661 /* shift not more than 16 bits */
6664 if (WORDS_BIG_ENDIAN)
6680 /* Negate if negative */
6684 if (WORDS_BIG_ENDIAN)
6686 for (k = 0; k < 4; k++)
6688 acc = (unsigned EMULONG) (~(*isave) & 0xffff) + carry;
6689 if (WORDS_BIG_ENDIAN)
6701 /* Longhand square root routine. */
6704 static int esqinited = 0;
6705 static unsigned short sqrndbit[NI];
6709 unsigned EMUSHORT *x, *y;
6711 unsigned EMUSHORT temp[NI], num[NI], sq[NI], xx[NI];
6713 int i, j, k, n, nlups;
6718 sqrndbit[NI - 2] = 1;
6721 /* Check for arg <= 0 */
6722 i = ecmp (x, ezero);
6727 mtherr ("esqrt", DOMAIN);
6743 /* Bring in the arg and renormalize if it is denormal. */
6745 m = (EMULONG) xx[1]; /* local long word exponent */
6749 /* Divide exponent by 2 */
6751 exp = (unsigned short) ((m / 2) + 0x3ffe);
6753 /* Adjust if exponent odd */
6763 n = 8; /* get 8 bits of result per inner loop */
6769 /* bring in next word of arg */
6771 num[NI - 1] = xx[j + 3];
6772 /* Do additional bit on last outer loop, for roundoff. */
6775 for (i = 0; i < n; i++)
6777 /* Next 2 bits of arg */
6780 /* Shift up answer */
6782 /* Make trial divisor */
6783 for (k = 0; k < NI; k++)
6786 eaddm (sqrndbit, temp);
6787 /* Subtract and insert answer bit if it goes in */
6788 if (ecmpm (temp, num) <= 0)
6798 /* Adjust for extra, roundoff loop done. */
6799 exp += (NBITS - 1) - rndprc;
6801 /* Sticky bit = 1 if the remainder is nonzero. */
6803 for (i = 3; i < NI; i++)
6806 /* Renormalize and round off. */
6807 emdnorm (sq, k, 0, exp, 64);
6811 #endif /* EMU_NON_COMPILE not defined */
6813 /* Return the binary precision of the significand for a given
6814 floating point mode. The mode can hold an integer value
6815 that many bits wide, without losing any bits. */
6818 significand_size (mode)
6819 enum machine_mode mode;
6822 /* Don't test the modes, but their sizes, lest this
6823 code won't work for BITS_PER_UNIT != 8 . */
6825 switch (GET_MODE_BITSIZE (mode))
6829 #if TARGET_FLOAT_FORMAT == C4X_FLOAT_FORMAT
6836 #if TARGET_FLOAT_FORMAT == IEEE_FLOAT_FORMAT
6839 #if TARGET_FLOAT_FORMAT == IBM_FLOAT_FORMAT
6842 #if TARGET_FLOAT_FORMAT == VAX_FLOAT_FORMAT
6845 #if TARGET_FLOAT_FORMAT == C4X_FLOAT_FORMAT