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, 95, 96, 97, 1998 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.
84 If LONG_DOUBLE_TYPE_SIZE = 64 (the default, unless tm.h defines it)
85 then `long double' and `double' are both implemented, but they
86 both mean DFmode. In this case, the software floating-point
87 support available here is activated by writing
88 #define REAL_ARITHMETIC
91 The case LONG_DOUBLE_TYPE_SIZE = 128 activates TFmode support
92 and may deactivate XFmode since `long double' is used to refer
95 The macros FLOAT_WORDS_BIG_ENDIAN, HOST_FLOAT_WORDS_BIG_ENDIAN,
96 contributed by Richard Earnshaw <Richard.Earnshaw@cl.cam.ac.uk>,
97 separate the floating point unit's endian-ness from that of
98 the integer addressing. This permits one to define a big-endian
99 FPU on a little-endian machine (e.g., ARM). An extension to
100 BYTES_BIG_ENDIAN may be required for some machines in the future.
101 These optional macros may be defined in tm.h. In real.h, they
102 default to WORDS_BIG_ENDIAN, etc., so there is no need to define
103 them for any normal host or target machine on which the floats
104 and the integers have the same endian-ness. */
107 /* The following converts gcc macros into the ones used by this file. */
109 /* REAL_ARITHMETIC defined means that macros in real.h are
110 defined to call emulator functions. */
111 #ifdef REAL_ARITHMETIC
113 #if TARGET_FLOAT_FORMAT == VAX_FLOAT_FORMAT
114 /* PDP-11, Pro350, VAX: */
116 #else /* it's not VAX */
117 #if TARGET_FLOAT_FORMAT == IBM_FLOAT_FORMAT
118 /* IBM System/370 style */
120 #else /* it's also not an IBM */
121 #if TARGET_FLOAT_FORMAT == C4X_FLOAT_FORMAT
122 /* TMS320C3x/C4x style */
124 #else /* it's also not a C4X */
125 #if TARGET_FLOAT_FORMAT == IEEE_FLOAT_FORMAT
127 #else /* it's not IEEE either */
128 /* UNKnown arithmetic. We don't support this and can't go on. */
129 unknown arithmetic type
131 #endif /* not IEEE */
136 #define REAL_WORDS_BIG_ENDIAN FLOAT_WORDS_BIG_ENDIAN
139 /* REAL_ARITHMETIC not defined means that the *host's* data
140 structure will be used. It may differ by endian-ness from the
141 target machine's structure and will get its ends swapped
142 accordingly (but not here). Probably only the decimal <-> binary
143 functions in this file will actually be used in this case. */
145 #if HOST_FLOAT_FORMAT == VAX_FLOAT_FORMAT
147 #else /* it's not VAX */
148 #if HOST_FLOAT_FORMAT == IBM_FLOAT_FORMAT
149 /* IBM System/370 style */
151 #else /* it's also not an IBM */
152 #if HOST_FLOAT_FORMAT == IEEE_FLOAT_FORMAT
154 #else /* it's not IEEE either */
155 unknown arithmetic type
157 #endif /* not IEEE */
161 #define REAL_WORDS_BIG_ENDIAN HOST_FLOAT_WORDS_BIG_ENDIAN
163 #endif /* REAL_ARITHMETIC not defined */
165 /* Define INFINITY for support of infinity.
166 Define NANS for support of Not-a-Number's (NaN's). */
167 #if !defined(DEC) && !defined(IBM) && !defined(C4X)
172 /* Support of NaNs requires support of infinity. */
179 /* Find a host integer type that is at least 16 bits wide,
180 and another type at least twice whatever that size is. */
182 #if HOST_BITS_PER_CHAR >= 16
183 #define EMUSHORT char
184 #define EMUSHORT_SIZE HOST_BITS_PER_CHAR
185 #define EMULONG_SIZE (2 * HOST_BITS_PER_CHAR)
187 #if HOST_BITS_PER_SHORT >= 16
188 #define EMUSHORT short
189 #define EMUSHORT_SIZE HOST_BITS_PER_SHORT
190 #define EMULONG_SIZE (2 * HOST_BITS_PER_SHORT)
192 #if HOST_BITS_PER_INT >= 16
194 #define EMUSHORT_SIZE HOST_BITS_PER_INT
195 #define EMULONG_SIZE (2 * HOST_BITS_PER_INT)
197 #if HOST_BITS_PER_LONG >= 16
198 #define EMUSHORT long
199 #define EMUSHORT_SIZE HOST_BITS_PER_LONG
200 #define EMULONG_SIZE (2 * HOST_BITS_PER_LONG)
202 /* You will have to modify this program to have a smaller unit size. */
203 #define EMU_NON_COMPILE
209 #if HOST_BITS_PER_SHORT >= EMULONG_SIZE
210 #define EMULONG short
212 #if HOST_BITS_PER_INT >= EMULONG_SIZE
215 #if HOST_BITS_PER_LONG >= EMULONG_SIZE
218 #if HOST_BITS_PER_LONGLONG >= EMULONG_SIZE
219 #define EMULONG long long int
221 /* You will have to modify this program to have a smaller unit size. */
222 #define EMU_NON_COMPILE
229 /* The host interface doesn't work if no 16-bit size exists. */
230 #if EMUSHORT_SIZE != 16
231 #define EMU_NON_COMPILE
234 /* OK to continue compilation. */
235 #ifndef EMU_NON_COMPILE
237 /* Construct macros to translate between REAL_VALUE_TYPE and e type.
238 In GET_REAL and PUT_REAL, r and e are pointers.
239 A REAL_VALUE_TYPE is guaranteed to occupy contiguous locations
240 in memory, with no holes. */
242 #if LONG_DOUBLE_TYPE_SIZE == 96
243 /* Number of 16 bit words in external e type format */
245 #define MAXDECEXP 4932
246 #define MINDECEXP -4956
247 #define GET_REAL(r,e) bcopy ((char *) r, (char *) e, 2*NE)
248 #define PUT_REAL(e,r) bcopy ((char *) e, (char *) r, 2*NE)
249 #else /* no XFmode */
250 #if LONG_DOUBLE_TYPE_SIZE == 128
252 #define MAXDECEXP 4932
253 #define MINDECEXP -4977
254 #define GET_REAL(r,e) bcopy ((char *) r, (char *) e, 2*NE)
255 #define PUT_REAL(e,r) bcopy ((char *) e, (char *) r, 2*NE)
258 #define MAXDECEXP 4932
259 #define MINDECEXP -4956
260 #ifdef REAL_ARITHMETIC
261 /* Emulator uses target format internally
262 but host stores it in host endian-ness. */
264 #define GET_REAL(r,e) \
266 if (HOST_FLOAT_WORDS_BIG_ENDIAN == REAL_WORDS_BIG_ENDIAN) \
267 e53toe ((unsigned EMUSHORT *) (r), (e)); \
270 unsigned EMUSHORT w[4]; \
271 w[3] = ((EMUSHORT *) r)[0]; \
272 w[2] = ((EMUSHORT *) r)[1]; \
273 w[1] = ((EMUSHORT *) r)[2]; \
274 w[0] = ((EMUSHORT *) r)[3]; \
279 #define PUT_REAL(e,r) \
281 if (HOST_FLOAT_WORDS_BIG_ENDIAN == REAL_WORDS_BIG_ENDIAN) \
282 etoe53 ((e), (unsigned EMUSHORT *) (r)); \
285 unsigned EMUSHORT w[4]; \
287 *((EMUSHORT *) r) = w[3]; \
288 *((EMUSHORT *) r + 1) = w[2]; \
289 *((EMUSHORT *) r + 2) = w[1]; \
290 *((EMUSHORT *) r + 3) = w[0]; \
294 #else /* not REAL_ARITHMETIC */
296 /* emulator uses host format */
297 #define GET_REAL(r,e) e53toe ((unsigned EMUSHORT *) (r), (e))
298 #define PUT_REAL(e,r) etoe53 ((e), (unsigned EMUSHORT *) (r))
300 #endif /* not REAL_ARITHMETIC */
301 #endif /* not TFmode */
302 #endif /* not XFmode */
305 /* Number of 16 bit words in internal format */
308 /* Array offset to exponent */
311 /* Array offset to high guard word */
314 /* Number of bits of precision */
315 #define NBITS ((NI-4)*16)
317 /* Maximum number of decimal digits in ASCII conversion
320 #define NDEC (NBITS*8/27)
322 /* The exponent of 1.0 */
323 #define EXONE (0x3fff)
325 extern int extra_warnings;
326 extern unsigned EMUSHORT ezero[], ehalf[], eone[], etwo[];
327 extern unsigned EMUSHORT elog2[], esqrt2[];
329 static void endian PROTO((unsigned EMUSHORT *, long *,
331 static void eclear PROTO((unsigned EMUSHORT *));
332 static void emov PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
334 static void eabs PROTO((unsigned EMUSHORT *));
336 static void eneg PROTO((unsigned EMUSHORT *));
337 static int eisneg PROTO((unsigned EMUSHORT *));
338 static int eisinf PROTO((unsigned EMUSHORT *));
339 static int eisnan PROTO((unsigned EMUSHORT *));
340 static void einfin PROTO((unsigned EMUSHORT *));
341 static void enan PROTO((unsigned EMUSHORT *, int));
342 static void emovi PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
343 static void emovo PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
344 static void ecleaz PROTO((unsigned EMUSHORT *));
345 static void ecleazs PROTO((unsigned EMUSHORT *));
346 static void emovz PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
347 static void einan PROTO((unsigned EMUSHORT *));
348 static int eiisnan PROTO((unsigned EMUSHORT *));
349 static int eiisneg PROTO((unsigned EMUSHORT *));
351 static void eiinfin PROTO((unsigned EMUSHORT *));
353 static int eiisinf PROTO((unsigned EMUSHORT *));
354 static int ecmpm PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
355 static void eshdn1 PROTO((unsigned EMUSHORT *));
356 static void eshup1 PROTO((unsigned EMUSHORT *));
357 static void eshdn8 PROTO((unsigned EMUSHORT *));
358 static void eshup8 PROTO((unsigned EMUSHORT *));
359 static void eshup6 PROTO((unsigned EMUSHORT *));
360 static void eshdn6 PROTO((unsigned EMUSHORT *));
361 static void eaddm PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
\f
362 static void esubm PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
363 static void m16m PROTO((unsigned int, unsigned short *,
365 static int edivm PROTO((unsigned short *, unsigned short *));
366 static int emulm PROTO((unsigned short *, unsigned short *));
367 static void emdnorm PROTO((unsigned EMUSHORT *, int, int, EMULONG, int));
368 static void esub PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *,
369 unsigned EMUSHORT *));
370 static void eadd PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *,
371 unsigned EMUSHORT *));
372 static void eadd1 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *,
373 unsigned EMUSHORT *));
374 static void ediv PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *,
375 unsigned EMUSHORT *));
376 static void emul PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *,
377 unsigned EMUSHORT *));
378 static void e53toe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
379 static void e64toe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
380 static void e113toe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
381 static void e24toe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
382 static void etoe113 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
383 static void toe113 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
384 static void etoe64 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
385 static void toe64 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
386 static void etoe53 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
387 static void toe53 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
388 static void etoe24 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
389 static void toe24 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
390 static int ecmp PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
392 static void eround PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
394 static void ltoe PROTO((HOST_WIDE_INT *, unsigned EMUSHORT *));
395 static void ultoe PROTO((unsigned HOST_WIDE_INT *, unsigned EMUSHORT *));
396 static void eifrac PROTO((unsigned EMUSHORT *, HOST_WIDE_INT *,
397 unsigned EMUSHORT *));
398 static void euifrac PROTO((unsigned EMUSHORT *, unsigned HOST_WIDE_INT *,
399 unsigned EMUSHORT *));
400 static int eshift PROTO((unsigned EMUSHORT *, int));
401 static int enormlz PROTO((unsigned EMUSHORT *));
403 static void e24toasc PROTO((unsigned EMUSHORT *, char *, int));
404 static void e53toasc PROTO((unsigned EMUSHORT *, char *, int));
405 static void e64toasc PROTO((unsigned EMUSHORT *, char *, int));
406 static void e113toasc PROTO((unsigned EMUSHORT *, char *, int));
408 static void etoasc PROTO((unsigned EMUSHORT *, char *, int));
409 static void asctoe24 PROTO((char *, unsigned EMUSHORT *));
410 static void asctoe53 PROTO((char *, unsigned EMUSHORT *));
411 static void asctoe64 PROTO((char *, unsigned EMUSHORT *));
412 static void asctoe113 PROTO((char *, unsigned EMUSHORT *));
413 static void asctoe PROTO((char *, unsigned EMUSHORT *));
414 static void asctoeg PROTO((char *, unsigned EMUSHORT *, int));
415 static void efloor PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
417 static void efrexp PROTO((unsigned EMUSHORT *, int *,
418 unsigned EMUSHORT *));
420 static void eldexp PROTO((unsigned EMUSHORT *, int, unsigned EMUSHORT *));
422 static void eremain PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *,
423 unsigned EMUSHORT *));
425 static void eiremain PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
426 static void mtherr PROTO((char *, int));
428 static void dectoe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
429 static void etodec PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
430 static void todec PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
433 static void ibmtoe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *,
435 static void etoibm PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *,
437 static void toibm PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *,
441 static void c4xtoe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *,
443 static void etoc4x PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *,
445 static void toc4x PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *,
448 static void make_nan PROTO((unsigned EMUSHORT *, int, enum machine_mode));
450 static void uditoe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
451 static void ditoe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
452 static void etoudi PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
453 static void etodi PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
454 static void esqrt PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
457 /* Copy 32-bit numbers obtained from array containing 16-bit numbers,
458 swapping ends if required, into output array of longs. The
459 result is normally passed to fprintf by the ASM_OUTPUT_ macros. */
463 unsigned EMUSHORT e[];
465 enum machine_mode mode;
469 if (REAL_WORDS_BIG_ENDIAN)
474 /* Swap halfwords in the fourth long. */
475 th = (unsigned long) e[6] & 0xffff;
476 t = (unsigned long) e[7] & 0xffff;
481 /* Swap halfwords in the third long. */
482 th = (unsigned long) e[4] & 0xffff;
483 t = (unsigned long) e[5] & 0xffff;
486 /* fall into the double case */
489 /* Swap halfwords in the second word. */
490 th = (unsigned long) e[2] & 0xffff;
491 t = (unsigned long) e[3] & 0xffff;
494 /* fall into the float case */
498 /* Swap halfwords in the first word. */
499 th = (unsigned long) e[0] & 0xffff;
500 t = (unsigned long) e[1] & 0xffff;
511 /* Pack the output array without swapping. */
516 /* Pack the fourth long. */
517 th = (unsigned long) e[7] & 0xffff;
518 t = (unsigned long) e[6] & 0xffff;
523 /* Pack the third long.
524 Each element of the input REAL_VALUE_TYPE array has 16 useful bits
526 th = (unsigned long) e[5] & 0xffff;
527 t = (unsigned long) e[4] & 0xffff;
530 /* fall into the double case */
533 /* Pack the second long */
534 th = (unsigned long) e[3] & 0xffff;
535 t = (unsigned long) e[2] & 0xffff;
538 /* fall into the float case */
542 /* Pack the first long */
543 th = (unsigned long) e[1] & 0xffff;
544 t = (unsigned long) e[0] & 0xffff;
556 /* This is the implementation of the REAL_ARITHMETIC macro. */
559 earith (value, icode, r1, r2)
560 REAL_VALUE_TYPE *value;
565 unsigned EMUSHORT d1[NE], d2[NE], v[NE];
571 /* Return NaN input back to the caller. */
574 PUT_REAL (d1, value);
579 PUT_REAL (d2, value);
583 code = (enum tree_code) icode;
591 esub (d2, d1, v); /* d1 - d2 */
599 #ifndef REAL_INFINITY
600 if (ecmp (d2, ezero) == 0)
603 enan (v, eisneg (d1) ^ eisneg (d2));
610 ediv (d2, d1, v); /* d1/d2 */
613 case MIN_EXPR: /* min (d1,d2) */
614 if (ecmp (d1, d2) < 0)
620 case MAX_EXPR: /* max (d1,d2) */
621 if (ecmp (d1, d2) > 0)
634 /* Truncate REAL_VALUE_TYPE toward zero to signed HOST_WIDE_INT.
635 implements REAL_VALUE_RNDZINT (x) (etrunci (x)). */
641 unsigned EMUSHORT f[NE], g[NE];
657 /* Truncate REAL_VALUE_TYPE toward zero to unsigned HOST_WIDE_INT;
658 implements REAL_VALUE_UNSIGNED_RNDZINT (x) (etruncui (x)). */
664 unsigned EMUSHORT f[NE], g[NE];
666 unsigned HOST_WIDE_INT l;
680 /* This is the REAL_VALUE_ATOF function. It converts a decimal string to
681 binary, rounding off as indicated by the machine_mode argument. Then it
682 promotes the rounded value to REAL_VALUE_TYPE. */
689 unsigned EMUSHORT tem[NE], e[NE];
732 /* Expansion of REAL_NEGATE. */
738 unsigned EMUSHORT e[NE];
748 /* Round real toward zero to HOST_WIDE_INT;
749 implements REAL_VALUE_FIX (x). */
755 unsigned EMUSHORT f[NE], g[NE];
762 warning ("conversion from NaN to int");
770 /* Round real toward zero to unsigned HOST_WIDE_INT
771 implements REAL_VALUE_UNSIGNED_FIX (x).
772 Negative input returns zero. */
774 unsigned HOST_WIDE_INT
778 unsigned EMUSHORT f[NE], g[NE];
779 unsigned HOST_WIDE_INT l;
785 warning ("conversion from NaN to unsigned int");
794 /* REAL_VALUE_FROM_INT macro. */
797 ereal_from_int (d, i, j, mode)
800 enum machine_mode mode;
802 unsigned EMUSHORT df[NE], dg[NE];
803 HOST_WIDE_INT low, high;
806 if (GET_MODE_CLASS (mode) != MODE_FLOAT)
813 /* complement and add 1 */
820 eldexp (eone, HOST_BITS_PER_WIDE_INT, df);
821 ultoe ((unsigned HOST_WIDE_INT *) &high, dg);
823 ultoe ((unsigned HOST_WIDE_INT *) &low, df);
828 /* A REAL_VALUE_TYPE may not be wide enough to hold the two HOST_WIDE_INTS.
829 Avoid double-rounding errors later by rounding off now from the
830 extra-wide internal format to the requested precision. */
831 switch (GET_MODE_BITSIZE (mode))
861 /* REAL_VALUE_FROM_UNSIGNED_INT macro. */
864 ereal_from_uint (d, i, j, mode)
866 unsigned HOST_WIDE_INT i, j;
867 enum machine_mode mode;
869 unsigned EMUSHORT df[NE], dg[NE];
870 unsigned HOST_WIDE_INT low, high;
872 if (GET_MODE_CLASS (mode) != MODE_FLOAT)
876 eldexp (eone, HOST_BITS_PER_WIDE_INT, df);
882 /* A REAL_VALUE_TYPE may not be wide enough to hold the two HOST_WIDE_INTS.
883 Avoid double-rounding errors later by rounding off now from the
884 extra-wide internal format to the requested precision. */
885 switch (GET_MODE_BITSIZE (mode))
915 /* REAL_VALUE_TO_INT macro. */
918 ereal_to_int (low, high, rr)
919 HOST_WIDE_INT *low, *high;
922 unsigned EMUSHORT d[NE], df[NE], dg[NE], dh[NE];
929 warning ("conversion from NaN to int");
935 /* convert positive value */
942 eldexp (eone, HOST_BITS_PER_WIDE_INT, df);
943 ediv (df, d, dg); /* dg = d / 2^32 is the high word */
944 euifrac (dg, (unsigned HOST_WIDE_INT *) high, dh);
945 emul (df, dh, dg); /* fractional part is the low word */
946 euifrac (dg, (unsigned HOST_WIDE_INT *)low, dh);
949 /* complement and add 1 */
959 /* REAL_VALUE_LDEXP macro. */
966 unsigned EMUSHORT e[NE], y[NE];
979 /* These routines are conditionally compiled because functions
980 of the same names may be defined in fold-const.c. */
982 #ifdef REAL_ARITHMETIC
984 /* Check for infinity in a REAL_VALUE_TYPE. */
990 unsigned EMUSHORT e[NE];
1000 /* Check whether a REAL_VALUE_TYPE item is a NaN. */
1006 unsigned EMUSHORT e[NE];
1010 return (eisnan (e));
1017 /* Check for a negative REAL_VALUE_TYPE number.
1018 This just checks the sign bit, so that -0 counts as negative. */
1024 return ereal_isneg (x);
1027 /* Expansion of REAL_VALUE_TRUNCATE.
1028 The result is in floating point, rounded to nearest or even. */
1031 real_value_truncate (mode, arg)
1032 enum machine_mode mode;
1033 REAL_VALUE_TYPE arg;
1035 unsigned EMUSHORT e[NE], t[NE];
1081 /* If an unsupported type was requested, presume that
1082 the machine files know something useful to do with
1083 the unmodified value. */
1092 /* Try to change R into its exact multiplicative inverse in machine mode
1093 MODE. Return nonzero function value if successful. */
1096 exact_real_inverse (mode, r)
1097 enum machine_mode mode;
1100 unsigned EMUSHORT e[NE], einv[NE];
1101 REAL_VALUE_TYPE rinv;
1106 /* Test for input in range. Don't transform IEEE special values. */
1107 if (eisinf (e) || eisnan (e) || (ecmp (e, ezero) == 0))
1110 /* Test for a power of 2: all significand bits zero except the MSB.
1111 We are assuming the target has binary (or hex) arithmetic. */
1112 if (e[NE - 2] != 0x8000)
1115 for (i = 0; i < NE - 2; i++)
1121 /* Compute the inverse and truncate it to the required mode. */
1122 ediv (e, eone, einv);
1123 PUT_REAL (einv, &rinv);
1124 rinv = real_value_truncate (mode, rinv);
1126 #ifdef CHECK_FLOAT_VALUE
1127 /* This check is not redundant. It may, for example, flush
1128 a supposedly IEEE denormal value to zero. */
1130 if (CHECK_FLOAT_VALUE (mode, rinv, i))
1133 GET_REAL (&rinv, einv);
1135 /* Check the bits again, because the truncation might have
1136 generated an arbitrary saturation value on overflow. */
1137 if (einv[NE - 2] != 0x8000)
1140 for (i = 0; i < NE - 2; i++)
1146 /* Fail if the computed inverse is out of range. */
1147 if (eisinf (einv) || eisnan (einv) || (ecmp (einv, ezero) == 0))
1150 /* Output the reciprocal and return success flag. */
1154 #endif /* REAL_ARITHMETIC defined */
1156 /* Used for debugging--print the value of R in human-readable format
1165 REAL_VALUE_TO_DECIMAL (r, "%.20g", dstr);
1166 fprintf (stderr, "%s", dstr);
1170 /* The following routines convert REAL_VALUE_TYPE to the various floating
1171 point formats that are meaningful to supported computers.
1173 The results are returned in 32-bit pieces, each piece stored in a `long'.
1174 This is so they can be printed by statements like
1176 fprintf (file, "%lx, %lx", L[0], L[1]);
1178 that will work on both narrow- and wide-word host computers. */
1180 /* Convert R to a 128-bit long double precision value. The output array L
1181 contains four 32-bit pieces of the result, in the order they would appear
1189 unsigned EMUSHORT e[NE];
1193 endian (e, l, TFmode);
1196 /* Convert R to a double extended precision value. The output array L
1197 contains three 32-bit pieces of the result, in the order they would
1198 appear in memory. */
1205 unsigned EMUSHORT e[NE];
1209 endian (e, l, XFmode);
1212 /* Convert R to a double precision value. The output array L contains two
1213 32-bit pieces of the result, in the order they would appear in memory. */
1220 unsigned EMUSHORT e[NE];
1224 endian (e, l, DFmode);
1227 /* Convert R to a single precision float value stored in the least-significant
1228 bits of a `long'. */
1234 unsigned EMUSHORT e[NE];
1239 endian (e, &l, SFmode);
1243 /* Convert X to a decimal ASCII string S for output to an assembly
1244 language file. Note, there is no standard way to spell infinity or
1245 a NaN, so these values may require special treatment in the tm.h
1249 ereal_to_decimal (x, s)
1253 unsigned EMUSHORT e[NE];
1259 /* Compare X and Y. Return 1 if X > Y, 0 if X == Y, -1 if X < Y,
1260 or -2 if either is a NaN. */
1264 REAL_VALUE_TYPE x, y;
1266 unsigned EMUSHORT ex[NE], ey[NE];
1270 return (ecmp (ex, ey));
1273 /* Return 1 if the sign bit of X is set, else return 0. */
1279 unsigned EMUSHORT ex[NE];
1282 return (eisneg (ex));
1285 /* End of REAL_ARITHMETIC interface */
1288 Extended precision IEEE binary floating point arithmetic routines
1290 Numbers are stored in C language as arrays of 16-bit unsigned
1291 short integers. The arguments of the routines are pointers to
1294 External e type data structure, similar to Intel 8087 chip
1295 temporary real format but possibly with a larger significand:
1297 NE-1 significand words (least significant word first,
1298 most significant bit is normally set)
1299 exponent (value = EXONE for 1.0,
1300 top bit is the sign)
1303 Internal exploded e-type data structure of a number (a "word" is 16 bits):
1305 ei[0] sign word (0 for positive, 0xffff for negative)
1306 ei[1] biased exponent (value = EXONE for the number 1.0)
1307 ei[2] high guard word (always zero after normalization)
1309 to ei[NI-2] significand (NI-4 significand words,
1310 most significant word first,
1311 most significant bit is set)
1312 ei[NI-1] low guard word (0x8000 bit is rounding place)
1316 Routines for external format e-type numbers
1318 asctoe (string, e) ASCII string to extended double e type
1319 asctoe64 (string, &d) ASCII string to long double
1320 asctoe53 (string, &d) ASCII string to double
1321 asctoe24 (string, &f) ASCII string to single
1322 asctoeg (string, e, prec) ASCII string to specified precision
1323 e24toe (&f, e) IEEE single precision to e type
1324 e53toe (&d, e) IEEE double precision to e type
1325 e64toe (&d, e) IEEE long double precision to e type
1326 e113toe (&d, e) 128-bit long double precision to e type
1328 eabs (e) absolute value
1330 eadd (a, b, c) c = b + a
1332 ecmp (a, b) Returns 1 if a > b, 0 if a == b,
1333 -1 if a < b, -2 if either a or b is a NaN.
1334 ediv (a, b, c) c = b / a
1335 efloor (a, b) truncate to integer, toward -infinity
1336 efrexp (a, exp, s) extract exponent and significand
1337 eifrac (e, &l, frac) e to HOST_WIDE_INT and e type fraction
1338 euifrac (e, &l, frac) e to unsigned HOST_WIDE_INT and e type fraction
1339 einfin (e) set e to infinity, leaving its sign alone
1340 eldexp (a, n, b) multiply by 2**n
1342 emul (a, b, c) c = b * a
1345 eround (a, b) b = nearest integer value to a
1347 esub (a, b, c) c = b - a
1349 e24toasc (&f, str, n) single to ASCII string, n digits after decimal
1350 e53toasc (&d, str, n) double to ASCII string, n digits after decimal
1351 e64toasc (&d, str, n) 80-bit long double to ASCII string
1352 e113toasc (&d, str, n) 128-bit long double to ASCII string
1354 etoasc (e, str, n) e to ASCII string, n digits after decimal
1355 etoe24 (e, &f) convert e type to IEEE single precision
1356 etoe53 (e, &d) convert e type to IEEE double precision
1357 etoe64 (e, &d) convert e type to IEEE long double precision
1358 ltoe (&l, e) HOST_WIDE_INT to e type
1359 ultoe (&l, e) unsigned HOST_WIDE_INT to e type
1360 eisneg (e) 1 if sign bit of e != 0, else 0
1361 eisinf (e) 1 if e has maximum exponent (non-IEEE)
1362 or is infinite (IEEE)
1363 eisnan (e) 1 if e is a NaN
1366 Routines for internal format exploded e-type numbers
1368 eaddm (ai, bi) add significands, bi = bi + ai
1370 ecleazs (ei) set ei = 0 but leave its sign alone
1371 ecmpm (ai, bi) compare significands, return 1, 0, or -1
1372 edivm (ai, bi) divide significands, bi = bi / ai
1373 emdnorm (ai,l,s,exp) normalize and round off
1374 emovi (a, ai) convert external a to internal ai
1375 emovo (ai, a) convert internal ai to external a
1376 emovz (ai, bi) bi = ai, low guard word of bi = 0
1377 emulm (ai, bi) multiply significands, bi = bi * ai
1378 enormlz (ei) left-justify the significand
1379 eshdn1 (ai) shift significand and guards down 1 bit
1380 eshdn8 (ai) shift down 8 bits
1381 eshdn6 (ai) shift down 16 bits
1382 eshift (ai, n) shift ai n bits up (or down if n < 0)
1383 eshup1 (ai) shift significand and guards up 1 bit
1384 eshup8 (ai) shift up 8 bits
1385 eshup6 (ai) shift up 16 bits
1386 esubm (ai, bi) subtract significands, bi = bi - ai
1387 eiisinf (ai) 1 if infinite
1388 eiisnan (ai) 1 if a NaN
1389 eiisneg (ai) 1 if sign bit of ai != 0, else 0
1390 einan (ai) set ai = NaN
1392 eiinfin (ai) set ai = infinity
1395 The result is always normalized and rounded to NI-4 word precision
1396 after each arithmetic operation.
1398 Exception flags are NOT fully supported.
1400 Signaling NaN's are NOT supported; they are treated the same
1403 Define INFINITY for support of infinity; otherwise a
1404 saturation arithmetic is implemented.
1406 Define NANS for support of Not-a-Number items; otherwise the
1407 arithmetic will never produce a NaN output, and might be confused
1409 If NaN's are supported, the output of `ecmp (a,b)' is -2 if
1410 either a or b is a NaN. This means asking `if (ecmp (a,b) < 0)'
1411 may not be legitimate. Use `if (ecmp (a,b) == -1)' for `less than'
1414 Denormals are always supported here where appropriate (e.g., not
1415 for conversion to DEC numbers). */
1417 /* Definitions for error codes that are passed to the common error handling
1420 For Digital Equipment PDP-11 and VAX computers, certain
1421 IBM systems, and others that use numbers with a 56-bit
1422 significand, the symbol DEC should be defined. In this
1423 mode, most floating point constants are given as arrays
1424 of octal integers to eliminate decimal to binary conversion
1425 errors that might be introduced by the compiler.
1427 For computers, such as IBM PC, that follow the IEEE
1428 Standard for Binary Floating Point Arithmetic (ANSI/IEEE
1429 Std 754-1985), the symbol IEEE should be defined.
1430 These numbers have 53-bit significands. In this mode, constants
1431 are provided as arrays of hexadecimal 16 bit integers.
1432 The endian-ness of generated values is controlled by
1433 REAL_WORDS_BIG_ENDIAN.
1435 To accommodate other types of computer arithmetic, all
1436 constants are also provided in a normal decimal radix
1437 which one can hope are correctly converted to a suitable
1438 format by the available C language compiler. To invoke
1439 this mode, the symbol UNK is defined.
1441 An important difference among these modes is a predefined
1442 set of machine arithmetic constants for each. The numbers
1443 MACHEP (the machine roundoff error), MAXNUM (largest number
1444 represented), and several other parameters are preset by
1445 the configuration symbol. Check the file const.c to
1446 ensure that these values are correct for your computer.
1448 For ANSI C compatibility, define ANSIC equal to 1. Currently
1449 this affects only the atan2 function and others that use it. */
1451 /* Constant definitions for math error conditions. */
1453 #define DOMAIN 1 /* argument domain error */
1454 #define SING 2 /* argument singularity */
1455 #define OVERFLOW 3 /* overflow range error */
1456 #define UNDERFLOW 4 /* underflow range error */
1457 #define TLOSS 5 /* total loss of precision */
1458 #define PLOSS 6 /* partial loss of precision */
1459 #define INVALID 7 /* NaN-producing operation */
1461 /* e type constants used by high precision check routines */
1463 #if LONG_DOUBLE_TYPE_SIZE == 128
1465 unsigned EMUSHORT ezero[NE] =
1466 {0x0000, 0x0000, 0x0000, 0x0000,
1467 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000,};
1468 extern unsigned EMUSHORT ezero[];
1471 unsigned EMUSHORT ehalf[NE] =
1472 {0x0000, 0x0000, 0x0000, 0x0000,
1473 0x0000, 0x0000, 0x0000, 0x0000, 0x8000, 0x3ffe,};
1474 extern unsigned EMUSHORT ehalf[];
1477 unsigned EMUSHORT eone[NE] =
1478 {0x0000, 0x0000, 0x0000, 0x0000,
1479 0x0000, 0x0000, 0x0000, 0x0000, 0x8000, 0x3fff,};
1480 extern unsigned EMUSHORT eone[];
1483 unsigned EMUSHORT etwo[NE] =
1484 {0x0000, 0x0000, 0x0000, 0x0000,
1485 0x0000, 0x0000, 0x0000, 0x0000, 0x8000, 0x4000,};
1486 extern unsigned EMUSHORT etwo[];
1489 unsigned EMUSHORT e32[NE] =
1490 {0x0000, 0x0000, 0x0000, 0x0000,
1491 0x0000, 0x0000, 0x0000, 0x0000, 0x8000, 0x4004,};
1492 extern unsigned EMUSHORT e32[];
1494 /* 6.93147180559945309417232121458176568075500134360255E-1 */
1495 unsigned EMUSHORT elog2[NE] =
1496 {0x40f3, 0xf6af, 0x03f2, 0xb398,
1497 0xc9e3, 0x79ab, 0150717, 0013767, 0130562, 0x3ffe,};
1498 extern unsigned EMUSHORT elog2[];
1500 /* 1.41421356237309504880168872420969807856967187537695E0 */
1501 unsigned EMUSHORT esqrt2[NE] =
1502 {0x1d6f, 0xbe9f, 0x754a, 0x89b3,
1503 0x597d, 0x6484, 0174736, 0171463, 0132404, 0x3fff,};
1504 extern unsigned EMUSHORT esqrt2[];
1506 /* 3.14159265358979323846264338327950288419716939937511E0 */
1507 unsigned EMUSHORT epi[NE] =
1508 {0x2902, 0x1cd1, 0x80dc, 0x628b,
1509 0xc4c6, 0xc234, 0020550, 0155242, 0144417, 0040000,};
1510 extern unsigned EMUSHORT epi[];
1513 /* LONG_DOUBLE_TYPE_SIZE is other than 128 */
1514 unsigned EMUSHORT ezero[NE] =
1515 {0, 0000000, 0000000, 0000000, 0000000, 0000000,};
1516 unsigned EMUSHORT ehalf[NE] =
1517 {0, 0000000, 0000000, 0000000, 0100000, 0x3ffe,};
1518 unsigned EMUSHORT eone[NE] =
1519 {0, 0000000, 0000000, 0000000, 0100000, 0x3fff,};
1520 unsigned EMUSHORT etwo[NE] =
1521 {0, 0000000, 0000000, 0000000, 0100000, 0040000,};
1522 unsigned EMUSHORT e32[NE] =
1523 {0, 0000000, 0000000, 0000000, 0100000, 0040004,};
1524 unsigned EMUSHORT elog2[NE] =
1525 {0xc9e4, 0x79ab, 0150717, 0013767, 0130562, 0x3ffe,};
1526 unsigned EMUSHORT esqrt2[NE] =
1527 {0x597e, 0x6484, 0174736, 0171463, 0132404, 0x3fff,};
1528 unsigned EMUSHORT epi[NE] =
1529 {0xc4c6, 0xc234, 0020550, 0155242, 0144417, 0040000,};
1532 /* Control register for rounding precision.
1533 This can be set to 113 (if NE=10), 80 (if NE=6), 64, 56, 53, or 24 bits. */
1538 /* Clear out entire e-type number X. */
1542 register unsigned EMUSHORT *x;
1546 for (i = 0; i < NE; i++)
1550 /* Move e-type number from A to B. */
1554 register unsigned EMUSHORT *a, *b;
1558 for (i = 0; i < NE; i++)
1564 /* Absolute value of e-type X. */
1568 unsigned EMUSHORT x[];
1570 /* sign is top bit of last word of external format */
1571 x[NE - 1] &= 0x7fff;
1575 /* Negate the e-type number X. */
1579 unsigned EMUSHORT x[];
1582 x[NE - 1] ^= 0x8000; /* Toggle the sign bit */
1585 /* Return 1 if sign bit of e-type number X is nonzero, else zero. */
1589 unsigned EMUSHORT x[];
1592 if (x[NE - 1] & 0x8000)
1598 /* Return 1 if e-type number X is infinity, else return zero. */
1602 unsigned EMUSHORT x[];
1609 if ((x[NE - 1] & 0x7fff) == 0x7fff)
1615 /* Check if e-type number is not a number. The bit pattern is one that we
1616 defined, so we know for sure how to detect it. */
1620 unsigned EMUSHORT x[];
1625 /* NaN has maximum exponent */
1626 if ((x[NE - 1] & 0x7fff) != 0x7fff)
1628 /* ... and non-zero significand field. */
1629 for (i = 0; i < NE - 1; i++)
1639 /* Fill e-type number X with infinity pattern (IEEE)
1640 or largest possible number (non-IEEE). */
1644 register unsigned EMUSHORT *x;
1649 for (i = 0; i < NE - 1; i++)
1653 for (i = 0; i < NE - 1; i++)
1681 /* Output an e-type NaN.
1682 This generates Intel's quiet NaN pattern for extended real.
1683 The exponent is 7fff, the leading mantissa word is c000. */
1687 register unsigned EMUSHORT *x;
1692 for (i = 0; i < NE - 2; i++)
1695 *x = (sign << 15) | 0x7fff;
1698 /* Move in an e-type number A, converting it to exploded e-type B. */
1702 unsigned EMUSHORT *a, *b;
1704 register unsigned EMUSHORT *p, *q;
1708 p = a + (NE - 1); /* point to last word of external number */
1709 /* get the sign bit */
1714 /* get the exponent */
1716 *q++ &= 0x7fff; /* delete the sign bit */
1718 if ((*(q - 1) & 0x7fff) == 0x7fff)
1724 for (i = 3; i < NI; i++)
1730 for (i = 2; i < NI; i++)
1736 /* clear high guard word */
1738 /* move in the significand */
1739 for (i = 0; i < NE - 1; i++)
1741 /* clear low guard word */
1745 /* Move out exploded e-type number A, converting it to e type B. */
1749 unsigned EMUSHORT *a, *b;
1751 register unsigned EMUSHORT *p, *q;
1752 unsigned EMUSHORT i;
1756 q = b + (NE - 1); /* point to output exponent */
1757 /* combine sign and exponent */
1760 *q-- = *p++ | 0x8000;
1764 if (*(p - 1) == 0x7fff)
1769 enan (b, eiisneg (a));
1777 /* skip over guard word */
1779 /* move the significand */
1780 for (j = 0; j < NE - 1; j++)
1784 /* Clear out exploded e-type number XI. */
1788 register unsigned EMUSHORT *xi;
1792 for (i = 0; i < NI; i++)
1796 /* Clear out exploded e-type XI, but don't touch the sign. */
1800 register unsigned EMUSHORT *xi;
1805 for (i = 0; i < NI - 1; i++)
1809 /* Move exploded e-type number from A to B. */
1813 register unsigned EMUSHORT *a, *b;
1817 for (i = 0; i < NI - 1; i++)
1819 /* clear low guard word */
1823 /* Generate exploded e-type NaN.
1824 The explicit pattern for this is maximum exponent and
1825 top two significant bits set. */
1829 unsigned EMUSHORT x[];
1837 /* Return nonzero if exploded e-type X is a NaN. */
1841 unsigned EMUSHORT x[];
1845 if ((x[E] & 0x7fff) == 0x7fff)
1847 for (i = M + 1; i < NI; i++)
1856 /* Return nonzero if sign of exploded e-type X is nonzero. */
1860 unsigned EMUSHORT x[];
1867 /* Fill exploded e-type X with infinity pattern.
1868 This has maximum exponent and significand all zeros. */
1872 unsigned EMUSHORT x[];
1880 /* Return nonzero if exploded e-type X is infinite. */
1884 unsigned EMUSHORT x[];
1891 if ((x[E] & 0x7fff) == 0x7fff)
1897 /* Compare significands of numbers in internal exploded e-type format.
1898 Guard words are included in the comparison.
1906 register unsigned EMUSHORT *a, *b;
1910 a += M; /* skip up to significand area */
1912 for (i = M; i < NI; i++)
1920 if (*(--a) > *(--b))
1926 /* Shift significand of exploded e-type X down by 1 bit. */
1930 register unsigned EMUSHORT *x;
1932 register unsigned EMUSHORT bits;
1935 x += M; /* point to significand area */
1938 for (i = M; i < NI; i++)
1950 /* Shift significand of exploded e-type X up by 1 bit. */
1954 register unsigned EMUSHORT *x;
1956 register unsigned EMUSHORT bits;
1962 for (i = M; i < NI; i++)
1975 /* Shift significand of exploded e-type X down by 8 bits. */
1979 register unsigned EMUSHORT *x;
1981 register unsigned EMUSHORT newbyt, oldbyt;
1986 for (i = M; i < NI; i++)
1996 /* Shift significand of exploded e-type X up by 8 bits. */
2000 register unsigned EMUSHORT *x;
2003 register unsigned EMUSHORT newbyt, oldbyt;
2008 for (i = M; i < NI; i++)
2018 /* Shift significand of exploded e-type X up by 16 bits. */
2022 register unsigned EMUSHORT *x;
2025 register unsigned EMUSHORT *p;
2030 for (i = M; i < NI - 1; i++)
2036 /* Shift significand of exploded e-type X down by 16 bits. */
2040 register unsigned EMUSHORT *x;
2043 register unsigned EMUSHORT *p;
2048 for (i = M; i < NI - 1; i++)
2054 /* Add significands of exploded e-type X and Y. X + Y replaces Y. */
2058 unsigned EMUSHORT *x, *y;
2060 register unsigned EMULONG a;
2067 for (i = M; i < NI; i++)
2069 a = (unsigned EMULONG) (*x) + (unsigned EMULONG) (*y) + carry;
2074 *y = (unsigned EMUSHORT) a;
2080 /* Subtract significands of exploded e-type X and Y. Y - X replaces Y. */
2084 unsigned EMUSHORT *x, *y;
2093 for (i = M; i < NI; i++)
2095 a = (unsigned EMULONG) (*y) - (unsigned EMULONG) (*x) - carry;
2100 *y = (unsigned EMUSHORT) a;
2107 static unsigned EMUSHORT equot[NI];
2111 /* Radix 2 shift-and-add versions of multiply and divide */
2114 /* Divide significands */
2118 unsigned EMUSHORT den[], num[];
2121 register unsigned EMUSHORT *p, *q;
2122 unsigned EMUSHORT j;
2128 for (i = M; i < NI; i++)
2133 /* Use faster compare and subtraction if denominator has only 15 bits of
2139 for (i = M + 3; i < NI; i++)
2144 if ((den[M + 1] & 1) != 0)
2152 for (i = 0; i < NBITS + 2; i++)
2170 /* The number of quotient bits to calculate is NBITS + 1 scaling guard
2171 bit + 1 roundoff bit. */
2176 for (i = 0; i < NBITS + 2; i++)
2178 if (ecmpm (den, num) <= 0)
2181 j = 1; /* quotient bit = 1 */
2195 /* test for nonzero remainder after roundoff bit */
2198 for (i = M; i < NI; i++)
2206 for (i = 0; i < NI; i++)
2212 /* Multiply significands */
2216 unsigned EMUSHORT a[], b[];
2218 unsigned EMUSHORT *p, *q;
2223 for (i = M; i < NI; i++)
2228 while (*p == 0) /* significand is not supposed to be zero */
2233 if ((*p & 0xff) == 0)
2241 for (i = 0; i < k; i++)
2245 /* remember if there were any nonzero bits shifted out */
2252 for (i = 0; i < NI; i++)
2255 /* return flag for lost nonzero bits */
2261 /* Radix 65536 versions of multiply and divide. */
2263 /* Multiply significand of e-type number B
2264 by 16-bit quantity A, return e-type result to C. */
2269 unsigned EMUSHORT b[], c[];
2271 register unsigned EMUSHORT *pp;
2272 register unsigned EMULONG carry;
2273 unsigned EMUSHORT *ps;
2274 unsigned EMUSHORT p[NI];
2275 unsigned EMULONG aa, m;
2284 for (i=M+1; i<NI; i++)
2294 m = (unsigned EMULONG) aa * *ps--;
2295 carry = (m & 0xffff) + *pp;
2296 *pp-- = (unsigned EMUSHORT)carry;
2297 carry = (carry >> 16) + (m >> 16) + *pp;
2298 *pp = (unsigned EMUSHORT)carry;
2299 *(pp-1) = carry >> 16;
2302 for (i=M; i<NI; i++)
2306 /* Divide significands of exploded e-types NUM / DEN. Neither the
2307 numerator NUM nor the denominator DEN is permitted to have its high guard
2312 unsigned EMUSHORT den[], num[];
2315 register unsigned EMUSHORT *p;
2316 unsigned EMULONG tnum;
2317 unsigned EMUSHORT j, tdenm, tquot;
2318 unsigned EMUSHORT tprod[NI+1];
2324 for (i=M; i<NI; i++)
2330 for (i=M; i<NI; i++)
2332 /* Find trial quotient digit (the radix is 65536). */
2333 tnum = (((unsigned EMULONG) num[M]) << 16) + num[M+1];
2335 /* Do not execute the divide instruction if it will overflow. */
2336 if ((tdenm * 0xffffL) < tnum)
2339 tquot = tnum / tdenm;
2340 /* Multiply denominator by trial quotient digit. */
2341 m16m ((unsigned int)tquot, den, tprod);
2342 /* The quotient digit may have been overestimated. */
2343 if (ecmpm (tprod, num) > 0)
2347 if (ecmpm (tprod, num) > 0)
2357 /* test for nonzero remainder after roundoff bit */
2360 for (i=M; i<NI; i++)
2367 for (i=0; i<NI; i++)
2373 /* Multiply significands of exploded e-type A and B, result in B. */
2377 unsigned EMUSHORT a[], b[];
2379 unsigned EMUSHORT *p, *q;
2380 unsigned EMUSHORT pprod[NI];
2381 unsigned EMUSHORT j;
2386 for (i=M; i<NI; i++)
2392 for (i=M+1; i<NI; i++)
2400 m16m ((unsigned int) *p--, b, pprod);
2401 eaddm(pprod, equot);
2407 for (i=0; i<NI; i++)
2410 /* return flag for lost nonzero bits */
2416 /* Normalize and round off.
2418 The internal format number to be rounded is S.
2419 Input LOST is 0 if the value is exact. This is the so-called sticky bit.
2421 Input SUBFLG indicates whether the number was obtained
2422 by a subtraction operation. In that case if LOST is nonzero
2423 then the number is slightly smaller than indicated.
2425 Input EXP is the biased exponent, which may be negative.
2426 the exponent field of S is ignored but is replaced by
2427 EXP as adjusted by normalization and rounding.
2429 Input RCNTRL is the rounding control. If it is nonzero, the
2430 returned value will be rounded to RNDPRC bits.
2432 For future reference: In order for emdnorm to round off denormal
2433 significands at the right point, the input exponent must be
2434 adjusted to be the actual value it would have after conversion to
2435 the final floating point type. This adjustment has been
2436 implemented for all type conversions (etoe53, etc.) and decimal
2437 conversions, but not for the arithmetic functions (eadd, etc.).
2438 Data types having standard 15-bit exponents are not affected by
2439 this, but SFmode and DFmode are affected. For example, ediv with
2440 rndprc = 24 will not round correctly to 24-bit precision if the
2441 result is denormal. */
2443 static int rlast = -1;
2445 static unsigned EMUSHORT rmsk = 0;
2446 static unsigned EMUSHORT rmbit = 0;
2447 static unsigned EMUSHORT rebit = 0;
2449 static unsigned EMUSHORT rbit[NI];
2452 emdnorm (s, lost, subflg, exp, rcntrl)
2453 unsigned EMUSHORT s[];
2460 unsigned EMUSHORT r;
2465 /* a blank significand could mean either zero or infinity. */
2478 if ((j > NBITS) && (exp < 32767))
2486 if (exp > (EMULONG) (-NBITS - 1))
2499 /* Round off, unless told not to by rcntrl. */
2502 /* Set up rounding parameters if the control register changed. */
2503 if (rndprc != rlast)
2510 rw = NI - 1; /* low guard word */
2533 /* For DEC or IBM arithmetic */
2550 /* For C4x arithmetic */
2571 /* Shift down 1 temporarily if the data structure has an implied
2572 most significant bit and the number is denormal.
2573 Intel long double denormals also lose one bit of precision. */
2574 if ((exp <= 0) && (rndprc != NBITS)
2575 && ((rndprc != 64) || ((rndprc == 64) && ! REAL_WORDS_BIG_ENDIAN)))
2577 lost |= s[NI - 1] & 1;
2580 /* Clear out all bits below the rounding bit,
2581 remembering in r if any were nonzero. */
2595 if ((r & rmbit) != 0)
2600 { /* round to even */
2601 if ((s[re] & rebit) == 0)
2613 /* Undo the temporary shift for denormal values. */
2614 if ((exp <= 0) && (rndprc != NBITS)
2615 && ((rndprc != 64) || ((rndprc == 64) && ! REAL_WORDS_BIG_ENDIAN)))
2620 { /* overflow on roundoff */
2633 for (i = 2; i < NI - 1; i++)
2636 warning ("floating point overflow");
2640 for (i = M + 1; i < NI - 1; i++)
2643 if ((rndprc < 64) || (rndprc == 113))
2658 s[1] = (unsigned EMUSHORT) exp;
2661 /* Subtract. C = B - A, all e type numbers. */
2663 static int subflg = 0;
2667 unsigned EMUSHORT *a, *b, *c;
2681 /* Infinity minus infinity is a NaN.
2682 Test for subtracting infinities of the same sign. */
2683 if (eisinf (a) && eisinf (b)
2684 && ((eisneg (a) ^ eisneg (b)) == 0))
2686 mtherr ("esub", INVALID);
2695 /* Add. C = A + B, all e type. */
2699 unsigned EMUSHORT *a, *b, *c;
2703 /* NaN plus anything is a NaN. */
2714 /* Infinity minus infinity is a NaN.
2715 Test for adding infinities of opposite signs. */
2716 if (eisinf (a) && eisinf (b)
2717 && ((eisneg (a) ^ eisneg (b)) != 0))
2719 mtherr ("esub", INVALID);
2728 /* Arithmetic common to both addition and subtraction. */
2732 unsigned EMUSHORT *a, *b, *c;
2734 unsigned EMUSHORT ai[NI], bi[NI], ci[NI];
2736 EMULONG lt, lta, ltb;
2757 /* compare exponents */
2762 { /* put the larger number in bi */
2772 if (lt < (EMULONG) (-NBITS - 1))
2773 goto done; /* answer same as larger addend */
2775 lost = eshift (ai, k); /* shift the smaller number down */
2779 /* exponents were the same, so must compare significands */
2782 { /* the numbers are identical in magnitude */
2783 /* if different signs, result is zero */
2789 /* if same sign, result is double */
2790 /* double denormalized tiny number */
2791 if ((bi[E] == 0) && ((bi[3] & 0x8000) == 0))
2796 /* add 1 to exponent unless both are zero! */
2797 for (j = 1; j < NI - 1; j++)
2813 bi[E] = (unsigned EMUSHORT) ltb;
2817 { /* put the larger number in bi */
2833 emdnorm (bi, lost, subflg, ltb, 64);
2839 /* Divide: C = B/A, all e type. */
2843 unsigned EMUSHORT *a, *b, *c;
2845 unsigned EMUSHORT ai[NI], bi[NI];
2847 EMULONG lt, lta, ltb;
2849 /* IEEE says if result is not a NaN, the sign is "-" if and only if
2850 operands have opposite signs -- but flush -0 to 0 later if not IEEE. */
2851 sign = eisneg(a) ^ eisneg(b);
2854 /* Return any NaN input. */
2865 /* Zero over zero, or infinity over infinity, is a NaN. */
2866 if (((ecmp (a, ezero) == 0) && (ecmp (b, ezero) == 0))
2867 || (eisinf (a) && eisinf (b)))
2869 mtherr ("ediv", INVALID);
2874 /* Infinity over anything else is infinity. */
2881 /* Anything else over infinity is zero. */
2893 { /* See if numerator is zero. */
2894 for (i = 1; i < NI - 1; i++)
2898 ltb -= enormlz (bi);
2908 { /* possible divide by zero */
2909 for (i = 1; i < NI - 1; i++)
2913 lta -= enormlz (ai);
2917 /* Divide by zero is not an invalid operation.
2918 It is a divide-by-zero operation! */
2920 mtherr ("ediv", SING);
2926 /* calculate exponent */
2927 lt = ltb - lta + EXONE;
2928 emdnorm (bi, i, 0, lt, 64);
2935 && (ecmp (c, ezero) != 0)
2938 *(c+(NE-1)) |= 0x8000;
2940 *(c+(NE-1)) &= ~0x8000;
2943 /* Multiply e-types A and B, return e-type product C. */
2947 unsigned EMUSHORT *a, *b, *c;
2949 unsigned EMUSHORT ai[NI], bi[NI];
2951 EMULONG lt, lta, ltb;
2953 /* IEEE says if result is not a NaN, the sign is "-" if and only if
2954 operands have opposite signs -- but flush -0 to 0 later if not IEEE. */
2955 sign = eisneg(a) ^ eisneg(b);
2958 /* NaN times anything is the same NaN. */
2969 /* Zero times infinity is a NaN. */
2970 if ((eisinf (a) && (ecmp (b, ezero) == 0))
2971 || (eisinf (b) && (ecmp (a, ezero) == 0)))
2973 mtherr ("emul", INVALID);
2978 /* Infinity times anything else is infinity. */
2980 if (eisinf (a) || eisinf (b))
2992 for (i = 1; i < NI - 1; i++)
2996 lta -= enormlz (ai);
3007 for (i = 1; i < NI - 1; i++)
3011 ltb -= enormlz (bi);
3020 /* Multiply significands */
3022 /* calculate exponent */
3023 lt = lta + ltb - (EXONE - 1);
3024 emdnorm (bi, j, 0, lt, 64);
3031 && (ecmp (c, ezero) != 0)
3034 *(c+(NE-1)) |= 0x8000;
3036 *(c+(NE-1)) &= ~0x8000;
3039 /* Convert double precision PE to e-type Y. */
3043 unsigned EMUSHORT *pe, *y;
3052 ibmtoe (pe, y, DFmode);
3057 c4xtoe (pe, y, HFmode);
3060 register unsigned EMUSHORT r;
3061 register unsigned EMUSHORT *e, *p;
3062 unsigned EMUSHORT yy[NI];
3066 denorm = 0; /* flag if denormalized number */
3068 if (! REAL_WORDS_BIG_ENDIAN)
3074 yy[M] = (r & 0x0f) | 0x10;
3075 r &= ~0x800f; /* strip sign and 4 significand bits */
3080 if (! REAL_WORDS_BIG_ENDIAN)
3082 if (((pe[3] & 0xf) != 0) || (pe[2] != 0)
3083 || (pe[1] != 0) || (pe[0] != 0))
3085 enan (y, yy[0] != 0);
3091 if (((pe[0] & 0xf) != 0) || (pe[1] != 0)
3092 || (pe[2] != 0) || (pe[3] != 0))
3094 enan (y, yy[0] != 0);
3105 #endif /* INFINITY */
3107 /* If zero exponent, then the significand is denormalized.
3108 So take back the understood high significand bit. */
3119 if (! REAL_WORDS_BIG_ENDIAN)
3136 /* If zero exponent, then normalize the significand. */
3137 if ((k = enormlz (yy)) > NBITS)
3140 yy[E] -= (unsigned EMUSHORT) (k - 1);
3143 #endif /* not C4X */
3144 #endif /* not IBM */
3145 #endif /* not DEC */
3148 /* Convert double extended precision float PE to e type Y. */
3152 unsigned EMUSHORT *pe, *y;
3154 unsigned EMUSHORT yy[NI];
3155 unsigned EMUSHORT *e, *p, *q;
3160 for (i = 0; i < NE - 5; i++)
3162 /* This precision is not ordinarily supported on DEC or IBM. */
3164 for (i = 0; i < 5; i++)
3168 p = &yy[0] + (NE - 1);
3171 for (i = 0; i < 5; i++)
3175 if (! REAL_WORDS_BIG_ENDIAN)
3177 for (i = 0; i < 5; i++)
3180 /* For denormal long double Intel format, shift significand up one
3181 -- but only if the top significand bit is zero. A top bit of 1
3182 is "pseudodenormal" when the exponent is zero. */
3183 if((yy[NE-1] & 0x7fff) == 0 && (yy[NE-2] & 0x8000) == 0)
3185 unsigned EMUSHORT temp[NI];
3195 p = &yy[0] + (NE - 1);
3196 #ifdef ARM_EXTENDED_IEEE_FORMAT
3197 /* For ARMs, the exponent is in the lowest 15 bits of the word. */
3198 *p-- = (e[0] & 0x8000) | (e[1] & 0x7ffff);
3204 for (i = 0; i < 4; i++)
3209 /* Point to the exponent field and check max exponent cases. */
3211 if ((*p & 0x7fff) == 0x7fff)
3214 if (! REAL_WORDS_BIG_ENDIAN)
3216 for (i = 0; i < 4; i++)
3218 if ((i != 3 && pe[i] != 0)
3219 /* Anything but 0x8000 here, including 0, is a NaN. */
3220 || (i == 3 && pe[i] != 0x8000))
3222 enan (y, (*p & 0x8000) != 0);
3229 #ifdef ARM_EXTENDED_IEEE_FORMAT
3230 for (i = 2; i <= 5; i++)
3234 enan (y, (*p & 0x8000) != 0);
3239 /* In Motorola extended precision format, the most significant
3240 bit of an infinity mantissa could be either 1 or 0. It is
3241 the lower order bits that tell whether the value is a NaN. */
3242 if ((pe[2] & 0x7fff) != 0)
3245 for (i = 3; i <= 5; i++)
3250 enan (y, (*p & 0x8000) != 0);
3254 #endif /* not ARM */
3263 #endif /* INFINITY */
3266 for (i = 0; i < NE; i++)
3270 /* Convert 128-bit long double precision float PE to e type Y. */
3274 unsigned EMUSHORT *pe, *y;
3276 register unsigned EMUSHORT r;
3277 unsigned EMUSHORT *e, *p;
3278 unsigned EMUSHORT yy[NI];
3285 if (! REAL_WORDS_BIG_ENDIAN)
3297 if (! REAL_WORDS_BIG_ENDIAN)
3299 for (i = 0; i < 7; i++)
3303 enan (y, yy[0] != 0);
3310 for (i = 1; i < 8; i++)
3314 enan (y, yy[0] != 0);
3326 #endif /* INFINITY */
3330 if (! REAL_WORDS_BIG_ENDIAN)
3332 for (i = 0; i < 7; i++)
3338 for (i = 0; i < 7; i++)
3342 /* If denormal, remove the implied bit; else shift down 1. */
3355 /* Convert single precision float PE to e type Y. */
3359 unsigned EMUSHORT *pe, *y;
3363 ibmtoe (pe, y, SFmode);
3369 c4xtoe (pe, y, QFmode);
3373 register unsigned EMUSHORT r;
3374 register unsigned EMUSHORT *e, *p;
3375 unsigned EMUSHORT yy[NI];
3379 denorm = 0; /* flag if denormalized number */
3382 if (! REAL_WORDS_BIG_ENDIAN)
3392 yy[M] = (r & 0x7f) | 0200;
3393 r &= ~0x807f; /* strip sign and 7 significand bits */
3398 if (REAL_WORDS_BIG_ENDIAN)
3400 if (((pe[0] & 0x7f) != 0) || (pe[1] != 0))
3402 enan (y, yy[0] != 0);
3408 if (((pe[1] & 0x7f) != 0) || (pe[0] != 0))
3410 enan (y, yy[0] != 0);
3421 #endif /* INFINITY */
3423 /* If zero exponent, then the significand is denormalized.
3424 So take back the understood high significand bit. */
3437 if (! REAL_WORDS_BIG_ENDIAN)
3447 { /* if zero exponent, then normalize the significand */
3448 if ((k = enormlz (yy)) > NBITS)
3451 yy[E] -= (unsigned EMUSHORT) (k - 1);
3454 #endif /* not C4X */
3455 #endif /* not IBM */
3458 /* Convert e-type X to IEEE 128-bit long double format E. */
3462 unsigned EMUSHORT *x, *e;
3464 unsigned EMUSHORT xi[NI];
3471 make_nan (e, eisneg (x), TFmode);
3476 exp = (EMULONG) xi[E];
3481 /* round off to nearest or even */
3484 emdnorm (xi, 0, 0, exp, 64);
3490 /* Convert exploded e-type X, that has already been rounded to
3491 113-bit precision, to IEEE 128-bit long double format Y. */
3495 unsigned EMUSHORT *a, *b;
3497 register unsigned EMUSHORT *p, *q;
3498 unsigned EMUSHORT i;
3503 make_nan (b, eiisneg (a), TFmode);
3508 if (REAL_WORDS_BIG_ENDIAN)
3511 q = b + 7; /* point to output exponent */
3513 /* If not denormal, delete the implied bit. */
3518 /* combine sign and exponent */
3520 if (REAL_WORDS_BIG_ENDIAN)
3523 *q++ = *p++ | 0x8000;
3530 *q-- = *p++ | 0x8000;
3534 /* skip over guard word */
3536 /* move the significand */
3537 if (REAL_WORDS_BIG_ENDIAN)
3539 for (i = 0; i < 7; i++)
3544 for (i = 0; i < 7; i++)
3549 /* Convert e-type X to IEEE double extended format E. */
3553 unsigned EMUSHORT *x, *e;
3555 unsigned EMUSHORT xi[NI];
3562 make_nan (e, eisneg (x), XFmode);
3567 /* adjust exponent for offset */
3568 exp = (EMULONG) xi[E];
3573 /* round off to nearest or even */
3576 emdnorm (xi, 0, 0, exp, 64);
3582 /* Convert exploded e-type X, that has already been rounded to
3583 64-bit precision, to IEEE double extended format Y. */
3587 unsigned EMUSHORT *a, *b;
3589 register unsigned EMUSHORT *p, *q;
3590 unsigned EMUSHORT i;
3595 make_nan (b, eiisneg (a), XFmode);
3599 /* Shift denormal long double Intel format significand down one bit. */
3600 if ((a[E] == 0) && ! REAL_WORDS_BIG_ENDIAN)
3610 if (REAL_WORDS_BIG_ENDIAN)
3614 q = b + 4; /* point to output exponent */
3615 #if LONG_DOUBLE_TYPE_SIZE == 96
3616 /* Clear the last two bytes of 12-byte Intel format */
3622 /* combine sign and exponent */
3626 *q++ = *p++ | 0x8000;
3633 *q-- = *p++ | 0x8000;
3638 if (REAL_WORDS_BIG_ENDIAN)
3640 #ifdef ARM_EXTENDED_IEEE_FORMAT
3641 /* The exponent is in the lowest 15 bits of the first word. */
3642 *q++ = i ? 0x8000 : 0;
3646 *q++ = *p++ | 0x8000;
3655 *q-- = *p++ | 0x8000;
3660 /* skip over guard word */
3662 /* move the significand */
3664 for (i = 0; i < 4; i++)
3668 for (i = 0; i < 4; i++)
3672 if (REAL_WORDS_BIG_ENDIAN)
3674 for (i = 0; i < 4; i++)
3682 /* Intel long double infinity significand. */
3690 for (i = 0; i < 4; i++)
3696 /* e type to double precision. */
3699 /* Convert e-type X to DEC-format double E. */
3703 unsigned EMUSHORT *x, *e;
3705 etodec (x, e); /* see etodec.c */
3708 /* Convert exploded e-type X, that has already been rounded to
3709 56-bit double precision, to DEC double Y. */
3713 unsigned EMUSHORT *x, *y;
3720 /* Convert e-type X to IBM 370-format double E. */
3724 unsigned EMUSHORT *x, *e;
3726 etoibm (x, e, DFmode);
3729 /* Convert exploded e-type X, that has already been rounded to
3730 56-bit precision, to IBM 370 double Y. */
3734 unsigned EMUSHORT *x, *y;
3736 toibm (x, y, DFmode);
3739 #else /* it's neither DEC nor IBM */
3741 /* Convert e-type X to C4X-format long double E. */
3745 unsigned EMUSHORT *x, *e;
3747 etoc4x (x, e, HFmode);
3750 /* Convert exploded e-type X, that has already been rounded to
3751 56-bit precision, to IBM 370 double Y. */
3755 unsigned EMUSHORT *x, *y;
3757 toc4x (x, y, HFmode);
3760 #else /* it's neither DEC nor IBM nor C4X */
3762 /* Convert e-type X to IEEE double E. */
3766 unsigned EMUSHORT *x, *e;
3768 unsigned EMUSHORT xi[NI];
3775 make_nan (e, eisneg (x), DFmode);
3780 /* adjust exponent for offsets */
3781 exp = (EMULONG) xi[E] - (EXONE - 0x3ff);
3786 /* round off to nearest or even */
3789 emdnorm (xi, 0, 0, exp, 64);
3795 /* Convert exploded e-type X, that has already been rounded to
3796 53-bit precision, to IEEE double Y. */
3800 unsigned EMUSHORT *x, *y;
3802 unsigned EMUSHORT i;
3803 unsigned EMUSHORT *p;
3808 make_nan (y, eiisneg (x), DFmode);
3814 if (! REAL_WORDS_BIG_ENDIAN)
3817 *y = 0; /* output high order */
3819 *y = 0x8000; /* output sign bit */
3822 if (i >= (unsigned int) 2047)
3824 /* Saturate at largest number less than infinity. */
3827 if (! REAL_WORDS_BIG_ENDIAN)
3841 *y |= (unsigned EMUSHORT) 0x7fef;
3842 if (! REAL_WORDS_BIG_ENDIAN)
3867 i |= *p++ & (unsigned EMUSHORT) 0x0f; /* *p = xi[M] */
3868 *y |= (unsigned EMUSHORT) i; /* high order output already has sign bit set */
3869 if (! REAL_WORDS_BIG_ENDIAN)
3884 #endif /* not C4X */
3885 #endif /* not IBM */
3886 #endif /* not DEC */
3890 /* e type to single precision. */
3893 /* Convert e-type X to IBM 370 float E. */
3897 unsigned EMUSHORT *x, *e;
3899 etoibm (x, e, SFmode);
3902 /* Convert exploded e-type X, that has already been rounded to
3903 float precision, to IBM 370 float Y. */
3907 unsigned EMUSHORT *x, *y;
3909 toibm (x, y, SFmode);
3915 /* Convert e-type X to C4X float E. */
3919 unsigned EMUSHORT *x, *e;
3921 etoc4x (x, e, QFmode);
3924 /* Convert exploded e-type X, that has already been rounded to
3925 float precision, to IBM 370 float Y. */
3929 unsigned EMUSHORT *x, *y;
3931 toc4x (x, y, QFmode);
3936 /* Convert e-type X to IEEE float E. DEC float is the same as IEEE float. */
3940 unsigned EMUSHORT *x, *e;
3943 unsigned EMUSHORT xi[NI];
3949 make_nan (e, eisneg (x), SFmode);
3954 /* adjust exponent for offsets */
3955 exp = (EMULONG) xi[E] - (EXONE - 0177);
3960 /* round off to nearest or even */
3963 emdnorm (xi, 0, 0, exp, 64);
3969 /* Convert exploded e-type X, that has already been rounded to
3970 float precision, to IEEE float Y. */
3974 unsigned EMUSHORT *x, *y;
3976 unsigned EMUSHORT i;
3977 unsigned EMUSHORT *p;
3982 make_nan (y, eiisneg (x), SFmode);
3988 if (! REAL_WORDS_BIG_ENDIAN)
3994 *y = 0; /* output high order */
3996 *y = 0x8000; /* output sign bit */
3999 /* Handle overflow cases. */
4003 *y |= (unsigned EMUSHORT) 0x7f80;
4008 if (! REAL_WORDS_BIG_ENDIAN)
4016 #else /* no INFINITY */
4017 *y |= (unsigned EMUSHORT) 0x7f7f;
4022 if (! REAL_WORDS_BIG_ENDIAN)
4033 #endif /* no INFINITY */
4045 i |= *p++ & (unsigned EMUSHORT) 0x7f; /* *p = xi[M] */
4046 /* High order output already has sign bit set. */
4052 if (! REAL_WORDS_BIG_ENDIAN)
4061 #endif /* not C4X */
4062 #endif /* not IBM */
4064 /* Compare two e type numbers.
4068 -2 if either a or b is a NaN. */
4072 unsigned EMUSHORT *a, *b;
4074 unsigned EMUSHORT ai[NI], bi[NI];
4075 register unsigned EMUSHORT *p, *q;
4080 if (eisnan (a) || eisnan (b))
4089 { /* the signs are different */
4091 for (i = 1; i < NI - 1; i++)
4105 /* both are the same sign */
4120 return (0); /* equality */
4124 if (*(--p) > *(--q))
4125 return (msign); /* p is bigger */
4127 return (-msign); /* p is littler */
4131 /* Find e-type nearest integer to X, as floor (X + 0.5). */
4135 unsigned EMUSHORT *x, *y;
4142 /* Convert HOST_WIDE_INT LP to e type Y. */
4147 unsigned EMUSHORT *y;
4149 unsigned EMUSHORT yi[NI];
4150 unsigned HOST_WIDE_INT ll;
4156 /* make it positive */
4157 ll = (unsigned HOST_WIDE_INT) (-(*lp));
4158 yi[0] = 0xffff; /* put correct sign in the e type number */
4162 ll = (unsigned HOST_WIDE_INT) (*lp);
4164 /* move the long integer to yi significand area */
4165 #if HOST_BITS_PER_WIDE_INT == 64
4166 yi[M] = (unsigned EMUSHORT) (ll >> 48);
4167 yi[M + 1] = (unsigned EMUSHORT) (ll >> 32);
4168 yi[M + 2] = (unsigned EMUSHORT) (ll >> 16);
4169 yi[M + 3] = (unsigned EMUSHORT) ll;
4170 yi[E] = EXONE + 47; /* exponent if normalize shift count were 0 */
4172 yi[M] = (unsigned EMUSHORT) (ll >> 16);
4173 yi[M + 1] = (unsigned EMUSHORT) ll;
4174 yi[E] = EXONE + 15; /* exponent if normalize shift count were 0 */
4177 if ((k = enormlz (yi)) > NBITS)/* normalize the significand */
4178 ecleaz (yi); /* it was zero */
4180 yi[E] -= (unsigned EMUSHORT) k;/* subtract shift count from exponent */
4181 emovo (yi, y); /* output the answer */
4184 /* Convert unsigned HOST_WIDE_INT LP to e type Y. */
4188 unsigned HOST_WIDE_INT *lp;
4189 unsigned EMUSHORT *y;
4191 unsigned EMUSHORT yi[NI];
4192 unsigned HOST_WIDE_INT ll;
4198 /* move the long integer to ayi significand area */
4199 #if HOST_BITS_PER_WIDE_INT == 64
4200 yi[M] = (unsigned EMUSHORT) (ll >> 48);
4201 yi[M + 1] = (unsigned EMUSHORT) (ll >> 32);
4202 yi[M + 2] = (unsigned EMUSHORT) (ll >> 16);
4203 yi[M + 3] = (unsigned EMUSHORT) ll;
4204 yi[E] = EXONE + 47; /* exponent if normalize shift count were 0 */
4206 yi[M] = (unsigned EMUSHORT) (ll >> 16);
4207 yi[M + 1] = (unsigned EMUSHORT) ll;
4208 yi[E] = EXONE + 15; /* exponent if normalize shift count were 0 */
4211 if ((k = enormlz (yi)) > NBITS)/* normalize the significand */
4212 ecleaz (yi); /* it was zero */
4214 yi[E] -= (unsigned EMUSHORT) k; /* subtract shift count from exponent */
4215 emovo (yi, y); /* output the answer */
4219 /* Find signed HOST_WIDE_INT integer I and floating point fractional
4220 part FRAC of e-type (packed internal format) floating point input X.
4221 The integer output I has the sign of the input, except that
4222 positive overflow is permitted if FIXUNS_TRUNC_LIKE_FIX_TRUNC.
4223 The output e-type fraction FRAC is the positive fractional
4228 unsigned EMUSHORT *x;
4230 unsigned EMUSHORT *frac;
4232 unsigned EMUSHORT xi[NI];
4234 unsigned HOST_WIDE_INT ll;
4237 k = (int) xi[E] - (EXONE - 1);
4240 /* if exponent <= 0, integer = 0 and real output is fraction */
4245 if (k > (HOST_BITS_PER_WIDE_INT - 1))
4247 /* long integer overflow: output large integer
4248 and correct fraction */
4250 *i = ((unsigned HOST_WIDE_INT) 1) << (HOST_BITS_PER_WIDE_INT - 1);
4253 #ifdef FIXUNS_TRUNC_LIKE_FIX_TRUNC
4254 /* In this case, let it overflow and convert as if unsigned. */
4255 euifrac (x, &ll, frac);
4256 *i = (HOST_WIDE_INT) ll;
4259 /* In other cases, return the largest positive integer. */
4260 *i = (((unsigned HOST_WIDE_INT) 1) << (HOST_BITS_PER_WIDE_INT - 1)) - 1;
4265 warning ("overflow on truncation to integer");
4269 /* Shift more than 16 bits: first shift up k-16 mod 16,
4270 then shift up by 16's. */
4271 j = k - ((k >> 4) << 4);
4278 ll = (ll << 16) | xi[M];
4280 while ((k -= 16) > 0);
4287 /* shift not more than 16 bits */
4289 *i = (HOST_WIDE_INT) xi[M] & 0xffff;
4296 if ((k = enormlz (xi)) > NBITS)
4299 xi[E] -= (unsigned EMUSHORT) k;
4305 /* Find unsigned HOST_WIDE_INT integer I and floating point fractional part
4306 FRAC of e-type X. A negative input yields integer output = 0 but
4307 correct fraction. */
4310 euifrac (x, i, frac)
4311 unsigned EMUSHORT *x;
4312 unsigned HOST_WIDE_INT *i;
4313 unsigned EMUSHORT *frac;
4315 unsigned HOST_WIDE_INT ll;
4316 unsigned EMUSHORT xi[NI];
4320 k = (int) xi[E] - (EXONE - 1);
4323 /* if exponent <= 0, integer = 0 and argument is fraction */
4328 if (k > HOST_BITS_PER_WIDE_INT)
4330 /* Long integer overflow: output large integer
4331 and correct fraction.
4332 Note, the BSD microvax compiler says that ~(0UL)
4333 is a syntax error. */
4337 warning ("overflow on truncation to unsigned integer");
4341 /* Shift more than 16 bits: first shift up k-16 mod 16,
4342 then shift up by 16's. */
4343 j = k - ((k >> 4) << 4);
4350 ll = (ll << 16) | xi[M];
4352 while ((k -= 16) > 0);
4357 /* shift not more than 16 bits */
4359 *i = (HOST_WIDE_INT) xi[M] & 0xffff;
4362 if (xi[0]) /* A negative value yields unsigned integer 0. */
4368 if ((k = enormlz (xi)) > NBITS)
4371 xi[E] -= (unsigned EMUSHORT) k;
4376 /* Shift the significand of exploded e-type X up or down by SC bits. */
4380 unsigned EMUSHORT *x;
4383 unsigned EMUSHORT lost;
4384 unsigned EMUSHORT *p;
4397 lost |= *p; /* remember lost bits */
4438 return ((int) lost);
4441 /* Shift normalize the significand area of exploded e-type X.
4442 Return the shift count (up = positive). */
4446 unsigned EMUSHORT x[];
4448 register unsigned EMUSHORT *p;
4457 return (0); /* already normalized */
4463 /* With guard word, there are NBITS+16 bits available.
4464 Return true if all are zero. */
4468 /* see if high byte is zero */
4469 while ((*p & 0xff00) == 0)
4474 /* now shift 1 bit at a time */
4475 while ((*p & 0x8000) == 0)
4481 mtherr ("enormlz", UNDERFLOW);
4487 /* Normalize by shifting down out of the high guard word
4488 of the significand */
4503 mtherr ("enormlz", OVERFLOW);
4510 /* Powers of ten used in decimal <-> binary conversions. */
4515 #if LONG_DOUBLE_TYPE_SIZE == 128
4516 static unsigned EMUSHORT etens[NTEN + 1][NE] =
4518 {0x6576, 0x4a92, 0x804a, 0x153f,
4519 0xc94c, 0x979a, 0x8a20, 0x5202, 0xc460, 0x7525,}, /* 10**4096 */
4520 {0x6a32, 0xce52, 0x329a, 0x28ce,
4521 0xa74d, 0x5de4, 0xc53d, 0x3b5d, 0x9e8b, 0x5a92,}, /* 10**2048 */
4522 {0x526c, 0x50ce, 0xf18b, 0x3d28,
4523 0x650d, 0x0c17, 0x8175, 0x7586, 0xc976, 0x4d48,},
4524 {0x9c66, 0x58f8, 0xbc50, 0x5c54,
4525 0xcc65, 0x91c6, 0xa60e, 0xa0ae, 0xe319, 0x46a3,},
4526 {0x851e, 0xeab7, 0x98fe, 0x901b,
4527 0xddbb, 0xde8d, 0x9df9, 0xebfb, 0xaa7e, 0x4351,},
4528 {0x0235, 0x0137, 0x36b1, 0x336c,
4529 0xc66f, 0x8cdf, 0x80e9, 0x47c9, 0x93ba, 0x41a8,},
4530 {0x50f8, 0x25fb, 0xc76b, 0x6b71,
4531 0x3cbf, 0xa6d5, 0xffcf, 0x1f49, 0xc278, 0x40d3,},
4532 {0x0000, 0x0000, 0x0000, 0x0000,
4533 0xf020, 0xb59d, 0x2b70, 0xada8, 0x9dc5, 0x4069,},
4534 {0x0000, 0x0000, 0x0000, 0x0000,
4535 0x0000, 0x0000, 0x0400, 0xc9bf, 0x8e1b, 0x4034,},
4536 {0x0000, 0x0000, 0x0000, 0x0000,
4537 0x0000, 0x0000, 0x0000, 0x2000, 0xbebc, 0x4019,},
4538 {0x0000, 0x0000, 0x0000, 0x0000,
4539 0x0000, 0x0000, 0x0000, 0x0000, 0x9c40, 0x400c,},
4540 {0x0000, 0x0000, 0x0000, 0x0000,
4541 0x0000, 0x0000, 0x0000, 0x0000, 0xc800, 0x4005,},
4542 {0x0000, 0x0000, 0x0000, 0x0000,
4543 0x0000, 0x0000, 0x0000, 0x0000, 0xa000, 0x4002,}, /* 10**1 */
4546 static unsigned EMUSHORT emtens[NTEN + 1][NE] =
4548 {0x2030, 0xcffc, 0xa1c3, 0x8123,
4549 0x2de3, 0x9fde, 0xd2ce, 0x04c8, 0xa6dd, 0x0ad8,}, /* 10**-4096 */
4550 {0x8264, 0xd2cb, 0xf2ea, 0x12d4,
4551 0x4925, 0x2de4, 0x3436, 0x534f, 0xceae, 0x256b,}, /* 10**-2048 */
4552 {0xf53f, 0xf698, 0x6bd3, 0x0158,
4553 0x87a6, 0xc0bd, 0xda57, 0x82a5, 0xa2a6, 0x32b5,},
4554 {0xe731, 0x04d4, 0xe3f2, 0xd332,
4555 0x7132, 0xd21c, 0xdb23, 0xee32, 0x9049, 0x395a,},
4556 {0xa23e, 0x5308, 0xfefb, 0x1155,
4557 0xfa91, 0x1939, 0x637a, 0x4325, 0xc031, 0x3cac,},
4558 {0xe26d, 0xdbde, 0xd05d, 0xb3f6,
4559 0xac7c, 0xe4a0, 0x64bc, 0x467c, 0xddd0, 0x3e55,},
4560 {0x2a20, 0x6224, 0x47b3, 0x98d7,
4561 0x3f23, 0xe9a5, 0xa539, 0xea27, 0xa87f, 0x3f2a,},
4562 {0x0b5b, 0x4af2, 0xa581, 0x18ed,
4563 0x67de, 0x94ba, 0x4539, 0x1ead, 0xcfb1, 0x3f94,},
4564 {0xbf71, 0xa9b3, 0x7989, 0xbe68,
4565 0x4c2e, 0xe15b, 0xc44d, 0x94be, 0xe695, 0x3fc9,},
4566 {0x3d4d, 0x7c3d, 0x36ba, 0x0d2b,
4567 0xfdc2, 0xcefc, 0x8461, 0x7711, 0xabcc, 0x3fe4,},
4568 {0xc155, 0xa4a8, 0x404e, 0x6113,
4569 0xd3c3, 0x652b, 0xe219, 0x1758, 0xd1b7, 0x3ff1,},
4570 {0xd70a, 0x70a3, 0x0a3d, 0xa3d7,
4571 0x3d70, 0xd70a, 0x70a3, 0x0a3d, 0xa3d7, 0x3ff8,},
4572 {0xcccd, 0xcccc, 0xcccc, 0xcccc,
4573 0xcccc, 0xcccc, 0xcccc, 0xcccc, 0xcccc, 0x3ffb,}, /* 10**-1 */
4576 /* LONG_DOUBLE_TYPE_SIZE is other than 128 */
4577 static unsigned EMUSHORT etens[NTEN + 1][NE] =
4579 {0xc94c, 0x979a, 0x8a20, 0x5202, 0xc460, 0x7525,}, /* 10**4096 */
4580 {0xa74d, 0x5de4, 0xc53d, 0x3b5d, 0x9e8b, 0x5a92,}, /* 10**2048 */
4581 {0x650d, 0x0c17, 0x8175, 0x7586, 0xc976, 0x4d48,},
4582 {0xcc65, 0x91c6, 0xa60e, 0xa0ae, 0xe319, 0x46a3,},
4583 {0xddbc, 0xde8d, 0x9df9, 0xebfb, 0xaa7e, 0x4351,},
4584 {0xc66f, 0x8cdf, 0x80e9, 0x47c9, 0x93ba, 0x41a8,},
4585 {0x3cbf, 0xa6d5, 0xffcf, 0x1f49, 0xc278, 0x40d3,},
4586 {0xf020, 0xb59d, 0x2b70, 0xada8, 0x9dc5, 0x4069,},
4587 {0x0000, 0x0000, 0x0400, 0xc9bf, 0x8e1b, 0x4034,},
4588 {0x0000, 0x0000, 0x0000, 0x2000, 0xbebc, 0x4019,},
4589 {0x0000, 0x0000, 0x0000, 0x0000, 0x9c40, 0x400c,},
4590 {0x0000, 0x0000, 0x0000, 0x0000, 0xc800, 0x4005,},
4591 {0x0000, 0x0000, 0x0000, 0x0000, 0xa000, 0x4002,}, /* 10**1 */
4594 static unsigned EMUSHORT emtens[NTEN + 1][NE] =
4596 {0x2de4, 0x9fde, 0xd2ce, 0x04c8, 0xa6dd, 0x0ad8,}, /* 10**-4096 */
4597 {0x4925, 0x2de4, 0x3436, 0x534f, 0xceae, 0x256b,}, /* 10**-2048 */
4598 {0x87a6, 0xc0bd, 0xda57, 0x82a5, 0xa2a6, 0x32b5,},
4599 {0x7133, 0xd21c, 0xdb23, 0xee32, 0x9049, 0x395a,},
4600 {0xfa91, 0x1939, 0x637a, 0x4325, 0xc031, 0x3cac,},
4601 {0xac7d, 0xe4a0, 0x64bc, 0x467c, 0xddd0, 0x3e55,},
4602 {0x3f24, 0xe9a5, 0xa539, 0xea27, 0xa87f, 0x3f2a,},
4603 {0x67de, 0x94ba, 0x4539, 0x1ead, 0xcfb1, 0x3f94,},
4604 {0x4c2f, 0xe15b, 0xc44d, 0x94be, 0xe695, 0x3fc9,},
4605 {0xfdc2, 0xcefc, 0x8461, 0x7711, 0xabcc, 0x3fe4,},
4606 {0xd3c3, 0x652b, 0xe219, 0x1758, 0xd1b7, 0x3ff1,},
4607 {0x3d71, 0xd70a, 0x70a3, 0x0a3d, 0xa3d7, 0x3ff8,},
4608 {0xcccd, 0xcccc, 0xcccc, 0xcccc, 0xcccc, 0x3ffb,}, /* 10**-1 */
4613 /* Convert float value X to ASCII string STRING with NDIG digits after
4614 the decimal point. */
4617 e24toasc (x, string, ndigs)
4618 unsigned EMUSHORT x[];
4622 unsigned EMUSHORT w[NI];
4625 etoasc (w, string, ndigs);
4628 /* Convert double value X to ASCII string STRING with NDIG digits after
4629 the decimal point. */
4632 e53toasc (x, string, ndigs)
4633 unsigned EMUSHORT x[];
4637 unsigned EMUSHORT w[NI];
4640 etoasc (w, string, ndigs);
4643 /* Convert double extended value X to ASCII string STRING with NDIG digits
4644 after the decimal point. */
4647 e64toasc (x, string, ndigs)
4648 unsigned EMUSHORT x[];
4652 unsigned EMUSHORT w[NI];
4655 etoasc (w, string, ndigs);
4658 /* Convert 128-bit long double value X to ASCII string STRING with NDIG digits
4659 after the decimal point. */
4662 e113toasc (x, string, ndigs)
4663 unsigned EMUSHORT x[];
4667 unsigned EMUSHORT w[NI];
4670 etoasc (w, string, ndigs);
4674 /* Convert e-type X to ASCII string STRING with NDIGS digits after
4675 the decimal point. */
4677 static char wstring[80]; /* working storage for ASCII output */
4680 etoasc (x, string, ndigs)
4681 unsigned EMUSHORT x[];
4686 unsigned EMUSHORT y[NI], t[NI], u[NI], w[NI];
4687 unsigned EMUSHORT *p, *r, *ten;
4688 unsigned EMUSHORT sign;
4689 int i, j, k, expon, rndsav;
4691 unsigned EMUSHORT m;
4702 sprintf (wstring, " NaN ");
4706 rndprc = NBITS; /* set to full precision */
4707 emov (x, y); /* retain external format */
4708 if (y[NE - 1] & 0x8000)
4711 y[NE - 1] &= 0x7fff;
4718 ten = &etens[NTEN][0];
4720 /* Test for zero exponent */
4723 for (k = 0; k < NE - 1; k++)
4726 goto tnzro; /* denormalized number */
4728 goto isone; /* valid all zeros */
4732 /* Test for infinity. */
4733 if (y[NE - 1] == 0x7fff)
4736 sprintf (wstring, " -Infinity ");
4738 sprintf (wstring, " Infinity ");
4742 /* Test for exponent nonzero but significand denormalized.
4743 * This is an error condition.
4745 if ((y[NE - 1] != 0) && ((y[NE - 2] & 0x8000) == 0))
4747 mtherr ("etoasc", DOMAIN);
4748 sprintf (wstring, "NaN");
4752 /* Compare to 1.0 */
4761 { /* Number is greater than 1 */
4762 /* Convert significand to an integer and strip trailing decimal zeros. */
4764 u[NE - 1] = EXONE + NBITS - 1;
4766 p = &etens[NTEN - 4][0];
4772 for (j = 0; j < NE - 1; j++)
4785 /* Rescale from integer significand */
4786 u[NE - 1] += y[NE - 1] - (unsigned int) (EXONE + NBITS - 1);
4788 /* Find power of 10 */
4792 /* An unordered compare result shouldn't happen here. */
4793 while (ecmp (ten, u) <= 0)
4795 if (ecmp (p, u) <= 0)
4808 { /* Number is less than 1.0 */
4809 /* Pad significand with trailing decimal zeros. */
4812 while ((y[NE - 2] & 0x8000) == 0)
4821 for (i = 0; i < NDEC + 1; i++)
4823 if ((w[NI - 1] & 0x7) != 0)
4825 /* multiply by 10 */
4838 if (eone[NE - 1] <= u[1])
4850 while (ecmp (eone, w) > 0)
4852 if (ecmp (p, w) >= 0)
4867 /* Find the first (leading) digit. */
4873 digit = equot[NI - 1];
4874 while ((digit == 0) && (ecmp (y, ezero) != 0))
4882 digit = equot[NI - 1];
4890 /* Examine number of digits requested by caller. */
4908 *s++ = (char)digit + '0';
4911 /* Generate digits after the decimal point. */
4912 for (k = 0; k <= ndigs; k++)
4914 /* multiply current number by 10, without normalizing */
4921 *s++ = (char) equot[NI - 1] + '0';
4923 digit = equot[NI - 1];
4926 /* round off the ASCII string */
4929 /* Test for critical rounding case in ASCII output. */
4933 if (ecmp (t, ezero) != 0)
4934 goto roun; /* round to nearest */
4935 if ((*(s - 1) & 1) == 0)
4936 goto doexp; /* round to even */
4938 /* Round up and propagate carry-outs */
4942 /* Carry out to most significant digit? */
4949 /* Most significant digit carries to 10? */
4957 /* Round up and carry out from less significant digits */
4969 sprintf (ss, "e+%d", expon);
4971 sprintf (ss, "e%d", expon);
4973 sprintf (ss, "e%d", expon);
4976 /* copy out the working string */
4979 while (*ss == ' ') /* strip possible leading space */
4981 while ((*s++ = *ss++) != '\0')
4986 /* Convert ASCII string to floating point.
4988 Numeric input is a free format decimal number of any length, with
4989 or without decimal point. Entering E after the number followed by an
4990 integer number causes the second number to be interpreted as a power of
4991 10 to be multiplied by the first number (i.e., "scientific" notation). */
4993 /* Convert ASCII string S to single precision float value Y. */
4998 unsigned EMUSHORT *y;
5004 /* Convert ASCII string S to double precision value Y. */
5009 unsigned EMUSHORT *y;
5011 #if defined(DEC) || defined(IBM)
5023 /* Convert ASCII string S to double extended value Y. */
5028 unsigned EMUSHORT *y;
5033 /* Convert ASCII string S to 128-bit long double Y. */
5038 unsigned EMUSHORT *y;
5040 asctoeg (s, y, 113);
5043 /* Convert ASCII string S to e type Y. */
5048 unsigned EMUSHORT *y;
5050 asctoeg (s, y, NBITS);
5053 /* Convert ASCII string SS to e type Y, with a specified rounding precision
5057 asctoeg (ss, y, oprec)
5059 unsigned EMUSHORT *y;
5062 unsigned EMUSHORT yy[NI], xt[NI], tt[NI];
5063 int esign, decflg, sgnflg, nexp, exp, prec, lost;
5064 int k, trail, c, rndsav;
5066 unsigned EMUSHORT nsign, *p;
5067 char *sp, *s, *lstr;
5069 /* Copy the input string. */
5070 lstr = (char *) alloca (strlen (ss) + 1);
5072 while (*s == ' ') /* skip leading spaces */
5075 while ((*sp++ = *s++) != '\0')
5080 rndprc = NBITS; /* Set to full precision */
5093 if ((k >= 0) && (k <= 9))
5095 /* Ignore leading zeros */
5096 if ((prec == 0) && (decflg == 0) && (k == 0))
5098 /* Identify and strip trailing zeros after the decimal point. */
5099 if ((trail == 0) && (decflg != 0))
5102 while ((*sp >= '0') && (*sp <= '9'))
5104 /* Check for syntax error */
5106 if ((c != 'e') && (c != 'E') && (c != '\0')
5107 && (c != '\n') && (c != '\r') && (c != ' ')
5118 /* If enough digits were given to more than fill up the yy register,
5119 continuing until overflow into the high guard word yy[2]
5120 guarantees that there will be a roundoff bit at the top
5121 of the low guard word after normalization. */
5126 nexp += 1; /* count digits after decimal point */
5127 eshup1 (yy); /* multiply current number by 10 */
5133 xt[NI - 2] = (unsigned EMUSHORT) k;
5138 /* Mark any lost non-zero digit. */
5140 /* Count lost digits before the decimal point. */
5155 case '.': /* decimal point */
5185 mtherr ("asctoe", DOMAIN);
5194 /* Exponent interpretation */
5196 /* 0.0eXXX is zero, regardless of XXX. Check for the 0.0. */
5197 for (k = 0; k < NI; k++)
5208 /* check for + or - */
5216 while ((*s >= '0') && (*s <= '9'))
5220 if (exp > -(MINDECEXP))
5230 if (exp > MAXDECEXP)
5234 yy[E] = 0x7fff; /* infinity */
5237 if (exp < MINDECEXP)
5246 /* Pad trailing zeros to minimize power of 10, per IEEE spec. */
5247 while ((nexp > 0) && (yy[2] == 0))
5259 if ((k = enormlz (yy)) > NBITS)
5264 lexp = (EXONE - 1 + NBITS) - k;
5265 emdnorm (yy, lost, 0, lexp, 64);
5267 /* Convert to external format:
5269 Multiply by 10**nexp. If precision is 64 bits,
5270 the maximum relative error incurred in forming 10**n
5271 for 0 <= n <= 324 is 8.2e-20, at 10**180.
5272 For 0 <= n <= 999, the peak relative error is 1.4e-19 at 10**947.
5273 For 0 >= n >= -999, it is -1.55e-19 at 10**-435. */
5288 /* Punt. Can't handle this without 2 divides. */
5289 emovi (etens[0], tt);
5296 p = &etens[NTEN][0];
5306 while (exp <= MAXP);
5324 /* Round and convert directly to the destination type */
5326 lexp -= EXONE - 0x3ff;
5328 else if (oprec == 24 || oprec == 32)
5329 lexp -= (EXONE - 0x7f);
5332 else if (oprec == 24 || oprec == 56)
5333 lexp -= EXONE - (0x41 << 2);
5335 else if (oprec == 24)
5336 lexp -= EXONE - 0177;
5340 else if (oprec == 56)
5341 lexp -= EXONE - 0201;
5344 emdnorm (yy, k, 0, lexp, 64);
5354 todec (yy, y); /* see etodec.c */
5359 toibm (yy, y, DFmode);
5364 toc4x (yy, y, HFmode);
5388 /* Return Y = largest integer not greater than X (truncated toward minus
5391 static unsigned EMUSHORT bmask[] =
5414 unsigned EMUSHORT x[], y[];
5416 register unsigned EMUSHORT *p;
5418 unsigned EMUSHORT f[NE];
5420 emov (x, f); /* leave in external format */
5421 expon = (int) f[NE - 1];
5422 e = (expon & 0x7fff) - (EXONE - 1);
5428 /* number of bits to clear out */
5440 /* clear the remaining bits */
5442 /* truncate negatives toward minus infinity */
5445 if ((unsigned EMUSHORT) expon & (unsigned EMUSHORT) 0x8000)
5447 for (i = 0; i < NE - 1; i++)
5460 /* Return S and EXP such that S * 2^EXP = X and .5 <= S < 1.
5461 For example, 1.1 = 0.55 * 2^1. */
5465 unsigned EMUSHORT x[];
5467 unsigned EMUSHORT s[];
5469 unsigned EMUSHORT xi[NI];
5473 /* Handle denormalized numbers properly using long integer exponent. */
5474 li = (EMULONG) ((EMUSHORT) xi[1]);
5482 *exp = (int) (li - 0x3ffe);
5486 /* Return e type Y = X * 2^PWR2. */
5490 unsigned EMUSHORT x[];
5492 unsigned EMUSHORT y[];
5494 unsigned EMUSHORT xi[NI];
5502 emdnorm (xi, i, i, li, 64);
5508 /* C = remainder after dividing B by A, all e type values.
5509 Least significant integer quotient bits left in EQUOT. */
5513 unsigned EMUSHORT a[], b[], c[];
5515 unsigned EMUSHORT den[NI], num[NI];
5519 || (ecmp (a, ezero) == 0)
5527 if (ecmp (a, ezero) == 0)
5529 mtherr ("eremain", SING);
5535 eiremain (den, num);
5536 /* Sign of remainder = sign of quotient */
5545 /* Return quotient of exploded e-types NUM / DEN in EQUOT,
5546 remainder in NUM. */
5550 unsigned EMUSHORT den[], num[];
5553 unsigned EMUSHORT j;
5556 ld -= enormlz (den);
5558 ln -= enormlz (num);
5562 if (ecmpm (den, num) <= 0)
5574 emdnorm (num, 0, 0, ln, 0);
5577 /* Report an error condition CODE encountered in function NAME.
5578 CODE is one of the following:
5580 Mnemonic Value Significance
5582 DOMAIN 1 argument domain error
5583 SING 2 function singularity
5584 OVERFLOW 3 overflow range error
5585 UNDERFLOW 4 underflow range error
5586 TLOSS 5 total loss of precision
5587 PLOSS 6 partial loss of precision
5588 INVALID 7 NaN - producing operation
5589 EDOM 33 Unix domain error code
5590 ERANGE 34 Unix range error code
5592 The order of appearance of the following messages is bound to the
5593 error codes defined above. */
5596 static char *ermsg[NMSGS] =
5598 "unknown", /* error code 0 */
5599 "domain", /* error code 1 */
5600 "singularity", /* et seq. */
5603 "total loss of precision",
5604 "partial loss of precision",
5618 /* The string passed by the calling program is supposed to be the
5619 name of the function in which the error occurred.
5620 The code argument selects which error message string will be printed. */
5622 if ((code <= 0) || (code >= NMSGS))
5624 sprintf (errstr, " %s %s error", name, ermsg[code]);
5627 /* Set global error message word */
5632 /* Convert DEC double precision D to e type E. */
5636 unsigned EMUSHORT *d;
5637 unsigned EMUSHORT *e;
5639 unsigned EMUSHORT y[NI];
5640 register unsigned EMUSHORT r, *p;
5642 ecleaz (y); /* start with a zero */
5643 p = y; /* point to our number */
5644 r = *d; /* get DEC exponent word */
5645 if (*d & (unsigned int) 0x8000)
5646 *p = 0xffff; /* fill in our sign */
5647 ++p; /* bump pointer to our exponent word */
5648 r &= 0x7fff; /* strip the sign bit */
5649 if (r == 0) /* answer = 0 if high order DEC word = 0 */
5653 r >>= 7; /* shift exponent word down 7 bits */
5654 r += EXONE - 0201; /* subtract DEC exponent offset */
5655 /* add our e type exponent offset */
5656 *p++ = r; /* to form our exponent */
5658 r = *d++; /* now do the high order mantissa */
5659 r &= 0177; /* strip off the DEC exponent and sign bits */
5660 r |= 0200; /* the DEC understood high order mantissa bit */
5661 *p++ = r; /* put result in our high guard word */
5663 *p++ = *d++; /* fill in the rest of our mantissa */
5667 eshdn8 (y); /* shift our mantissa down 8 bits */
5672 /* Convert e type X to DEC double precision D. */
5676 unsigned EMUSHORT *x, *d;
5678 unsigned EMUSHORT xi[NI];
5683 /* Adjust exponent for offsets. */
5684 exp = (EMULONG) xi[E] - (EXONE - 0201);
5685 /* Round off to nearest or even. */
5688 emdnorm (xi, 0, 0, exp, 64);
5693 /* Convert exploded e-type X, that has already been rounded to
5694 56-bit precision, to DEC format double Y. */
5698 unsigned EMUSHORT *x, *y;
5700 unsigned EMUSHORT i;
5701 unsigned EMUSHORT *p;
5740 /* Convert IBM single/double precision to e type. */
5744 unsigned EMUSHORT *d;
5745 unsigned EMUSHORT *e;
5746 enum machine_mode mode;
5748 unsigned EMUSHORT y[NI];
5749 register unsigned EMUSHORT r, *p;
5752 ecleaz (y); /* start with a zero */
5753 p = y; /* point to our number */
5754 r = *d; /* get IBM exponent word */
5755 if (*d & (unsigned int) 0x8000)
5756 *p = 0xffff; /* fill in our sign */
5757 ++p; /* bump pointer to our exponent word */
5758 r &= 0x7f00; /* strip the sign bit */
5759 r >>= 6; /* shift exponent word down 6 bits */
5760 /* in fact shift by 8 right and 2 left */
5761 r += EXONE - (0x41 << 2); /* subtract IBM exponent offset */
5762 /* add our e type exponent offset */
5763 *p++ = r; /* to form our exponent */
5765 *p++ = *d++ & 0xff; /* now do the high order mantissa */
5766 /* strip off the IBM exponent and sign bits */
5767 if (mode != SFmode) /* there are only 2 words in SFmode */
5769 *p++ = *d++; /* fill in the rest of our mantissa */
5774 if (y[M] == 0 && y[M+1] == 0 && y[M+2] == 0 && y[M+3] == 0)
5777 y[E] -= 5 + enormlz (y); /* now normalise the mantissa */
5778 /* handle change in RADIX */
5784 /* Convert e type to IBM single/double precision. */
5788 unsigned EMUSHORT *x, *d;
5789 enum machine_mode mode;
5791 unsigned EMUSHORT xi[NI];
5796 exp = (EMULONG) xi[E] - (EXONE - (0x41 << 2)); /* adjust exponent for offsets */
5797 /* round off to nearest or even */
5800 emdnorm (xi, 0, 0, exp, 64);
5802 toibm (xi, d, mode);
5807 unsigned EMUSHORT *x, *y;
5808 enum machine_mode mode;
5810 unsigned EMUSHORT i;
5811 unsigned EMUSHORT *p;
5861 /* Convert C4X single/double precision to e type. */
5865 unsigned EMUSHORT *d;
5866 unsigned EMUSHORT *e;
5867 enum machine_mode mode;
5869 unsigned EMUSHORT y[NI];
5877 /* Short-circuit the zero case. */
5878 if ((d[0] == 0x8000)
5880 && ((mode == QFmode) || ((d[2] == 0x0000) && (d[3] == 0x0000))))
5891 ecleaz (y); /* start with a zero */
5892 r = d[0]; /* get sign/exponent part */
5893 if (r & (unsigned int) 0x0080)
5895 y[0] = 0xffff; /* fill in our sign */
5903 r >>= 8; /* Shift exponent word down 8 bits. */
5904 if (r & 0x80) /* Make the exponent negative if it is. */
5906 r = r | (~0 & ~0xff);
5911 /* Now do the high order mantissa. We don't "or" on the high bit
5912 because it is 2 (not 1) and is handled a little differently
5917 if (mode != QFmode) /* There are only 2 words in QFmode. */
5919 y[M+2] = d[2]; /* Fill in the rest of our mantissa. */
5929 /* Now do the two's complement on the data. */
5931 carry = 1; /* Initially add 1 for the two's complement. */
5932 for (i=size + M; i > M; i--)
5934 if (carry && (y[i] == 0x0000))
5936 /* We overflowed into the next word, carry is the same. */
5937 y[i] = carry ? 0x0000 : 0xffff;
5941 /* No overflow, just invert and add carry. */
5942 y[i] = ((~y[i]) + carry) & 0xffff;
5957 /* Add our e type exponent offset to form our exponent. */
5961 /* Now do the high order mantissa strip off the exponent and sign
5962 bits and add the high 1 bit. */
5963 y[M] = d[0] & 0x7f | 0x80;
5966 if (mode != QFmode) /* There are only 2 words in QFmode. */
5968 y[M+2] = d[2]; /* Fill in the rest of our mantissa. */
5978 /* Convert e type to C4X single/double precision. */
5982 unsigned EMUSHORT *x, *d;
5983 enum machine_mode mode;
5985 unsigned EMUSHORT xi[NI];
5991 /* Adjust exponent for offsets. */
5992 exp = (EMULONG) xi[E] - (EXONE - 0x7f);
5994 /* Round off to nearest or even. */
5996 rndprc = mode == QFmode ? 24 : 32;
5997 emdnorm (xi, 0, 0, exp, 64);
5999 toc4x (xi, d, mode);
6004 unsigned EMUSHORT *x, *y;
6005 enum machine_mode mode;
6012 /* Short-circuit the zero case */
6013 if ((x[0] == 0) /* Zero exponent and sign */
6015 && (x[M] == 0) /* The rest is for zero mantissa */
6017 /* Only check for double if necessary */
6018 && ((mode == QFmode) || ((x[M+2] == 0) && (x[M+3] == 0))))
6020 /* We have a zero. Put it into the output and return. */
6033 /* Negative number require a two's complement conversion of the
6039 i = ((int) x[1]) - 0x7f;
6041 /* Now add 1 to the inverted data to do the two's complement. */
6051 x[v] = carry ? 0x0000 : 0xffff;
6055 x[v] = ((~x[v]) + carry) & 0xffff;
6061 /* The following is a special case. The C4X negative float requires
6062 a zero in the high bit (because the format is (2 - x) x 2^m), so
6063 if a one is in that bit, we have to shift left one to get rid
6064 of it. This only occurs if the number is -1 x 2^m. */
6065 if (x[M+1] & 0x8000)
6067 /* This is the case of -1 x 2^m, we have to rid ourselves of the
6068 high sign bit and shift the exponent. */
6075 i = ((int) x[1]) - 0x7f;
6078 if ((i < -128) || (i > 127))
6093 y[0] |= ((i & 0xff) << 8);
6097 y[0] |= x[M] & 0x7f;
6107 /* Output a binary NaN bit pattern in the target machine's format. */
6109 /* If special NaN bit patterns are required, define them in tm.h
6110 as arrays of unsigned 16-bit shorts. Otherwise, use the default
6116 unsigned EMUSHORT TFbignan[8] =
6117 {0x7fff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff};
6118 unsigned EMUSHORT TFlittlenan[8] = {0, 0, 0, 0, 0, 0, 0x8000, 0xffff};
6126 unsigned EMUSHORT XFbignan[6] =
6127 {0x7fff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff};
6128 unsigned EMUSHORT XFlittlenan[6] = {0, 0, 0, 0xc000, 0xffff, 0};
6136 unsigned EMUSHORT DFbignan[4] = {0x7fff, 0xffff, 0xffff, 0xffff};
6137 unsigned EMUSHORT DFlittlenan[4] = {0, 0, 0, 0xfff8};
6145 unsigned EMUSHORT SFbignan[2] = {0x7fff, 0xffff};
6146 unsigned EMUSHORT SFlittlenan[2] = {0, 0xffc0};
6152 make_nan (nan, sign, mode)
6153 unsigned EMUSHORT *nan;
6155 enum machine_mode mode;
6158 unsigned EMUSHORT *p;
6162 /* Possibly the `reserved operand' patterns on a VAX can be
6163 used like NaN's, but probably not in the same way as IEEE. */
6164 #if !defined(DEC) && !defined(IBM) && !defined(C4X)
6167 if (REAL_WORDS_BIG_ENDIAN)
6175 if (REAL_WORDS_BIG_ENDIAN)
6183 if (REAL_WORDS_BIG_ENDIAN)
6192 if (REAL_WORDS_BIG_ENDIAN)
6202 if (REAL_WORDS_BIG_ENDIAN)
6203 *nan++ = (sign << 15) | *p++;
6206 if (! REAL_WORDS_BIG_ENDIAN)
6207 *nan = (sign << 15) | *p;
6210 /* This is the inverse of the function `etarsingle' invoked by
6211 REAL_VALUE_TO_TARGET_SINGLE. */
6214 ereal_unto_float (f)
6218 unsigned EMUSHORT s[2];
6219 unsigned EMUSHORT e[NE];
6221 /* Convert 32 bit integer to array of 16 bit pieces in target machine order.
6222 This is the inverse operation to what the function `endian' does. */
6223 if (REAL_WORDS_BIG_ENDIAN)
6225 s[0] = (unsigned EMUSHORT) (f >> 16);
6226 s[1] = (unsigned EMUSHORT) f;
6230 s[0] = (unsigned EMUSHORT) f;
6231 s[1] = (unsigned EMUSHORT) (f >> 16);
6233 /* Convert and promote the target float to E-type. */
6235 /* Output E-type to REAL_VALUE_TYPE. */
6241 /* This is the inverse of the function `etardouble' invoked by
6242 REAL_VALUE_TO_TARGET_DOUBLE. */
6245 ereal_unto_double (d)
6249 unsigned EMUSHORT s[4];
6250 unsigned EMUSHORT e[NE];
6252 /* Convert array of HOST_WIDE_INT to equivalent array of 16-bit pieces. */
6253 if (REAL_WORDS_BIG_ENDIAN)
6255 s[0] = (unsigned EMUSHORT) (d[0] >> 16);
6256 s[1] = (unsigned EMUSHORT) d[0];
6257 s[2] = (unsigned EMUSHORT) (d[1] >> 16);
6258 s[3] = (unsigned EMUSHORT) d[1];
6262 /* Target float words are little-endian. */
6263 s[0] = (unsigned EMUSHORT) d[0];
6264 s[1] = (unsigned EMUSHORT) (d[0] >> 16);
6265 s[2] = (unsigned EMUSHORT) d[1];
6266 s[3] = (unsigned EMUSHORT) (d[1] >> 16);
6268 /* Convert target double to E-type. */
6270 /* Output E-type to REAL_VALUE_TYPE. */
6276 /* Convert an SFmode target `float' value to a REAL_VALUE_TYPE.
6277 This is somewhat like ereal_unto_float, but the input types
6278 for these are different. */
6281 ereal_from_float (f)
6285 unsigned EMUSHORT s[2];
6286 unsigned EMUSHORT e[NE];
6288 /* Convert 32 bit integer to array of 16 bit pieces in target machine order.
6289 This is the inverse operation to what the function `endian' does. */
6290 if (REAL_WORDS_BIG_ENDIAN)
6292 s[0] = (unsigned EMUSHORT) (f >> 16);
6293 s[1] = (unsigned EMUSHORT) f;
6297 s[0] = (unsigned EMUSHORT) f;
6298 s[1] = (unsigned EMUSHORT) (f >> 16);
6300 /* Convert and promote the target float to E-type. */
6302 /* Output E-type to REAL_VALUE_TYPE. */
6308 /* Convert a DFmode target `double' value to a REAL_VALUE_TYPE.
6309 This is somewhat like ereal_unto_double, but the input types
6310 for these are different.
6312 The DFmode is stored as an array of HOST_WIDE_INT in the target's
6313 data format, with no holes in the bit packing. The first element
6314 of the input array holds the bits that would come first in the
6315 target computer's memory. */
6318 ereal_from_double (d)
6322 unsigned EMUSHORT s[4];
6323 unsigned EMUSHORT e[NE];
6325 /* Convert array of HOST_WIDE_INT to equivalent array of 16-bit pieces. */
6326 if (REAL_WORDS_BIG_ENDIAN)
6328 s[0] = (unsigned EMUSHORT) (d[0] >> 16);
6329 s[1] = (unsigned EMUSHORT) d[0];
6330 #if HOST_BITS_PER_WIDE_INT == 32
6331 s[2] = (unsigned EMUSHORT) (d[1] >> 16);
6332 s[3] = (unsigned EMUSHORT) d[1];
6334 /* In this case the entire target double is contained in the
6335 first array element. The second element of the input is
6337 s[2] = (unsigned EMUSHORT) (d[0] >> 48);
6338 s[3] = (unsigned EMUSHORT) (d[0] >> 32);
6343 /* Target float words are little-endian. */
6344 s[0] = (unsigned EMUSHORT) d[0];
6345 s[1] = (unsigned EMUSHORT) (d[0] >> 16);
6346 #if HOST_BITS_PER_WIDE_INT == 32
6347 s[2] = (unsigned EMUSHORT) d[1];
6348 s[3] = (unsigned EMUSHORT) (d[1] >> 16);
6350 s[2] = (unsigned EMUSHORT) (d[0] >> 32);
6351 s[3] = (unsigned EMUSHORT) (d[0] >> 48);
6354 /* Convert target double to E-type. */
6356 /* Output E-type to REAL_VALUE_TYPE. */
6363 /* Convert target computer unsigned 64-bit integer to e-type.
6364 The endian-ness of DImode follows the convention for integers,
6365 so we use WORDS_BIG_ENDIAN here, not REAL_WORDS_BIG_ENDIAN. */
6369 unsigned EMUSHORT *di; /* Address of the 64-bit int. */
6370 unsigned EMUSHORT *e;
6372 unsigned EMUSHORT yi[NI];
6376 if (WORDS_BIG_ENDIAN)
6378 for (k = M; k < M + 4; k++)
6383 for (k = M + 3; k >= M; k--)
6386 yi[E] = EXONE + 47; /* exponent if normalize shift count were 0 */
6387 if ((k = enormlz (yi)) > NBITS)/* normalize the significand */
6388 ecleaz (yi); /* it was zero */
6390 yi[E] -= (unsigned EMUSHORT) k;/* subtract shift count from exponent */
6394 /* Convert target computer signed 64-bit integer to e-type. */
6398 unsigned EMUSHORT *di; /* Address of the 64-bit int. */
6399 unsigned EMUSHORT *e;
6401 unsigned EMULONG acc;
6402 unsigned EMUSHORT yi[NI];
6403 unsigned EMUSHORT carry;
6407 if (WORDS_BIG_ENDIAN)
6409 for (k = M; k < M + 4; k++)
6414 for (k = M + 3; k >= M; k--)
6417 /* Take absolute value */
6423 for (k = M + 3; k >= M; k--)
6425 acc = (unsigned EMULONG) (~yi[k] & 0xffff) + carry;
6432 yi[E] = EXONE + 47; /* exponent if normalize shift count were 0 */
6433 if ((k = enormlz (yi)) > NBITS)/* normalize the significand */
6434 ecleaz (yi); /* it was zero */
6436 yi[E] -= (unsigned EMUSHORT) k;/* subtract shift count from exponent */
6443 /* Convert e-type to unsigned 64-bit int. */
6447 unsigned EMUSHORT *x;
6448 unsigned EMUSHORT *i;
6450 unsigned EMUSHORT xi[NI];
6459 k = (int) xi[E] - (EXONE - 1);
6462 for (j = 0; j < 4; j++)
6468 for (j = 0; j < 4; j++)
6471 warning ("overflow on truncation to integer");
6476 /* Shift more than 16 bits: first shift up k-16 mod 16,
6477 then shift up by 16's. */
6478 j = k - ((k >> 4) << 4);
6482 if (WORDS_BIG_ENDIAN)
6493 if (WORDS_BIG_ENDIAN)
6498 while ((k -= 16) > 0);
6502 /* shift not more than 16 bits */
6507 if (WORDS_BIG_ENDIAN)
6526 /* Convert e-type to signed 64-bit int. */
6530 unsigned EMUSHORT *x;
6531 unsigned EMUSHORT *i;
6533 unsigned EMULONG acc;
6534 unsigned EMUSHORT xi[NI];
6535 unsigned EMUSHORT carry;
6536 unsigned EMUSHORT *isave;
6540 k = (int) xi[E] - (EXONE - 1);
6543 for (j = 0; j < 4; j++)
6549 for (j = 0; j < 4; j++)
6552 warning ("overflow on truncation to integer");
6558 /* Shift more than 16 bits: first shift up k-16 mod 16,
6559 then shift up by 16's. */
6560 j = k - ((k >> 4) << 4);
6564 if (WORDS_BIG_ENDIAN)
6575 if (WORDS_BIG_ENDIAN)
6580 while ((k -= 16) > 0);
6584 /* shift not more than 16 bits */
6587 if (WORDS_BIG_ENDIAN)
6603 /* Negate if negative */
6607 if (WORDS_BIG_ENDIAN)
6609 for (k = 0; k < 4; k++)
6611 acc = (unsigned EMULONG) (~(*isave) & 0xffff) + carry;
6612 if (WORDS_BIG_ENDIAN)
6624 /* Longhand square root routine. */
6627 static int esqinited = 0;
6628 static unsigned short sqrndbit[NI];
6632 unsigned EMUSHORT *x, *y;
6634 unsigned EMUSHORT temp[NI], num[NI], sq[NI], xx[NI];
6636 int i, j, k, n, nlups;
6641 sqrndbit[NI - 2] = 1;
6644 /* Check for arg <= 0 */
6645 i = ecmp (x, ezero);
6650 mtherr ("esqrt", DOMAIN);
6666 /* Bring in the arg and renormalize if it is denormal. */
6668 m = (EMULONG) xx[1]; /* local long word exponent */
6672 /* Divide exponent by 2 */
6674 exp = (unsigned short) ((m / 2) + 0x3ffe);
6676 /* Adjust if exponent odd */
6686 n = 8; /* get 8 bits of result per inner loop */
6692 /* bring in next word of arg */
6694 num[NI - 1] = xx[j + 3];
6695 /* Do additional bit on last outer loop, for roundoff. */
6698 for (i = 0; i < n; i++)
6700 /* Next 2 bits of arg */
6703 /* Shift up answer */
6705 /* Make trial divisor */
6706 for (k = 0; k < NI; k++)
6709 eaddm (sqrndbit, temp);
6710 /* Subtract and insert answer bit if it goes in */
6711 if (ecmpm (temp, num) <= 0)
6721 /* Adjust for extra, roundoff loop done. */
6722 exp += (NBITS - 1) - rndprc;
6724 /* Sticky bit = 1 if the remainder is nonzero. */
6726 for (i = 3; i < NI; i++)
6729 /* Renormalize and round off. */
6730 emdnorm (sq, k, 0, exp, 64);
6734 #endif /* EMU_NON_COMPILE not defined */
6736 /* Return the binary precision of the significand for a given
6737 floating point mode. The mode can hold an integer value
6738 that many bits wide, without losing any bits. */
6741 significand_size (mode)
6742 enum machine_mode mode;
6745 /* Don't test the modes, but their sizes, lest this
6746 code won't work for BITS_PER_UNIT != 8 . */
6748 switch (GET_MODE_BITSIZE (mode))
6752 #if TARGET_FLOAT_FORMAT == C4X_FLOAT_FORMAT
6759 #if TARGET_FLOAT_FORMAT == IEEE_FLOAT_FORMAT
6762 #if TARGET_FLOAT_FORMAT == IBM_FLOAT_FORMAT
6765 #if TARGET_FLOAT_FORMAT == VAX_FLOAT_FORMAT
6768 #if TARGET_FLOAT_FORMAT == C4X_FLOAT_FORMAT