/* real.c - implementation of REAL_ARITHMETIC, REAL_VALUE_ATOF,
and support for XFmode IEEE extended real floating point arithmetic.
- Copyright (C) 1993, 1994 Free Software Foundation, Inc.
+ Copyright (C) 1993, 1994, 1995, 1996, 1997 Free Software Foundation, Inc.
Contributed by Stephen L. Moshier (moshier@world.std.com).
This file is part of GNU CC.
You should have received a copy of the GNU General Public License
along with GNU CC; see the file COPYING. If not, write to
-the Free Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA. */
+the Free Software Foundation, 59 Temple Place - Suite 330,
+Boston, MA 02111-1307, USA. */
+#include "config.h"
#include <stdio.h>
#include <errno.h>
-#include "config.h"
#include "tree.h"
#ifndef errno
its decimal conversion functions can be used if desired (see
real.h).
-The first part of this file interfaces gcc to ieee.c, which is a
-floating point arithmetic suite that was not written with gcc in
-mind. The interface is followed by ieee.c itself and related
-items. Avoid changing ieee.c unless you have suitable test
-programs available. A special version of the PARANOIA floating
-point arithmetic tester, modified for this purpose, can be found
-on usc.edu : /pub/C-numanal/ieeetest.zoo. Some tutorial
-information on ieee.c is given in my book: S. L. Moshier,
-_Methods and Programs for Mathematical Functions_, Prentice-Hall
-or Simon & Schuster Int'l, 1989. A library of XFmode elementary
-transcendental functions can be obtained by ftp from
-research.att.com: netlib/cephes/ldouble.shar.Z */
+The first part of this file interfaces gcc to a floating point
+arithmetic suite that was not written with gcc in mind. Avoid
+changing the low-level arithmetic routines unless you have suitable
+test programs available. A special version of the PARANOIA floating
+point arithmetic tester, modified for this purpose, can be found on
+usc.edu: /pub/C-numanal/ieeetest.zoo. Other tests, and libraries of
+XFmode and TFmode transcendental functions, can be obtained by ftp from
+netlib.att.com: netlib/cephes. */
\f
/* Type of computer arithmetic.
- Only one of DEC, IBM, MIEEE, IBMPC, or UNK should get defined.
-
- `MIEEE' refers generically to big-endian IEEE floating-point data
- structure. This definition should work in SFmode `float' type and
- DFmode `double' type on virtually all big-endian IEEE machines.
- If LONG_DOUBLE_TYPE_SIZE has been defined to be 96, then MIEEE
- also invokes the particular XFmode (`long double' type) data
- structure used by the Motorola 680x0 series processors.
-
- `IBMPC' refers generally to little-endian IEEE machines. In this
- case, if LONG_DOUBLE_TYPE_SIZE has been defined to be 96, then
- IBMPC also invokes the particular XFmode `long double' data
- structure used by the Intel 80x86 series processors.
+ Only one of DEC, IBM, IEEE, or UNK should get defined.
+
+ `IEEE', when REAL_WORDS_BIG_ENDIAN is non-zero, refers generically
+ to big-endian IEEE floating-point data structure. This definition
+ should work in SFmode `float' type and DFmode `double' type on
+ virtually all big-endian IEEE machines. If LONG_DOUBLE_TYPE_SIZE
+ has been defined to be 96, then IEEE also invokes the particular
+ XFmode (`long double' type) data structure used by the Motorola
+ 680x0 series processors.
+
+ `IEEE', when REAL_WORDS_BIG_ENDIAN is zero, refers generally to
+ little-endian IEEE machines. In this case, if LONG_DOUBLE_TYPE_SIZE
+ has been defined to be 96, then IEEE also invokes the particular
+ XFmode `long double' data structure used by the Intel 80x86 series
+ processors.
`DEC' refers specifically to the Digital Equipment Corp PDP-11
and VAX floating point data structure. This model currently
#define IBM 1
#else /* it's also not an IBM */
#if TARGET_FLOAT_FORMAT == IEEE_FLOAT_FORMAT
-#if FLOAT_WORDS_BIG_ENDIAN
-/* Motorola IEEE, high order words come first (Sun workstation): */
-#define MIEEE 1
-#else /* not big-endian */
-/* Intel IEEE, low order words come first:
- */
-#define IBMPC 1
-#endif /* big-endian */
+#define IEEE
#else /* it's not IEEE either */
-/* UNKnown arithmetic. We don't support this and can't go on. */
+/* UNKnown arithmetic. We don't support this and can't go on. */
unknown arithmetic type
#define UNK 1
#endif /* not IEEE */
#endif /* not IBM */
#endif /* not VAX */
+#define REAL_WORDS_BIG_ENDIAN FLOAT_WORDS_BIG_ENDIAN
+
#else
/* REAL_ARITHMETIC not defined means that the *host's* data
structure will be used. It may differ by endian-ness from the
#define IBM 1
#else /* it's also not an IBM */
#if HOST_FLOAT_FORMAT == IEEE_FLOAT_FORMAT
-#if HOST_FLOAT_WORDS_BIG_ENDIAN
-#define MIEEE 1
-#else /* not big-endian */
-#define IBMPC 1
-#endif /* big-endian */
+#define IEEE
#else /* it's not IEEE either */
unknown arithmetic type
#define UNK 1
#endif /* not IBM */
#endif /* not VAX */
+#define REAL_WORDS_BIG_ENDIAN HOST_FLOAT_WORDS_BIG_ENDIAN
+
#endif /* REAL_ARITHMETIC not defined */
/* Define INFINITY for support of infinity.
#define NANS
#endif
-/* Support of NaNs requires support of infinity. */
+/* Support of NaNs requires support of infinity. */
#ifdef NANS
#ifndef INFINITY
#define INFINITY
#endif
\f
/* Find a host integer type that is at least 16 bits wide,
- and another type at least twice whatever that size is. */
+ and another type at least twice whatever that size is. */
#if HOST_BITS_PER_CHAR >= 16
#define EMUSHORT char
#define EMUSHORT_SIZE HOST_BITS_PER_LONG
#define EMULONG_SIZE (2 * HOST_BITS_PER_LONG)
#else
-/* You will have to modify this program to have a smaller unit size. */
+/* You will have to modify this program to have a smaller unit size. */
#define EMU_NON_COMPILE
#endif
#endif
#if HOST_BITS_PER_LONG >= EMULONG_SIZE
#define EMULONG long
#else
-#if HOST_BITS_PER_LONG_LONG >= EMULONG_SIZE
+#if HOST_BITS_PER_LONGLONG >= EMULONG_SIZE
#define EMULONG long long int
#else
-/* You will have to modify this program to have a smaller unit size. */
+/* You will have to modify this program to have a smaller unit size. */
#define EMU_NON_COMPILE
#endif
#endif
#endif
-/* The host interface doesn't work if no 16-bit size exists. */
+/* The host interface doesn't work if no 16-bit size exists. */
#if EMUSHORT_SIZE != 16
#define EMU_NON_COMPILE
#endif
-/* OK to continue compilation. */
+/* OK to continue compilation. */
#ifndef EMU_NON_COMPILE
/* Construct macros to translate between REAL_VALUE_TYPE and e type.
#define NE 6
#define MAXDECEXP 4932
#define MINDECEXP -4956
-#define GET_REAL(r,e) bcopy (r, e, 2*NE)
-#define PUT_REAL(e,r) bcopy (e, r, 2*NE)
+#define GET_REAL(r,e) bcopy ((char *) r, (char *) e, 2*NE)
+#define PUT_REAL(e,r) bcopy ((char *) e, (char *) r, 2*NE)
#else /* no XFmode */
#if LONG_DOUBLE_TYPE_SIZE == 128
#define NE 10
#define MAXDECEXP 4932
#define MINDECEXP -4977
-#define GET_REAL(r,e) bcopy (r, e, 2*NE)
-#define PUT_REAL(e,r) bcopy (e, r, 2*NE)
+#define GET_REAL(r,e) bcopy ((char *) r, (char *) e, 2*NE)
+#define PUT_REAL(e,r) bcopy ((char *) e, (char *) r, 2*NE)
#else
#define NE 6
#define MAXDECEXP 4932
#define MINDECEXP -4956
#ifdef REAL_ARITHMETIC
/* Emulator uses target format internally
- but host stores it in host endian-ness. */
-
-#if HOST_FLOAT_WORDS_BIG_ENDIAN == FLOAT_WORDS_BIG_ENDIAN
-#define GET_REAL(r,e) e53toe ((unsigned EMUSHORT*) (r), (e))
-#define PUT_REAL(e,r) etoe53 ((e), (unsigned EMUSHORT *) (r))
-
-#else /* endian-ness differs */
-/* emulator uses target endian-ness internally */
-#define GET_REAL(r,e) \
-do { unsigned EMUSHORT w[4]; \
- w[3] = ((EMUSHORT *) r)[0]; \
- w[2] = ((EMUSHORT *) r)[1]; \
- w[1] = ((EMUSHORT *) r)[2]; \
- w[0] = ((EMUSHORT *) r)[3]; \
- e53toe (w, (e)); } while (0)
-
-#define PUT_REAL(e,r) \
-do { unsigned EMUSHORT w[4]; \
- etoe53 ((e), w); \
- *((EMUSHORT *) r) = w[3]; \
- *((EMUSHORT *) r + 1) = w[2]; \
- *((EMUSHORT *) r + 2) = w[1]; \
- *((EMUSHORT *) r + 3) = w[0]; } while (0)
-
-#endif /* endian-ness differs */
+ but host stores it in host endian-ness. */
+
+#define GET_REAL(r,e) \
+do { \
+ if (HOST_FLOAT_WORDS_BIG_ENDIAN == REAL_WORDS_BIG_ENDIAN) \
+ e53toe ((unsigned EMUSHORT *) (r), (e)); \
+ else \
+ { \
+ unsigned EMUSHORT w[4]; \
+ w[3] = ((EMUSHORT *) r)[0]; \
+ w[2] = ((EMUSHORT *) r)[1]; \
+ w[1] = ((EMUSHORT *) r)[2]; \
+ w[0] = ((EMUSHORT *) r)[3]; \
+ e53toe (w, (e)); \
+ } \
+ } while (0)
+
+#define PUT_REAL(e,r) \
+do { \
+ if (HOST_FLOAT_WORDS_BIG_ENDIAN == REAL_WORDS_BIG_ENDIAN) \
+ etoe53 ((e), (unsigned EMUSHORT *) (r)); \
+ else \
+ { \
+ unsigned EMUSHORT w[4]; \
+ etoe53 ((e), w); \
+ *((EMUSHORT *) r) = w[3]; \
+ *((EMUSHORT *) r + 1) = w[2]; \
+ *((EMUSHORT *) r + 2) = w[1]; \
+ *((EMUSHORT *) r + 3) = w[0]; \
+ } \
+ } while (0)
#else /* not REAL_ARITHMETIC */
unsigned EMUSHORT *));
static void eiremain PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
static void mtherr PROTO((char *, int));
+#ifdef DEC
static void dectoe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
static void etodec PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
static void todec PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
+#endif
+#ifdef IBM
static void ibmtoe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *,
enum machine_mode));
static void etoibm PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *,
enum machine_mode));
static void toibm PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *,
enum machine_mode));
+#endif
static void make_nan PROTO((unsigned EMUSHORT *, int, enum machine_mode));
static void uditoe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
static void ditoe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *));
{
unsigned long th, t;
-#if FLOAT_WORDS_BIG_ENDIAN
- switch (mode)
+ if (REAL_WORDS_BIG_ENDIAN)
{
+ switch (mode)
+ {
- case TFmode:
- /* Swap halfwords in the fourth long. */
- th = (unsigned long) e[6] & 0xffff;
- t = (unsigned long) e[7] & 0xffff;
- t |= th << 16;
- x[3] = (long) t;
-
- case XFmode:
-
- /* Swap halfwords in the third long. */
- th = (unsigned long) e[4] & 0xffff;
- t = (unsigned long) e[5] & 0xffff;
- t |= th << 16;
- x[2] = (long) t;
- /* fall into the double case */
-
- case DFmode:
-
- /* swap halfwords in the second word */
- th = (unsigned long) e[2] & 0xffff;
- t = (unsigned long) e[3] & 0xffff;
- t |= th << 16;
- x[1] = (long) t;
- /* fall into the float case */
-
- case HFmode:
- case SFmode:
-
- /* swap halfwords in the first word */
- th = (unsigned long) e[0] & 0xffff;
- t = (unsigned long) e[1] & 0xffff;
- t |= th << 16;
- x[0] = t;
- break;
+ case TFmode:
+ /* Swap halfwords in the fourth long. */
+ th = (unsigned long) e[6] & 0xffff;
+ t = (unsigned long) e[7] & 0xffff;
+ t |= th << 16;
+ x[3] = (long) t;
+
+ case XFmode:
+
+ /* Swap halfwords in the third long. */
+ th = (unsigned long) e[4] & 0xffff;
+ t = (unsigned long) e[5] & 0xffff;
+ t |= th << 16;
+ x[2] = (long) t;
+ /* fall into the double case */
+
+ case DFmode:
+
+ /* swap halfwords in the second word */
+ th = (unsigned long) e[2] & 0xffff;
+ t = (unsigned long) e[3] & 0xffff;
+ t |= th << 16;
+ x[1] = (long) t;
+ /* fall into the float case */
+
+ case HFmode:
+ case SFmode:
+
+ /* swap halfwords in the first word */
+ th = (unsigned long) e[0] & 0xffff;
+ t = (unsigned long) e[1] & 0xffff;
+ t |= th << 16;
+ x[0] = (long) t;
+ break;
- default:
- abort ();
+ default:
+ abort ();
+ }
}
-
-#else
-
- /* Pack the output array without swapping. */
-
- switch (mode)
+ else
{
+ /* Pack the output array without swapping. */
- case TFmode:
-
- /* Pack the fourth long. */
- th = (unsigned long) e[7] & 0xffff;
- t = (unsigned long) e[6] & 0xffff;
- t |= th << 16;
- x[3] = (long) t;
-
- case XFmode:
-
- /* Pack the third long.
- Each element of the input REAL_VALUE_TYPE array has 16 useful bits
- in it. */
- th = (unsigned long) e[5] & 0xffff;
- t = (unsigned long) e[4] & 0xffff;
- t |= th << 16;
- x[2] = (long) t;
- /* fall into the double case */
-
- case DFmode:
-
- /* pack the second long */
- th = (unsigned long) e[3] & 0xffff;
- t = (unsigned long) e[2] & 0xffff;
- t |= th << 16;
- x[1] = (long) t;
- /* fall into the float case */
-
- case HFmode:
- case SFmode:
+ switch (mode)
+ {
- /* pack the first long */
- th = (unsigned long) e[1] & 0xffff;
- t = (unsigned long) e[0] & 0xffff;
- t |= th << 16;
- x[0] = t;
- break;
+ case TFmode:
+
+ /* Pack the fourth long. */
+ th = (unsigned long) e[7] & 0xffff;
+ t = (unsigned long) e[6] & 0xffff;
+ t |= th << 16;
+ x[3] = (long) t;
+
+ case XFmode:
+
+ /* Pack the third long.
+ Each element of the input REAL_VALUE_TYPE array has 16 useful bits
+ in it. */
+ th = (unsigned long) e[5] & 0xffff;
+ t = (unsigned long) e[4] & 0xffff;
+ t |= th << 16;
+ x[2] = (long) t;
+ /* fall into the double case */
+
+ case DFmode:
+
+ /* pack the second long */
+ th = (unsigned long) e[3] & 0xffff;
+ t = (unsigned long) e[2] & 0xffff;
+ t |= th << 16;
+ x[1] = (long) t;
+ /* fall into the float case */
+
+ case HFmode:
+ case SFmode:
+
+ /* pack the first long */
+ th = (unsigned long) e[1] & 0xffff;
+ t = (unsigned long) e[0] & 0xffff;
+ t |= th << 16;
+ x[0] = (long) t;
+ break;
- default:
- abort ();
+ default:
+ abort ();
+ }
}
-
-#endif
}
GET_REAL (r1, d1);
GET_REAL (r2, d2);
#ifdef NANS
-/* Return NaN input back to the caller. */
+/* Return NaN input back to the caller. */
if (eisnan (d1))
{
PUT_REAL (d1, value);
/* REAL_VALUE_FROM_INT macro. */
void
-ereal_from_int (d, i, j)
+ereal_from_int (d, i, j, mode)
REAL_VALUE_TYPE *d;
HOST_WIDE_INT i, j;
+ enum machine_mode mode;
{
unsigned EMUSHORT df[NE], dg[NE];
HOST_WIDE_INT low, high;
int sign;
+ if (GET_MODE_CLASS (mode) != MODE_FLOAT)
+ abort ();
sign = 0;
low = i;
if ((high = j) < 0)
eadd (df, dg, dg);
if (sign)
eneg (dg);
+
+ /* A REAL_VALUE_TYPE may not be wide enough to hold the two HOST_WIDE_INTS.
+ Avoid double-rounding errors later by rounding off now from the
+ extra-wide internal format to the requested precision. */
+ switch (GET_MODE_BITSIZE (mode))
+ {
+ case 32:
+ etoe24 (dg, df);
+ e24toe (df, dg);
+ break;
+
+ case 64:
+ etoe53 (dg, df);
+ e53toe (df, dg);
+ break;
+
+ case 96:
+ etoe64 (dg, df);
+ e64toe (df, dg);
+ break;
+
+ case 128:
+ etoe113 (dg, df);
+ e113toe (df, dg);
+ break;
+
+ default:
+ abort ();
+ }
+
PUT_REAL (dg, d);
}
/* REAL_VALUE_FROM_UNSIGNED_INT macro. */
void
-ereal_from_uint (d, i, j)
+ereal_from_uint (d, i, j, mode)
REAL_VALUE_TYPE *d;
unsigned HOST_WIDE_INT i, j;
+ enum machine_mode mode;
{
unsigned EMUSHORT df[NE], dg[NE];
unsigned HOST_WIDE_INT low, high;
+ if (GET_MODE_CLASS (mode) != MODE_FLOAT)
+ abort ();
low = i;
high = j;
eldexp (eone, HOST_BITS_PER_WIDE_INT, df);
emul (dg, df, dg);
ultoe (&low, df);
eadd (df, dg, dg);
+
+ /* A REAL_VALUE_TYPE may not be wide enough to hold the two HOST_WIDE_INTS.
+ Avoid double-rounding errors later by rounding off now from the
+ extra-wide internal format to the requested precision. */
+ switch (GET_MODE_BITSIZE (mode))
+ {
+ case 32:
+ etoe24 (dg, df);
+ e24toe (df, dg);
+ break;
+
+ case 64:
+ etoe53 (dg, df);
+ e53toe (df, dg);
+ break;
+
+ case 96:
+ etoe64 (dg, df);
+ e64toe (df, dg);
+ break;
+
+ case 128:
+ etoe113 (dg, df);
+ e113toe (df, dg);
+ break;
+
+ default:
+ abort ();
+ }
+
PUT_REAL (dg, d);
}
#ifdef REAL_ARITHMETIC
-/* Check for infinity in a REAL_VALUE_TYPE. */
+/* Check for infinity in a REAL_VALUE_TYPE. */
int
target_isinf (x)
#endif
}
-
-/* Check whether a REAL_VALUE_TYPE item is a NaN. */
+/* Check whether a REAL_VALUE_TYPE item is a NaN. */
int
target_isnan (x)
/* Check for a negative REAL_VALUE_TYPE number.
- This just checks the sign bit, so that -0 counts as negative. */
+ This just checks the sign bit, so that -0 counts as negative. */
int
target_negative (x)
return (r);
}
+/* Try to change R into its exact multiplicative inverse in machine mode
+ MODE. Return nonzero function value if successful. */
+
+int
+exact_real_inverse (mode, r)
+ enum machine_mode mode;
+ REAL_VALUE_TYPE *r;
+{
+ unsigned EMUSHORT e[NE], einv[NE];
+ REAL_VALUE_TYPE rinv;
+ int i;
+
+ GET_REAL (r, e);
+
+ /* Test for input in range. Don't transform IEEE special values. */
+ if (eisinf (e) || eisnan (e) || (ecmp (e, ezero) == 0))
+ return 0;
+
+ /* Test for a power of 2: all significand bits zero except the MSB.
+ We are assuming the target has binary (or hex) arithmetic. */
+ if (e[NE - 2] != 0x8000)
+ return 0;
+
+ for (i = 0; i < NE - 2; i++)
+ {
+ if (e[i] != 0)
+ return 0;
+ }
+
+ /* Compute the inverse and truncate it to the required mode. */
+ ediv (e, eone, einv);
+ PUT_REAL (einv, &rinv);
+ rinv = real_value_truncate (mode, rinv);
+
+#ifdef CHECK_FLOAT_VALUE
+ /* This check is not redundant. It may, for example, flush
+ a supposedly IEEE denormal value to zero. */
+ i = 0;
+ if (CHECK_FLOAT_VALUE (mode, rinv, i))
+ return 0;
+#endif
+ GET_REAL (&rinv, einv);
+
+ /* Check the bits again, because the truncation might have
+ generated an arbitrary saturation value on overflow. */
+ if (einv[NE - 2] != 0x8000)
+ return 0;
+
+ for (i = 0; i < NE - 2; i++)
+ {
+ if (einv[i] != 0)
+ return 0;
+ }
+
+ /* Fail if the computed inverse is out of range. */
+ if (eisinf (einv) || eisnan (einv) || (ecmp (einv, ezero) == 0))
+ return 0;
+
+ /* Output the reciprocal and return success flag. */
+ PUT_REAL (einv, r);
+ return 1;
+}
#endif /* REAL_ARITHMETIC defined */
/* Used for debugging--print the value of R in human-readable format
}
\f
-/* Target values are arrays of host longs. A long is guaranteed
- to be at least 32 bits wide. */
+/* The following routines convert REAL_VALUE_TYPE to the various floating
+ point formats that are meaningful to supported computers.
+
+ The results are returned in 32-bit pieces, each piece stored in a `long'.
+ This is so they can be printed by statements like
+
+ fprintf (file, "%lx, %lx", L[0], L[1]);
-/* 128-bit long double */
+ that will work on both narrow- and wide-word host computers. */
+
+/* Convert R to a 128-bit long double precision value. The output array L
+ contains four 32-bit pieces of the result, in the order they would appear
+ in memory. */
void
etartdouble (r, l)
endian (e, l, TFmode);
}
-/* 80-bit long double */
+/* Convert R to a double extended precision value. The output array L
+ contains three 32-bit pieces of the result, in the order they would
+ appear in memory. */
void
etarldouble (r, l)
endian (e, l, XFmode);
}
+/* Convert R to a double precision value. The output array L contains two
+ 32-bit pieces of the result, in the order they would appear in memory. */
+
void
etardouble (r, l)
REAL_VALUE_TYPE r;
endian (e, l, DFmode);
}
+/* Convert R to a single precision float value stored in the least-significant
+ bits of a `long'. */
+
long
etarsingle (r)
REAL_VALUE_TYPE r;
return ((long) l);
}
+/* Convert X to a decimal ASCII string S for output to an assembly
+ language file. Note, there is no standard way to spell infinity or
+ a NaN, so these values may require special treatment in the tm.h
+ macros. */
+
void
ereal_to_decimal (x, s)
REAL_VALUE_TYPE x;
etoasc (e, s, 20);
}
+/* Compare X and Y. Return 1 if X > Y, 0 if X == Y, -1 if X < Y,
+ or -2 if either is a NaN. */
+
int
ereal_cmp (x, y)
REAL_VALUE_TYPE x, y;
return (ecmp (ex, ey));
}
+/* Return 1 if the sign bit of X is set, else return 0. */
+
int
ereal_isneg (x)
REAL_VALUE_TYPE x;
short integers. The arguments of the routines are pointers to
the arrays.
- External e type data structure, simulates Intel 8087 chip
+ External e type data structure, similar to Intel 8087 chip
temporary real format but possibly with a larger significand:
NE-1 significand words (least significant word first,
top bit is the sign)
- Internal data structure of a number (a "word" is 16 bits):
+ Internal exploded e-type data structure of a number (a "word" is 16 bits):
ei[0] sign word (0 for positive, 0xffff for negative)
ei[1] biased exponent (value = EXONE for the number 1.0)
- Routines for external format numbers
+ Routines for external format e-type numbers
asctoe (string, e) ASCII string to extended double e type
asctoe64 (string, &d) ASCII string to long double
eisnan (e) 1 if e is a NaN
- Routines for internal format numbers
+ Routines for internal format exploded e-type numbers
eaddm (ai, bi) add significands, bi = bi + ai
ecleaz (ei) ei = 0
For computers, such as IBM PC, that follow the IEEE
Standard for Binary Floating Point Arithmetic (ANSI/IEEE
- Std 754-1985), the symbol IBMPC or MIEEE should be defined.
+ Std 754-1985), the symbol IEEE should be defined.
These numbers have 53-bit significands. In this mode, constants
are provided as arrays of hexadecimal 16 bit integers.
+ The endian-ness of generated values is controlled by
+ REAL_WORDS_BIG_ENDIAN.
To accommodate other types of computer arithmetic, all
constants are also provided in a normal decimal radix
ensure that these values are correct for your computer.
For ANSI C compatibility, define ANSIC equal to 1. Currently
- this affects only the atan2 function and others that use it. */
+ this affects only the atan2 function and others that use it. */
/* Constant definitions for math error conditions. */
{0xc4c6, 0xc234, 0020550, 0155242, 0144417, 0040000,};
#endif
-
-
/* Control register for rounding precision.
This can be set to 113 (if NE=10), 80 (if NE=6), 64, 56, 53, or 24 bits. */
int rndprc = NBITS;
extern int rndprc;
-/* Clear out entire external format number. */
+/* Clear out entire e-type number X. */
static void
eclear (x)
*x++ = 0;
}
-
-
-/* Move external format number from a to b. */
+/* Move e-type number from A to B. */
static void
emov (a, b)
}
-/* Absolute value of external format number. */
+/* Absolute value of e-type X. */
static void
eabs (x)
x[NE - 1] &= 0x7fff;
}
-/* Negate external format number. */
+/* Negate the e-type number X. */
static void
eneg (x)
x[NE - 1] ^= 0x8000; /* Toggle the sign bit */
}
-
-
-/* Return 1 if sign bit of external format number is nonzero, else zero. */
+/* Return 1 if sign bit of e-type number X is nonzero, else zero. */
static int
eisneg (x)
return (0);
}
-
-/* Return 1 if external format number is infinity, else return zero. */
+/* Return 1 if e-type number X is infinity, else return zero. */
static int
eisinf (x)
return (0);
}
-
/* Check if e-type number is not a number. The bit pattern is one that we
defined, so we know for sure how to detect it. */
/* NaN has maximum exponent */
if ((x[NE - 1] & 0x7fff) != 0x7fff)
return (0);
- /* ... and non-zero significand field. */
+ /* ... and non-zero significand field. */
for (i = 0; i < NE - 1; i++)
{
if (*x++ != 0)
return (0);
}
-/* Fill external format number with infinity pattern (IEEE)
- or largest possible number (non-IEEE). */
+/* Fill e-type number X with infinity pattern (IEEE)
+ or largest possible number (non-IEEE). */
static void
einfin (x)
#endif
}
-
/* Output an e-type NaN.
This generates Intel's quiet NaN pattern for extended real.
The exponent is 7fff, the leading mantissa word is c000. */
*x = (sign << 15) | 0x7fff;
}
-
-/* Move in external format number, converting it to internal format. */
+/* Move in an e-type number A, converting it to exploded e-type B. */
static void
emovi (a, b)
*q = 0;
}
-
-/* Move internal format number out, converting it to external format. */
+/* Move out exploded e-type number A, converting it to e type B. */
static void
emovo (a, b)
*q-- = *p++;
}
-/* Clear out internal format number. */
+/* Clear out exploded e-type number XI. */
static void
ecleaz (xi)
*xi++ = 0;
}
-
-/* Same, but don't touch the sign. */
+/* Clear out exploded e-type XI, but don't touch the sign. */
static void
ecleazs (xi)
*xi++ = 0;
}
-
-
-/* Move internal format number from a to b. */
+/* Move exploded e-type number from A to B. */
static void
emovz (a, b)
*b = 0;
}
-/* Generate internal format NaN.
+/* Generate exploded e-type NaN.
The explicit pattern for this is maximum exponent and
top two significant bits set. */
x[M + 1] = 0xc000;
}
-/* Return nonzero if internal format number is a NaN. */
+/* Return nonzero if exploded e-type X is a NaN. */
static int
eiisnan (x)
return (0);
}
-/* Return nonzero if sign of internal format number is nonzero. */
+/* Return nonzero if sign of exploded e-type X is nonzero. */
static int
eiisneg (x)
return x[0] != 0;
}
-/* Fill internal format number with infinity pattern.
+/* Fill exploded e-type X with infinity pattern.
This has maximum exponent and significand all zeros. */
static void
x[E] = 0x7fff;
}
-/* Return nonzero if internal format number is infinite. */
+/* Return nonzero if exploded e-type X is infinite. */
static int
eiisinf (x)
}
-/* Compare significands of numbers in internal format.
+/* Compare significands of numbers in internal exploded e-type format.
Guard words are included in the comparison.
Returns +1 if a > b
return (-1);
}
-
-/* Shift significand down by 1 bit. */
+/* Shift significand of exploded e-type X down by 1 bit. */
static void
eshdn1 (x)
}
}
-
-
-/* Shift significand up by 1 bit. */
+/* Shift significand of exploded e-type X up by 1 bit. */
static void
eshup1 (x)
}
-/* Shift significand down by 8 bits. */
+/* Shift significand of exploded e-type X down by 8 bits. */
static void
eshdn8 (x)
}
}
-/* Shift significand up by 8 bits. */
+/* Shift significand of exploded e-type X up by 8 bits. */
static void
eshup8 (x)
}
}
-/* Shift significand up by 16 bits. */
+/* Shift significand of exploded e-type X up by 16 bits. */
static void
eshup6 (x)
*p = 0;
}
-/* Shift significand down by 16 bits. */
+/* Shift significand of exploded e-type X down by 16 bits. */
static void
eshdn6 (x)
*(--p) = 0;
}
-\f
-/* Add significands. x + y replaces y. */
+
+/* Add significands of exploded e-type X and Y. X + Y replaces Y. */
static void
eaddm (x, y)
}
}
-/* Subtract significands. y - x replaces y. */
+/* Subtract significands of exploded e-type X and Y. Y - X replaces Y. */
static void
esubm (x, y)
/* Multiply significands */
+
int
emulm (a, b)
unsigned EMUSHORT a[], b[];
#else
-/* Radix 65536 versions of multiply and divide */
+/* Radix 65536 versions of multiply and divide. */
-
-/* Multiply significand of e-type number b
- by 16-bit quantity a, e-type result to c. */
+/* Multiply significand of e-type number B
+ by 16-bit quantity A, return e-type result to C. */
static void
m16m (a, b, c)
unsigned int a;
- unsigned short b[], c[];
+ unsigned EMUSHORT b[], c[];
{
- register unsigned short *pp;
- register unsigned long carry;
- unsigned short *ps;
- unsigned short p[NI];
- unsigned long aa, m;
+ register unsigned EMUSHORT *pp;
+ register unsigned EMULONG carry;
+ unsigned EMUSHORT *ps;
+ unsigned EMUSHORT p[NI];
+ unsigned EMULONG aa, m;
int i;
aa = a;
}
else
{
- m = (unsigned long) aa * *ps--;
+ m = (unsigned EMULONG) aa * *ps--;
carry = (m & 0xffff) + *pp;
- *pp-- = (unsigned short)carry;
+ *pp-- = (unsigned EMUSHORT)carry;
carry = (carry >> 16) + (m >> 16) + *pp;
- *pp = (unsigned short)carry;
+ *pp = (unsigned EMUSHORT)carry;
*(pp-1) = carry >> 16;
}
}
c[i] = p[i];
}
-
-/* Divide significands. Neither the numerator nor the denominator
- is permitted to have its high guard word nonzero. */
+/* Divide significands of exploded e-types NUM / DEN. Neither the
+ numerator NUM nor the denominator DEN is permitted to have its high guard
+ word nonzero. */
static int
edivm (den, num)
- unsigned short den[], num[];
+ unsigned EMUSHORT den[], num[];
{
int i;
- register unsigned short *p;
- unsigned long tnum;
- unsigned short j, tdenm, tquot;
- unsigned short tprod[NI+1];
+ register unsigned EMUSHORT *p;
+ unsigned EMULONG tnum;
+ unsigned EMUSHORT j, tdenm, tquot;
+ unsigned EMUSHORT tprod[NI+1];
p = &equot[0];
*p++ = num[0];
tdenm = den[M+1];
for (i=M; i<NI; i++)
{
- /* Find trial quotient digit (the radix is 65536). */
- tnum = (((unsigned long) num[M]) << 16) + num[M+1];
+ /* Find trial quotient digit (the radix is 65536). */
+ tnum = (((unsigned EMULONG) num[M]) << 16) + num[M+1];
- /* Do not execute the divide instruction if it will overflow. */
+ /* Do not execute the divide instruction if it will overflow. */
if ((tdenm * 0xffffL) < tnum)
tquot = 0xffff;
else
tquot = tnum / tdenm;
- /* Multiply denominator by trial quotient digit. */
+ /* Multiply denominator by trial quotient digit. */
m16m ((unsigned int)tquot, den, tprod);
- /* The quotient digit may have been overestimated. */
+ /* The quotient digit may have been overestimated. */
if (ecmpm (tprod, num) > 0)
{
tquot -= 1;
return ((int)j);
}
+/* Multiply significands of exploded e-type A and B, result in B. */
-
-/* Multiply significands */
static int
emulm (a, b)
- unsigned short a[], b[];
+ unsigned EMUSHORT a[], b[];
{
- unsigned short *p, *q;
- unsigned short pprod[NI];
- unsigned short j;
+ unsigned EMUSHORT *p, *q;
+ unsigned EMUSHORT pprod[NI];
+ unsigned EMUSHORT j;
int i;
equot[0] = b[0];
/* Normalize and round off.
- The internal format number to be rounded is "s".
- Input "lost" indicates whether or not the number is exact.
- This is the so-called sticky bit.
+ The internal format number to be rounded is S.
+ Input LOST is 0 if the value is exact. This is the so-called sticky bit.
- Input "subflg" indicates whether the number was obtained
- by a subtraction operation. In that case if lost is nonzero
+ Input SUBFLG indicates whether the number was obtained
+ by a subtraction operation. In that case if LOST is nonzero
then the number is slightly smaller than indicated.
- Input "exp" is the biased exponent, which may be negative.
- the exponent field of "s" is ignored but is replaced by
- "exp" as adjusted by normalization and rounding.
+ Input EXP is the biased exponent, which may be negative.
+ the exponent field of S is ignored but is replaced by
+ EXP as adjusted by normalization and rounding.
- Input "rcntrl" is the rounding control.
+ Input RCNTRL is the rounding control. If it is nonzero, the
+ returned value will be rounded to RNDPRC bits.
For future reference: In order for emdnorm to round off denormal
significands at the right point, the input exponent must be
/* Normalize */
j = enormlz (s);
- /* a blank significand could mean either zero or infinity. */
+ /* a blank significand could mean either zero or infinity. */
#ifndef INFINITY
if (j > NBITS)
{
return;
}
}
- /* Round off, unless told not to by rcntrl. */
+ /* Round off, unless told not to by rcntrl. */
if (rcntrl == 0)
goto mdfin;
- /* Set up rounding parameters if the control register changed. */
+ /* Set up rounding parameters if the control register changed. */
if (rndprc != rlast)
{
ecleaz (rbit);
}
/* Shift down 1 temporarily if the data structure has an implied
- most significant bit and the number is denormal. */
- if ((exp <= 0) && (rndprc != 64) && (rndprc != NBITS))
+ most significant bit and the number is denormal.
+ Intel long double denormals also lose one bit of precision. */
+ if ((exp <= 0) && (rndprc != NBITS)
+ && ((rndprc != 64) || ((rndprc == 64) && ! REAL_WORDS_BIG_ENDIAN)))
{
lost |= s[NI - 1] & 1;
eshdn1 (s);
eaddm (rbit, s);
}
mddone:
- if ((exp <= 0) && (rndprc != 64) && (rndprc != NBITS))
+/* Undo the temporary shift for denormal values. */
+ if ((exp <= 0) && (rndprc != NBITS)
+ && ((rndprc != 64) || ((rndprc == 64) && ! REAL_WORDS_BIG_ENDIAN)))
{
eshup1 (s);
}
s[1] = (unsigned EMUSHORT) exp;
}
-
-
-/* Subtract external format numbers. */
+/* Subtract. C = B - A, all e type numbers. */
static int subflg = 0;
return;
}
/* Infinity minus infinity is a NaN.
- Test for subtracting infinities of the same sign. */
+ Test for subtracting infinities of the same sign. */
if (eisinf (a) && eisinf (b)
&& ((eisneg (a) ^ eisneg (b)) == 0))
{
eadd1 (a, b, c);
}
-
-/* Add. */
+/* Add. C = A + B, all e type. */
static void
eadd (a, b, c)
{
#ifdef NANS
-/* NaN plus anything is a NaN. */
+/* NaN plus anything is a NaN. */
if (eisnan (a))
{
emov (a, c);
return;
}
/* Infinity minus infinity is a NaN.
- Test for adding infinities of opposite signs. */
+ Test for adding infinities of opposite signs. */
if (eisinf (a) && eisinf (b)
&& ((eisneg (a) ^ eisneg (b)) != 0))
{
eadd1 (a, b, c);
}
+/* Arithmetic common to both addition and subtraction. */
+
static void
eadd1 (a, b, c)
unsigned EMUSHORT *a, *b, *c;
return;
}
/* if same sign, result is double */
- /* double denomalized tiny number */
+ /* double denormalized tiny number */
if ((bi[E] == 0) && ((bi[3] & 0x8000) == 0))
{
eshup1 (bi);
{
if (bi[j] != 0)
{
- /* This could overflow, but let emovo take care of that. */
ltb += 1;
+ if (ltb >= 0x7fff)
+ {
+ eclear (c);
+ if (ai[0] != 0)
+ eneg (c);
+ einfin (c);
+ return;
+ }
break;
}
}
emovo (bi, c);
}
-
-
-/* Divide. */
+/* Divide: C = B/A, all e type. */
static void
ediv (a, b, c)
unsigned EMUSHORT *a, *b, *c;
{
unsigned EMUSHORT ai[NI], bi[NI];
- int i;
+ int i, sign;
EMULONG lt, lta, ltb;
+/* IEEE says if result is not a NaN, the sign is "-" if and only if
+ operands have opposite signs -- but flush -0 to 0 later if not IEEE. */
+ sign = eisneg(a) ^ eisneg(b);
+
#ifdef NANS
-/* Return any NaN input. */
+/* Return any NaN input. */
if (eisnan (a))
{
emov (a, c);
emov (b, c);
return;
}
-/* Zero over zero, or infinity over infinity, is a NaN. */
+/* Zero over zero, or infinity over infinity, is a NaN. */
if (((ecmp (a, ezero) == 0) && (ecmp (b, ezero) == 0))
|| (eisinf (a) && eisinf (b)))
{
mtherr ("ediv", INVALID);
- enan (c, eisneg (a) ^ eisneg (b));
+ enan (c, sign);
return;
}
#endif
-/* Infinity over anything else is infinity. */
+/* Infinity over anything else is infinity. */
#ifdef INFINITY
if (eisinf (b))
{
- if (eisneg (a) ^ eisneg (b))
- *(c + (NE - 1)) = 0x8000;
- else
- *(c + (NE - 1)) = 0;
einfin (c);
- return;
+ goto divsign;
}
-/* Anything else over infinity is zero. */
+/* Anything else over infinity is zero. */
if (eisinf (a))
{
eclear (c);
- return;
+ goto divsign;
}
#endif
emovi (a, ai);
lta = ai[E];
ltb = bi[E];
if (bi[E] == 0)
- { /* See if numerator is zero. */
+ { /* See if numerator is zero. */
for (i = 1; i < NI - 1; i++)
{
if (bi[i] != 0)
}
}
eclear (c);
- return;
+ goto divsign;
}
dnzro1:
goto dnzro2;
}
}
- if (ai[0] == bi[0])
- *(c + (NE - 1)) = 0;
- else
- *(c + (NE - 1)) = 0x8000;
/* Divide by zero is not an invalid operation.
It is a divide-by-zero operation! */
einfin (c);
mtherr ("ediv", SING);
- return;
+ goto divsign;
}
dnzro2:
/* calculate exponent */
lt = ltb - lta + EXONE;
emdnorm (bi, i, 0, lt, 64);
- /* set the sign */
- if (ai[0] == bi[0])
- bi[0] = 0;
- else
- bi[0] = 0Xffff;
emovo (bi, c);
-}
+ divsign:
+ if (sign
+#ifndef IEEE
+ && (ecmp (c, ezero) != 0)
+#endif
+ )
+ *(c+(NE-1)) |= 0x8000;
+ else
+ *(c+(NE-1)) &= ~0x8000;
+}
-/* Multiply. */
+/* Multiply e-types A and B, return e-type product C. */
static void
emul (a, b, c)
unsigned EMUSHORT *a, *b, *c;
{
unsigned EMUSHORT ai[NI], bi[NI];
- int i, j;
+ int i, j, sign;
EMULONG lt, lta, ltb;
+/* IEEE says if result is not a NaN, the sign is "-" if and only if
+ operands have opposite signs -- but flush -0 to 0 later if not IEEE. */
+ sign = eisneg(a) ^ eisneg(b);
+
#ifdef NANS
-/* NaN times anything is the same NaN. */
+/* NaN times anything is the same NaN. */
if (eisnan (a))
{
emov (a, c);
emov (b, c);
return;
}
-/* Zero times infinity is a NaN. */
+/* Zero times infinity is a NaN. */
if ((eisinf (a) && (ecmp (b, ezero) == 0))
|| (eisinf (b) && (ecmp (a, ezero) == 0)))
{
mtherr ("emul", INVALID);
- enan (c, eisneg (a) ^ eisneg (b));
+ enan (c, sign);
return;
}
#endif
-/* Infinity times anything else is infinity. */
+/* Infinity times anything else is infinity. */
#ifdef INFINITY
if (eisinf (a) || eisinf (b))
{
- if (eisneg (a) ^ eisneg (b))
- *(c + (NE - 1)) = 0x8000;
- else
- *(c + (NE - 1)) = 0;
einfin (c);
- return;
+ goto mulsign;
}
#endif
emovi (a, ai);
}
}
eclear (c);
- return;
+ goto mulsign;
}
mnzer1:
}
}
eclear (c);
- return;
+ goto mulsign;
}
mnzer2:
/* calculate exponent */
lt = lta + ltb - (EXONE - 1);
emdnorm (bi, j, 0, lt, 64);
- /* calculate sign of product */
- if (ai[0] == bi[0])
- bi[0] = 0;
- else
- bi[0] = 0xffff;
emovo (bi, c);
-}
-
+ mulsign:
+ if (sign
+#ifndef IEEE
+ && (ecmp (c, ezero) != 0)
+#endif
+ )
+ *(c+(NE-1)) |= 0x8000;
+ else
+ *(c+(NE-1)) &= ~0x8000;
+}
-/* Convert IEEE double precision to e type. */
+/* Convert double precision PE to e-type Y. */
static void
e53toe (pe, y)
{
#ifdef DEC
- dectoe (pe, y); /* see etodec.c */
+ dectoe (pe, y);
#else
#ifdef IBM
e = pe;
denorm = 0; /* flag if denormalized number */
ecleaz (yy);
-#ifdef IBMPC
- e += 3;
-#endif
+ if (! REAL_WORDS_BIG_ENDIAN)
+ e += 3;
r = *e;
yy[0] = 0;
if (r & 0x8000)
if (r == 0x7ff0)
{
#ifdef NANS
-#ifdef IBMPC
- if (((pe[3] & 0xf) != 0) || (pe[2] != 0)
- || (pe[1] != 0) || (pe[0] != 0))
+ if (! REAL_WORDS_BIG_ENDIAN)
{
- enan (y, yy[0] != 0);
- return;
+ if (((pe[3] & 0xf) != 0) || (pe[2] != 0)
+ || (pe[1] != 0) || (pe[0] != 0))
+ {
+ enan (y, yy[0] != 0);
+ return;
+ }
}
-#else
- if (((pe[0] & 0xf) != 0) || (pe[1] != 0)
- || (pe[2] != 0) || (pe[3] != 0))
+ else
{
- enan (y, yy[0] != 0);
- return;
+ if (((pe[0] & 0xf) != 0) || (pe[1] != 0)
+ || (pe[2] != 0) || (pe[3] != 0))
+ {
+ enan (y, yy[0] != 0);
+ return;
+ }
}
-#endif
#endif /* NANS */
eclear (y);
einfin (y);
#endif /* INFINITY */
r >>= 4;
/* If zero exponent, then the significand is denormalized.
- So take back the understood high significand bit. */
+ So take back the understood high significand bit. */
if (r == 0)
{
r += EXONE - 01777;
yy[E] = r;
p = &yy[M + 1];
-#ifdef IBMPC
- *p++ = *(--e);
- *p++ = *(--e);
- *p++ = *(--e);
-#endif
-#ifdef MIEEE
- ++e;
- *p++ = *e++;
- *p++ = *e++;
- *p++ = *e++;
+#ifdef IEEE
+ if (! REAL_WORDS_BIG_ENDIAN)
+ {
+ *p++ = *(--e);
+ *p++ = *(--e);
+ *p++ = *(--e);
+ }
+ else
+ {
+ ++e;
+ *p++ = *e++;
+ *p++ = *e++;
+ *p++ = *e++;
+ }
#endif
eshift (yy, -5);
if (denorm)
#endif /* not DEC */
}
+/* Convert double extended precision float PE to e type Y. */
+
static void
e64toe (pe, y)
unsigned EMUSHORT *pe, *y;
p = yy;
for (i = 0; i < NE - 5; i++)
*p++ = 0;
-#ifdef IBMPC
- for (i = 0; i < 5; i++)
- *p++ = *e++;
-#endif
-/* This precision is not ordinarily supported on DEC or IBM. */
+/* This precision is not ordinarily supported on DEC or IBM. */
#ifdef DEC
for (i = 0; i < 5; i++)
*p++ = *e++;
for (i = 0; i < 5; i++)
*p-- = *e++;
#endif
-#ifdef MIEEE
- p = &yy[0] + (NE - 1);
- *p-- = *e++;
- ++e;
- for (i = 0; i < 4; i++)
- *p-- = *e++;
+#ifdef IEEE
+ if (! REAL_WORDS_BIG_ENDIAN)
+ {
+ for (i = 0; i < 5; i++)
+ *p++ = *e++;
+
+ /* For denormal long double Intel format, shift significand up one
+ -- but only if the top significand bit is zero. A top bit of 1
+ is "pseudodenormal" when the exponent is zero. */
+ if((yy[NE-1] & 0x7fff) == 0 && (yy[NE-2] & 0x8000) == 0)
+ {
+ unsigned EMUSHORT temp[NI];
+
+ emovi(yy, temp);
+ eshup1(temp);
+ emovo(temp,y);
+ return;
+ }
+ }
+ else
+ {
+ p = &yy[0] + (NE - 1);
+#ifdef ARM_EXTENDED_IEEE_FORMAT
+ /* For ARMs, the exponent is in the lowest 15 bits of the word. */
+ *p-- = (e[0] & 0x8000) | (e[1] & 0x7ffff);
+ e += 2;
+#else
+ *p-- = *e++;
+ ++e;
+#endif
+ for (i = 0; i < 4; i++)
+ *p-- = *e++;
+ }
#endif
- p = yy;
- q = y;
#ifdef INFINITY
- if (*p == 0x7fff)
+ /* Point to the exponent field and check max exponent cases. */
+ p = &yy[NE - 1];
+ if ((*p & 0x7fff) == 0x7fff)
{
#ifdef NANS
-#ifdef IBMPC
- for (i = 0; i < 4; i++)
+ if (! REAL_WORDS_BIG_ENDIAN)
{
- if (pe[i] != 0)
+ for (i = 0; i < 4; i++)
{
- enan (y, (*p & 0x8000) != 0);
- return;
+ if ((i != 3 && pe[i] != 0)
+ /* Anything but 0x8000 here, including 0, is a NaN. */
+ || (i == 3 && pe[i] != 0x8000))
+ {
+ enan (y, (*p & 0x8000) != 0);
+ return;
+ }
}
}
-#else
- for (i = 1; i <= 4; i++)
+ else
{
- if (pe[i] != 0)
+#ifdef ARM_EXTENDED_IEEE_FORMAT
+ for (i = 2; i <= 5; i++)
{
- enan (y, (*p & 0x8000) != 0);
- return;
+ if (pe[i] != 0)
+ {
+ enan (y, (*p & 0x8000) != 0);
+ return;
+ }
+ }
+#else /* not ARM */
+ /* In Motorola extended precision format, the most significant
+ bit of an infinity mantissa could be either 1 or 0. It is
+ the lower order bits that tell whether the value is a NaN. */
+ if ((pe[2] & 0x7fff) != 0)
+ goto bigend_nan;
+
+ for (i = 3; i <= 5; i++)
+ {
+ if (pe[i] != 0)
+ {
+bigend_nan:
+ enan (y, (*p & 0x8000) != 0);
+ return;
+ }
}
+#endif /* not ARM */
}
-#endif
#endif /* NANS */
eclear (y);
einfin (y);
return;
}
#endif /* INFINITY */
+ p = yy;
+ q = y;
for (i = 0; i < NE; i++)
*q++ = *p++;
}
+/* Convert 128-bit long double precision float PE to e type Y. */
static void
e113toe (pe, y)
e = pe;
denorm = 0;
ecleaz (yy);
-#ifdef IBMPC
- e += 7;
+#ifdef IEEE
+ if (! REAL_WORDS_BIG_ENDIAN)
+ e += 7;
#endif
r = *e;
yy[0] = 0;
if (r == 0x7fff)
{
#ifdef NANS
-#ifdef IBMPC
- for (i = 0; i < 7; i++)
+ if (! REAL_WORDS_BIG_ENDIAN)
{
- if (pe[i] != 0)
+ for (i = 0; i < 7; i++)
{
- enan (y, yy[0] != 0);
- return;
- }
+ if (pe[i] != 0)
+ {
+ enan (y, yy[0] != 0);
+ return;
+ }
+ }
}
-#else
- for (i = 1; i < 8; i++)
+ else
{
- if (pe[i] != 0)
+ for (i = 1; i < 8; i++)
{
- enan (y, yy[0] != 0);
- return;
+ if (pe[i] != 0)
+ {
+ enan (y, yy[0] != 0);
+ return;
+ }
}
}
-#endif
#endif /* NANS */
eclear (y);
einfin (y);
#endif /* INFINITY */
yy[E] = r;
p = &yy[M + 1];
-#ifdef IBMPC
- for (i = 0; i < 7; i++)
- *p++ = *(--e);
-#endif
-#ifdef MIEEE
- ++e;
- for (i = 0; i < 7; i++)
- *p++ = *e++;
+#ifdef IEEE
+ if (! REAL_WORDS_BIG_ENDIAN)
+ {
+ for (i = 0; i < 7; i++)
+ *p++ = *(--e);
+ }
+ else
+ {
+ ++e;
+ for (i = 0; i < 7; i++)
+ *p++ = *e++;
+ }
#endif
-/* If denormal, remove the implied bit; else shift down 1. */
+/* If denormal, remove the implied bit; else shift down 1. */
if (r == 0)
{
yy[M] = 0;
emovo (yy, y);
}
-
-/* Convert IEEE single precision to e type. */
+/* Convert single precision float PE to e type Y. */
static void
e24toe (pe, y)
e = pe;
denorm = 0; /* flag if denormalized number */
ecleaz (yy);
-#ifdef IBMPC
- e += 1;
+#ifdef IEEE
+ if (! REAL_WORDS_BIG_ENDIAN)
+ e += 1;
#endif
#ifdef DEC
e += 1;
if (r == 0x7f80)
{
#ifdef NANS
-#ifdef MIEEE
- if (((pe[0] & 0x7f) != 0) || (pe[1] != 0))
+ if (REAL_WORDS_BIG_ENDIAN)
{
- enan (y, yy[0] != 0);
- return;
+ if (((pe[0] & 0x7f) != 0) || (pe[1] != 0))
+ {
+ enan (y, yy[0] != 0);
+ return;
+ }
}
-#else
- if (((pe[1] & 0x7f) != 0) || (pe[0] != 0))
+ else
{
- enan (y, yy[0] != 0);
- return;
+ if (((pe[1] & 0x7f) != 0) || (pe[0] != 0))
+ {
+ enan (y, yy[0] != 0);
+ return;
+ }
}
-#endif
#endif /* NANS */
eclear (y);
einfin (y);
#endif /* INFINITY */
r >>= 7;
/* If zero exponent, then the significand is denormalized.
- So take back the understood high significand bit. */
+ So take back the understood high significand bit. */
if (r == 0)
{
denorm = 1;
r += EXONE - 0177;
yy[E] = r;
p = &yy[M + 1];
-#ifdef IBMPC
- *p++ = *(--e);
-#endif
#ifdef DEC
*p++ = *(--e);
#endif
-#ifdef MIEEE
- ++e;
- *p++ = *e++;
+#ifdef IEEE
+ if (! REAL_WORDS_BIG_ENDIAN)
+ *p++ = *(--e);
+ else
+ {
+ ++e;
+ *p++ = *e++;
+ }
#endif
eshift (yy, -8);
if (denorm)
#endif /* not IBM */
}
+/* Convert e-type X to IEEE 128-bit long double format E. */
static void
etoe113 (x, e)
toe113 (xi, e);
}
-/* Move out internal format to ieee long double */
+/* Convert exploded e-type X, that has already been rounded to
+ 113-bit precision, to IEEE 128-bit long double format Y. */
static void
toe113 (a, b)
}
#endif
p = a;
-#ifdef MIEEE
- q = b;
-#else
- q = b + 7; /* point to output exponent */
-#endif
+ if (REAL_WORDS_BIG_ENDIAN)
+ q = b;
+ else
+ q = b + 7; /* point to output exponent */
- /* If not denormal, delete the implied bit. */
+ /* If not denormal, delete the implied bit. */
if (a[E] != 0)
{
eshup1 (a);
}
/* combine sign and exponent */
i = *p++;
-#ifdef MIEEE
- if (i)
- *q++ = *p++ | 0x8000;
- else
- *q++ = *p++;
-#else
- if (i)
- *q-- = *p++ | 0x8000;
+ if (REAL_WORDS_BIG_ENDIAN)
+ {
+ if (i)
+ *q++ = *p++ | 0x8000;
+ else
+ *q++ = *p++;
+ }
else
- *q-- = *p++;
-#endif
+ {
+ if (i)
+ *q-- = *p++ | 0x8000;
+ else
+ *q-- = *p++;
+ }
/* skip over guard word */
++p;
/* move the significand */
-#ifdef MIEEE
- for (i = 0; i < 7; i++)
- *q++ = *p++;
-#else
- for (i = 0; i < 7; i++)
- *q-- = *p++;
-#endif
+ if (REAL_WORDS_BIG_ENDIAN)
+ {
+ for (i = 0; i < 7; i++)
+ *q++ = *p++;
+ }
+ else
+ {
+ for (i = 0; i < 7; i++)
+ *q-- = *p++;
+ }
}
+/* Convert e-type X to IEEE double extended format E. */
+
static void
etoe64 (x, e)
unsigned EMUSHORT *x, *e;
toe64 (xi, e);
}
-
-/* Move out internal format to ieee long double. */
+/* Convert exploded e-type X, that has already been rounded to
+ 64-bit precision, to IEEE double extended format Y. */
static void
toe64 (a, b)
return;
}
#endif
+ /* Shift denormal long double Intel format significand down one bit. */
+ if ((a[E] == 0) && ! REAL_WORDS_BIG_ENDIAN)
+ eshdn1 (a);
p = a;
-#if defined(MIEEE) || defined(IBM)
+#ifdef IBM
q = b;
-#else
- q = b + 4; /* point to output exponent */
+#endif
+#ifdef DEC
+ q = b + 4;
+#endif
+#ifdef IEEE
+ if (REAL_WORDS_BIG_ENDIAN)
+ q = b;
+ else
+ {
+ q = b + 4; /* point to output exponent */
#if LONG_DOUBLE_TYPE_SIZE == 96
- /* Clear the last two bytes of 12-byte Intel format */
- *(q+1) = 0;
+ /* Clear the last two bytes of 12-byte Intel format */
+ *(q+1) = 0;
#endif
+ }
#endif
/* combine sign and exponent */
i = *p++;
-#if defined(MIEEE) || defined(IBM)
+#ifdef IBM
if (i)
*q++ = *p++ | 0x8000;
else
*q++ = *p++;
*q++ = 0;
-#else
+#endif
+#ifdef DEC
if (i)
*q-- = *p++ | 0x8000;
else
*q-- = *p++;
#endif
+#ifdef IEEE
+ if (REAL_WORDS_BIG_ENDIAN)
+ {
+#ifdef ARM_EXTENDED_IEEE_FORMAT
+ /* The exponent is in the lowest 15 bits of the first word. */
+ *q++ = i ? 0x8000 : 0;
+ *q++ = *p++;
+#else
+ if (i)
+ *q++ = *p++ | 0x8000;
+ else
+ *q++ = *p++;
+ *q++ = 0;
+#endif
+ }
+ else
+ {
+ if (i)
+ *q-- = *p++ | 0x8000;
+ else
+ *q-- = *p++;
+ }
+#endif
/* skip over guard word */
++p;
/* move the significand */
-#if defined(MIEEE) || defined(IBM)
+#ifdef IBM
for (i = 0; i < 4; i++)
*q++ = *p++;
-#else
+#endif
+#ifdef DEC
for (i = 0; i < 4; i++)
*q-- = *p++;
#endif
+#ifdef IEEE
+ if (REAL_WORDS_BIG_ENDIAN)
+ {
+ for (i = 0; i < 4; i++)
+ *q++ = *p++;
+ }
+ else
+ {
+#ifdef INFINITY
+ if (eiisinf (a))
+ {
+ /* Intel long double infinity significand. */
+ *q-- = 0x8000;
+ *q-- = 0;
+ *q-- = 0;
+ *q = 0;
+ return;
+ }
+#endif
+ for (i = 0; i < 4; i++)
+ *q-- = *p++;
+ }
+#endif
}
-
-/* e type to IEEE double precision. */
+/* e type to double precision. */
#ifdef DEC
+/* Convert e-type X to DEC-format double E. */
static void
etoe53 (x, e)
etodec (x, e); /* see etodec.c */
}
+/* Convert exploded e-type X, that has already been rounded to
+ 56-bit double precision, to DEC double Y. */
+
static void
toe53 (x, y)
unsigned EMUSHORT *x, *y;
#else
#ifdef IBM
+/* Convert e-type X to IBM 370-format double E. */
static void
etoe53 (x, e)
etoibm (x, e, DFmode);
}
+/* Convert exploded e-type X, that has already been rounded to
+ 56-bit precision, to IBM 370 double Y. */
+
static void
toe53 (x, y)
unsigned EMUSHORT *x, *y;
#else /* it's neither DEC nor IBM */
+/* Convert e-type X to IEEE double E. */
+
static void
etoe53 (x, e)
unsigned EMUSHORT *x, *e;
toe53 (xi, e);
}
+/* Convert exploded e-type X, that has already been rounded to
+ 53-bit precision, to IEEE double Y. */
static void
toe53 (x, y)
}
#endif
p = &x[0];
-#ifdef IBMPC
- y += 3;
+#ifdef IEEE
+ if (! REAL_WORDS_BIG_ENDIAN)
+ y += 3;
#endif
*y = 0; /* output high order */
if (*p++)
i = *p++;
if (i >= (unsigned int) 2047)
- { /* Saturate at largest number less than infinity. */
+ {
+ /* Saturate at largest number less than infinity. */
#ifdef INFINITY
*y |= 0x7ff0;
-#ifdef IBMPC
- *(--y) = 0;
- *(--y) = 0;
- *(--y) = 0;
-#endif
-#ifdef MIEEE
- ++y;
- *y++ = 0;
- *y++ = 0;
- *y++ = 0;
-#endif
+ if (! REAL_WORDS_BIG_ENDIAN)
+ {
+ *(--y) = 0;
+ *(--y) = 0;
+ *(--y) = 0;
+ }
+ else
+ {
+ ++y;
+ *y++ = 0;
+ *y++ = 0;
+ *y++ = 0;
+ }
#else
*y |= (unsigned EMUSHORT) 0x7fef;
-#ifdef IBMPC
- *(--y) = 0xffff;
- *(--y) = 0xffff;
- *(--y) = 0xffff;
-#endif
-#ifdef MIEEE
- ++y;
- *y++ = 0xffff;
- *y++ = 0xffff;
- *y++ = 0xffff;
-#endif
+ if (! REAL_WORDS_BIG_ENDIAN)
+ {
+ *(--y) = 0xffff;
+ *(--y) = 0xffff;
+ *(--y) = 0xffff;
+ }
+ else
+ {
+ ++y;
+ *y++ = 0xffff;
+ *y++ = 0xffff;
+ *y++ = 0xffff;
+ }
#endif
return;
}
}
i |= *p++ & (unsigned EMUSHORT) 0x0f; /* *p = xi[M] */
*y |= (unsigned EMUSHORT) i; /* high order output already has sign bit set */
-#ifdef IBMPC
- *(--y) = *p++;
- *(--y) = *p++;
- *(--y) = *p;
-#endif
-#ifdef MIEEE
- ++y;
- *y++ = *p++;
- *y++ = *p++;
- *y++ = *p++;
-#endif
+ if (! REAL_WORDS_BIG_ENDIAN)
+ {
+ *(--y) = *p++;
+ *(--y) = *p++;
+ *(--y) = *p;
+ }
+ else
+ {
+ ++y;
+ *y++ = *p++;
+ *y++ = *p++;
+ *y++ = *p++;
+ }
}
#endif /* not IBM */
-/* e type to IEEE single precision. */
+/* e type to single precision. */
#ifdef IBM
+/* Convert e-type X to IBM 370 float E. */
static void
etoe24 (x, e)
etoibm (x, e, SFmode);
}
+/* Convert exploded e-type X, that has already been rounded to
+ float precision, to IBM 370 float Y. */
+
static void
toe24 (x, y)
unsigned EMUSHORT *x, *y;
}
#else
+/* Convert e-type X to IEEE float E. DEC float is the same as IEEE float. */
static void
etoe24 (x, e)
toe24 (xi, e);
}
+/* Convert exploded e-type X, that has already been rounded to
+ float precision, to IEEE float Y. */
+
static void
toe24 (x, y)
unsigned EMUSHORT *x, *y;
}
#endif
p = &x[0];
-#ifdef IBMPC
- y += 1;
+#ifdef IEEE
+ if (! REAL_WORDS_BIG_ENDIAN)
+ y += 1;
#endif
#ifdef DEC
y += 1;
*y = 0x8000; /* output sign bit */
i = *p++;
-/* Handle overflow cases. */
+/* Handle overflow cases. */
if (i >= 255)
{
#ifdef INFINITY
*y |= (unsigned EMUSHORT) 0x7f80;
-#ifdef IBMPC
- *(--y) = 0;
-#endif
#ifdef DEC
*(--y) = 0;
#endif
-#ifdef MIEEE
- ++y;
- *y = 0;
+#ifdef IEEE
+ if (! REAL_WORDS_BIG_ENDIAN)
+ *(--y) = 0;
+ else
+ {
+ ++y;
+ *y = 0;
+ }
#endif
#else /* no INFINITY */
*y |= (unsigned EMUSHORT) 0x7f7f;
-#ifdef IBMPC
- *(--y) = 0xffff;
-#endif
#ifdef DEC
*(--y) = 0xffff;
#endif
-#ifdef MIEEE
- ++y;
- *y = 0xffff;
+#ifdef IEEE
+ if (! REAL_WORDS_BIG_ENDIAN)
+ *(--y) = 0xffff;
+ else
+ {
+ ++y;
+ *y = 0xffff;
+ }
#endif
#ifdef ERANGE
errno = ERANGE;
eshift (x, 8);
}
i |= *p++ & (unsigned EMUSHORT) 0x7f; /* *p = xi[M] */
- *y |= i; /* high order output already has sign bit set */
-#ifdef IBMPC
- *(--y) = *p;
-#endif
+ /* High order output already has sign bit set. */
+ *y |= i;
#ifdef DEC
*(--y) = *p;
#endif
-#ifdef MIEEE
- ++y;
- *y = *p;
+#ifdef IEEE
+ if (! REAL_WORDS_BIG_ENDIAN)
+ *(--y) = *p;
+ else
+ {
+ ++y;
+ *y = *p;
+ }
#endif
}
#endif /* not IBM */
return (0); /* equality */
-
-
diff:
if (*(--p) > *(--q))
return (-msign); /* p is littler */
}
-
-
-
-/* Find nearest integer to x = floor (x + 0.5). */
+/* Find e-type nearest integer to X, as floor (X + 0.5). */
static void
eround (x, y)
efloor (y, y);
}
-
-
-
-/* Convert HOST_WIDE_INT to e type. */
+/* Convert HOST_WIDE_INT LP to e type Y. */
static void
ltoe (lp, y)
emovo (yi, y); /* output the answer */
}
-/* Convert unsigned HOST_WIDE_INT to e type. */
+/* Convert unsigned HOST_WIDE_INT LP to e type Y. */
static void
ultoe (lp, y)
}
-/* Find signed HOST_WIDE_INT integer and floating point fractional
- parts of e-type (packed internal format) floating point input X.
+/* Find signed HOST_WIDE_INT integer I and floating point fractional
+ part FRAC of e-type (packed internal format) floating point input X.
The integer output I has the sign of the input, except that
positive overflow is permitted if FIXUNS_TRUNC_LIKE_FIX_TRUNC.
The output e-type fraction FRAC is the positive fractional
}
-/* Find unsigned HOST_WIDE_INT integer and floating point fractional parts.
- A negative e type input yields integer output = 0
- but correct fraction. */
+/* Find unsigned HOST_WIDE_INT integer I and floating point fractional part
+ FRAC of e-type X. A negative input yields integer output = 0 but
+ correct fraction. */
static void
euifrac (x, i, frac)
*i = (HOST_WIDE_INT) xi[M] & 0xffff;
}
- if (xi[0]) /* A negative value yields unsigned integer 0. */
+ if (xi[0]) /* A negative value yields unsigned integer 0. */
*i = 0L;
xi[0] = 0;
emovo (xi, frac);
}
-
-
-/* Shift significand area up or down by the number of bits given by SC. */
+/* Shift the significand of exploded e-type X up or down by SC bits. */
static int
eshift (x, sc)
return ((int) lost);
}
-
-
-/* Shift normalize the significand area pointed to by argument.
- Shift count (up = positive) is returned. */
+/* Shift normalize the significand area of exploded e-type X.
+ Return the shift count (up = positive). */
static int
enormlz (x)
return (sc);
}
-
-
-
-/* Convert e type number to decimal format ASCII string.
- The constants are for 64 bit precision. */
+/* Powers of ten used in decimal <-> binary conversions. */
#define NTEN 12
#define MAXP 4096
};
#endif
+/* Convert float value X to ASCII string STRING with NDIG digits after
+ the decimal point. */
+
static void
e24toasc (x, string, ndigs)
unsigned EMUSHORT x[];
etoasc (w, string, ndigs);
}
+/* Convert double value X to ASCII string STRING with NDIG digits after
+ the decimal point. */
static void
e53toasc (x, string, ndigs)
etoasc (w, string, ndigs);
}
+/* Convert double extended value X to ASCII string STRING with NDIG digits
+ after the decimal point. */
static void
e64toasc (x, string, ndigs)
etoasc (w, string, ndigs);
}
+/* Convert 128-bit long double value X to ASCII string STRING with NDIG digits
+ after the decimal point. */
+
static void
e113toasc (x, string, ndigs)
unsigned EMUSHORT x[];
etoasc (w, string, ndigs);
}
+/* Convert e-type X to ASCII string STRING with NDIGS digits after
+ the decimal point. */
static char wstring[80]; /* working storage for ASCII output */
if (y[k] != 0)
goto tnzro; /* denormalized number */
}
- goto isone; /* legal all zeros */
+ goto isone; /* valid all zeros */
}
tnzro:
- /* Test for infinity. */
+ /* Test for infinity. */
if (y[NE - 1] == 0x7fff)
{
if (sign)
if (i < 0)
{ /* Number is greater than 1 */
- /* Convert significand to an integer and strip trailing decimal zeros. */
+ /* Convert significand to an integer and strip trailing decimal zeros. */
emov (y, u);
u[NE - 1] = EXONE + NBITS - 1;
emov (eone, t);
m = MAXP;
p = &etens[0][0];
- /* An unordered compare result shouldn't happen here. */
+ /* An unordered compare result shouldn't happen here. */
while (ecmp (ten, u) <= 0)
{
if (ecmp (p, u) <= 0)
}
else
{ /* Number is less than 1.0 */
- /* Pad significand with trailing decimal zeros. */
+ /* Pad significand with trailing decimal zeros. */
if (y[NE - 1] == 0)
{
while ((y[NE - 2] & 0x8000) == 0)
ediv (t, eone, t);
}
isone:
- /* Find the first (leading) digit. */
+ /* Find the first (leading) digit. */
emovi (t, w);
emovz (w, t);
emovi (y, w);
*s++ = '-';
else
*s++ = ' ';
- /* Examine number of digits requested by caller. */
+ /* Examine number of digits requested by caller. */
if (ndigs < 0)
ndigs = 0;
if (ndigs > NDEC)
*s++ = (char)digit + '0';
*s++ = '.';
}
- /* Generate digits after the decimal point. */
+ /* Generate digits after the decimal point. */
for (k = 0; k <= ndigs; k++)
{
/* multiply current number by 10, without normalizing */
/* round off the ASCII string */
if (digit > 4)
{
- /* Test for critical rounding case in ASCII output. */
+ /* Test for critical rounding case in ASCII output. */
if (digit == 5)
{
emovo (y, t);
}
-/* Convert ASCII string to quadruple precision floating point
+/* Convert ASCII string to floating point.
- Numeric input is free field decimal number with max of 15 digits with or
- without decimal point entered as ASCII from teletype. Entering E after
- the number followed by a second number causes the second number to be
- interpreted as a power of 10 to be multiplied by the first number
- (i.e., "scientific" notation). */
+ Numeric input is a free format decimal number of any length, with
+ or without decimal point. Entering E after the number followed by an
+ integer number causes the second number to be interpreted as a power of
+ 10 to be multiplied by the first number (i.e., "scientific" notation). */
-/* ASCII to single */
+/* Convert ASCII string S to single precision float value Y. */
static void
asctoe24 (s, y)
}
-/* ASCII to double */
+/* Convert ASCII string S to double precision value Y. */
static void
asctoe53 (s, y)
}
-/* ASCII to long double */
+/* Convert ASCII string S to double extended value Y. */
static void
asctoe64 (s, y)
asctoeg (s, y, 64);
}
-/* ASCII to 128-bit long double */
+/* Convert ASCII string S to 128-bit long double Y. */
static void
asctoe113 (s, y)
asctoeg (s, y, 113);
}
-/* ASCII to super double */
+/* Convert ASCII string S to e type Y. */
static void
asctoe (s, y)
asctoeg (s, y, NBITS);
}
-
-/* ASCII to e type, with specified rounding precision = oprec. */
+/* Convert ASCII string SS to e type Y, with a specified rounding precision
+ of OPREC bits. */
static void
asctoeg (ss, y, oprec)
unsigned EMUSHORT nsign, *p;
char *sp, *s, *lstr;
- /* Copy the input string. */
+ /* Copy the input string. */
lstr = (char *) alloca (strlen (ss) + 1);
s = ss;
while (*s == ' ') /* skip leading spaces */
/* Ignore leading zeros */
if ((prec == 0) && (decflg == 0) && (k == 0))
goto donchr;
- /* Identify and strip trailing zeros after the decimal point. */
+ /* Identify and strip trailing zeros after the decimal point. */
if ((trail == 0) && (decflg != 0))
{
sp = s;
/* Exponent interpretation */
expnt:
+ /* 0.0eXXX is zero, regardless of XXX. Check for the 0.0. */
+ for (k = 0; k < NI; k++)
+ {
+ if (yy[k] != 0)
+ goto read_expnt;
+ }
+ goto aexit;
+read_expnt:
esign = 1;
exp = 0;
++s;
daldone:
nexp = exp - nexp;
- /* Pad trailing zeros to minimize power of 10, per IEEE spec. */
+ /* Pad trailing zeros to minimize power of 10, per IEEE spec. */
while ((nexp > 0) && (yy[2] == 0))
{
emovz (yy, xt);
esign = -1;
if (nexp > 4096)
{
- /* Punt. Can't handle this without 2 divides. */
+ /* Punt. Can't handle this without 2 divides. */
emovi (etens[0], tt);
lexp -= tt[E];
k = edivm (tt, yy);
-/* y = largest integer not greater than x (truncated toward minus infinity) */
+/* Return Y = largest integer not greater than X (truncated toward minus
+ infinity). */
static unsigned EMUSHORT bmask[] =
{
}
-/* Returns s and exp such that s * 2**exp = x and .5 <= s < 1.
- For example, 1.1 = 0.55 * 2**1
- Handles denormalized numbers properly using long integer exp. */
+/* Return S and EXP such that S * 2^EXP = X and .5 <= S < 1.
+ For example, 1.1 = 0.55 * 2^1. */
static void
efrexp (x, exp, s)
EMULONG li;
emovi (x, xi);
+ /* Handle denormalized numbers properly using long integer exponent. */
li = (EMULONG) ((EMUSHORT) xi[1]);
if (li == 0)
*exp = (int) (li - 0x3ffe);
}
-
-
-/* Return y = x * 2**pwr2. */
+/* Return e type Y = X * 2^PWR2. */
static void
eldexp (x, pwr2, y)
}
-/* c = remainder after dividing b by a
- Least significant integer quotient bits left in equot[]. */
+/* C = remainder after dividing B by A, all e type values.
+ Least significant integer quotient bits left in EQUOT. */
static void
eremain (a, b, c)
emovo (num, c);
}
+/* Return quotient of exploded e-types NUM / DEN in EQUOT,
+ remainder in NUM. */
+
static void
eiremain (den, num)
unsigned EMUSHORT den[], num[];
j = 1;
}
else
- {
j = 0;
- }
eshup1 (equot);
equot[NI - 1] |= j;
eshup1 (num);
emdnorm (num, 0, 0, ln, 0);
}
-/* This routine may be called to report one of the following
- error conditions (in the include file mconf.h).
+/* Report an error condition CODE encountered in function NAME.
+ CODE is one of the following:
Mnemonic Value Significance
EDOM 33 Unix domain error code
ERANGE 34 Unix range error code
- The default version of the file prints the function name,
- passed to it by the pointer fctnam, followed by the
- error condition. The display is directed to the standard
- output device. The routine then returns to the calling
- program. Users may wish to modify the program to abort by
- calling exit under severe error conditions such as domain
- errors.
-
- Since all error conditions pass control to this function,
- the display may be easily changed, eliminated, or directed
- to an error logging device. */
-
-/* Note: the order of appearance of the following messages is bound to the
+ The order of appearance of the following messages is bound to the
error codes defined above. */
#define NMSGS 8
{
char errstr[80];
- /* Display string passed by calling program, which is supposed to be the
+ /* The string passed by the calling program is supposed to be the
name of the function in which the error occurred.
-
- Display error message defined by the code argument. */
+ The code argument selects which error message string will be printed. */
if ((code <= 0) || (code >= NMSGS))
code = 0;
}
#ifdef DEC
-/* Convert DEC double precision to e type. */
+/* Convert DEC double precision D to e type E. */
static void
dectoe (d, e)
emovo (y, e);
}
-
-
-/*
-; convert e type to DEC double precision
-; double d;
-; EMUSHORT e[NE];
-; etodec (e, &d);
-*/
+/* Convert e type X to DEC double precision D. */
static void
etodec (x, d)
int rndsav;
emovi (x, xi);
- exp = (EMULONG) xi[E] - (EXONE - 0201); /* adjust exponent for offsets */
-/* round off to nearest or even */
+ /* Adjust exponent for offsets. */
+ exp = (EMULONG) xi[E] - (EXONE - 0201);
+ /* Round off to nearest or even. */
rndsav = rndprc;
rndprc = 56;
emdnorm (xi, 0, 0, exp, 64);
todec (xi, d);
}
+/* Convert exploded e-type X, that has already been rounded to
+ 56-bit precision, to DEC format double Y. */
+
static void
todec (x, y)
unsigned EMUSHORT *x, *y;
/* If special NaN bit patterns are required, define them in tm.h
as arrays of unsigned 16-bit shorts. Otherwise, use the default
- patterns here. */
+ patterns here. */
#ifdef TFMODE_NAN
TFMODE_NAN;
#else
-#ifdef MIEEE
-unsigned EMUSHORT TFnan[8] =
+#ifdef IEEE
+unsigned EMUSHORT TFbignan[8] =
{0x7fff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff};
-#endif
-#ifdef IBMPC
-unsigned EMUSHORT TFnan[8] = {0, 0, 0, 0, 0, 0, 0x8000, 0xffff};
+unsigned EMUSHORT TFlittlenan[8] = {0, 0, 0, 0, 0, 0, 0x8000, 0xffff};
#endif
#endif
#ifdef XFMODE_NAN
XFMODE_NAN;
#else
-#ifdef MIEEE
-unsigned EMUSHORT XFnan[6] = {0x7fff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff};
-#endif
-#ifdef IBMPC
-unsigned EMUSHORT XFnan[6] = {0, 0, 0, 0xc000, 0xffff, 0};
+#ifdef IEEE
+unsigned EMUSHORT XFbignan[6] =
+ {0x7fff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff};
+unsigned EMUSHORT XFlittlenan[6] = {0, 0, 0, 0xc000, 0xffff, 0};
#endif
#endif
#ifdef DFMODE_NAN
DFMODE_NAN;
#else
-#ifdef MIEEE
-unsigned EMUSHORT DFnan[4] = {0x7fff, 0xffff, 0xffff, 0xffff};
-#endif
-#ifdef IBMPC
-unsigned EMUSHORT DFnan[4] = {0, 0, 0, 0xfff8};
+#ifdef IEEE
+unsigned EMUSHORT DFbignan[4] = {0x7fff, 0xffff, 0xffff, 0xffff};
+unsigned EMUSHORT DFlittlenan[4] = {0, 0, 0, 0xfff8};
#endif
#endif
#ifdef SFMODE_NAN
SFMODE_NAN;
#else
-#ifdef MIEEE
-unsigned EMUSHORT SFnan[2] = {0x7fff, 0xffff};
-#endif
-#ifdef IBMPC
-unsigned EMUSHORT SFnan[2] = {0, 0xffc0};
+#ifdef IEEE
+unsigned EMUSHORT SFbignan[2] = {0x7fff, 0xffff};
+unsigned EMUSHORT SFlittlenan[2] = {0, 0xffc0};
#endif
#endif
switch (mode)
{
/* Possibly the `reserved operand' patterns on a VAX can be
- used like NaN's, but probably not in the same way as IEEE. */
+ used like NaN's, but probably not in the same way as IEEE. */
#if !defined(DEC) && !defined(IBM)
case TFmode:
n = 8;
- p = TFnan;
+ if (REAL_WORDS_BIG_ENDIAN)
+ p = TFbignan;
+ else
+ p = TFlittlenan;
break;
case XFmode:
n = 6;
- p = XFnan;
+ if (REAL_WORDS_BIG_ENDIAN)
+ p = XFbignan;
+ else
+ p = XFlittlenan;
break;
case DFmode:
n = 4;
- p = DFnan;
+ if (REAL_WORDS_BIG_ENDIAN)
+ p = DFbignan;
+ else
+ p = DFlittlenan;
break;
case HFmode:
case SFmode:
n = 2;
- p = SFnan;
+ if (REAL_WORDS_BIG_ENDIAN)
+ p = SFbignan;
+ else
+ p = SFlittlenan;
break;
#endif
default:
abort ();
}
-#ifdef MIEEE
- *nan++ = (sign << 15) | *p++;
-#endif
+ if (REAL_WORDS_BIG_ENDIAN)
+ *nan++ = (sign << 15) | *p++;
while (--n != 0)
*nan++ = *p++;
-#ifndef MIEEE
- *nan = (sign << 15) | *p;
-#endif
+ if (! REAL_WORDS_BIG_ENDIAN)
+ *nan = (sign << 15) | *p;
}
-/* Convert an SFmode target `float' value to a REAL_VALUE_TYPE.
- This is the inverse of the function `etarsingle' invoked by
+/* This is the inverse of the function `etarsingle' invoked by
REAL_VALUE_TO_TARGET_SINGLE. */
REAL_VALUE_TYPE
-ereal_from_float (f)
- HOST_WIDE_INT f;
+ereal_unto_float (f)
+ long f;
{
REAL_VALUE_TYPE r;
unsigned EMUSHORT s[2];
/* Convert 32 bit integer to array of 16 bit pieces in target machine order.
This is the inverse operation to what the function `endian' does. */
-#if FLOAT_WORDS_BIG_ENDIAN
- s[0] = (unsigned EMUSHORT) (f >> 16);
- s[1] = (unsigned EMUSHORT) f;
-#else
- s[0] = (unsigned EMUSHORT) f;
- s[1] = (unsigned EMUSHORT) (f >> 16);
-#endif
+ if (REAL_WORDS_BIG_ENDIAN)
+ {
+ s[0] = (unsigned EMUSHORT) (f >> 16);
+ s[1] = (unsigned EMUSHORT) f;
+ }
+ else
+ {
+ s[0] = (unsigned EMUSHORT) f;
+ s[1] = (unsigned EMUSHORT) (f >> 16);
+ }
/* Convert and promote the target float to E-type. */
e24toe (s, e);
/* Output E-type to REAL_VALUE_TYPE. */
}
+/* This is the inverse of the function `etardouble' invoked by
+ REAL_VALUE_TO_TARGET_DOUBLE. */
+
+REAL_VALUE_TYPE
+ereal_unto_double (d)
+ long d[];
+{
+ REAL_VALUE_TYPE r;
+ unsigned EMUSHORT s[4];
+ unsigned EMUSHORT e[NE];
+
+ /* Convert array of HOST_WIDE_INT to equivalent array of 16-bit pieces. */
+ if (REAL_WORDS_BIG_ENDIAN)
+ {
+ s[0] = (unsigned EMUSHORT) (d[0] >> 16);
+ s[1] = (unsigned EMUSHORT) d[0];
+ s[2] = (unsigned EMUSHORT) (d[1] >> 16);
+ s[3] = (unsigned EMUSHORT) d[1];
+ }
+ else
+ {
+ /* Target float words are little-endian. */
+ s[0] = (unsigned EMUSHORT) d[0];
+ s[1] = (unsigned EMUSHORT) (d[0] >> 16);
+ s[2] = (unsigned EMUSHORT) d[1];
+ s[3] = (unsigned EMUSHORT) (d[1] >> 16);
+ }
+ /* Convert target double to E-type. */
+ e53toe (s, e);
+ /* Output E-type to REAL_VALUE_TYPE. */
+ PUT_REAL (e, &r);
+ return r;
+}
+
+
+/* Convert an SFmode target `float' value to a REAL_VALUE_TYPE.
+ This is somewhat like ereal_unto_float, but the input types
+ for these are different. */
+
+REAL_VALUE_TYPE
+ereal_from_float (f)
+ HOST_WIDE_INT f;
+{
+ REAL_VALUE_TYPE r;
+ unsigned EMUSHORT s[2];
+ unsigned EMUSHORT e[NE];
+
+ /* Convert 32 bit integer to array of 16 bit pieces in target machine order.
+ This is the inverse operation to what the function `endian' does. */
+ if (REAL_WORDS_BIG_ENDIAN)
+ {
+ s[0] = (unsigned EMUSHORT) (f >> 16);
+ s[1] = (unsigned EMUSHORT) f;
+ }
+ else
+ {
+ s[0] = (unsigned EMUSHORT) f;
+ s[1] = (unsigned EMUSHORT) (f >> 16);
+ }
+ /* Convert and promote the target float to E-type. */
+ e24toe (s, e);
+ /* Output E-type to REAL_VALUE_TYPE. */
+ PUT_REAL (e, &r);
+ return r;
+}
+
+
/* Convert a DFmode target `double' value to a REAL_VALUE_TYPE.
- This is the inverse of the function `etardouble' invoked by
- REAL_VALUE_TO_TARGET_DOUBLE.
+ This is somewhat like ereal_unto_double, but the input types
+ for these are different.
The DFmode is stored as an array of HOST_WIDE_INT in the target's
data format, with no holes in the bit packing. The first element
unsigned EMUSHORT e[NE];
/* Convert array of HOST_WIDE_INT to equivalent array of 16-bit pieces. */
-#if FLOAT_WORDS_BIG_ENDIAN
- s[0] = (unsigned EMUSHORT) (d[0] >> 16);
- s[1] = (unsigned EMUSHORT) d[0];
+ if (REAL_WORDS_BIG_ENDIAN)
+ {
+ s[0] = (unsigned EMUSHORT) (d[0] >> 16);
+ s[1] = (unsigned EMUSHORT) d[0];
#if HOST_BITS_PER_WIDE_INT == 32
- s[2] = (unsigned EMUSHORT) (d[1] >> 16);
- s[3] = (unsigned EMUSHORT) d[1];
+ s[2] = (unsigned EMUSHORT) (d[1] >> 16);
+ s[3] = (unsigned EMUSHORT) d[1];
#else
- /* In this case the entire target double is contained in the
- first array element. The second element of the input is ignored. */
- s[2] = (unsigned EMUSHORT) (d[0] >> 48);
- s[3] = (unsigned EMUSHORT) (d[0] >> 32);
+ /* In this case the entire target double is contained in the
+ first array element. The second element of the input is
+ ignored. */
+ s[2] = (unsigned EMUSHORT) (d[0] >> 48);
+ s[3] = (unsigned EMUSHORT) (d[0] >> 32);
#endif
-#else
-/* Target float words are little-endian. */
- s[0] = (unsigned EMUSHORT) d[0];
- s[1] = (unsigned EMUSHORT) (d[0] >> 16);
+ }
+ else
+ {
+ /* Target float words are little-endian. */
+ s[0] = (unsigned EMUSHORT) d[0];
+ s[1] = (unsigned EMUSHORT) (d[0] >> 16);
#if HOST_BITS_PER_WIDE_INT == 32
- s[2] = (unsigned EMUSHORT) d[1];
- s[3] = (unsigned EMUSHORT) (d[1] >> 16);
+ s[2] = (unsigned EMUSHORT) d[1];
+ s[3] = (unsigned EMUSHORT) (d[1] >> 16);
#else
- s[2] = (unsigned EMUSHORT) (d[0] >> 32);
- s[3] = (unsigned EMUSHORT) (d[0] >> 48);
-#endif
+ s[2] = (unsigned EMUSHORT) (d[0] >> 32);
+ s[3] = (unsigned EMUSHORT) (d[0] >> 48);
#endif
- /* Convert target double to E-type. */
+ }
+ /* Convert target double to E-type. */
e53toe (s, e);
- /* Output E-type to REAL_VALUE_TYPE. */
+ /* Output E-type to REAL_VALUE_TYPE. */
PUT_REAL (e, &r);
return r;
}
/* Convert target computer unsigned 64-bit integer to e-type.
The endian-ness of DImode follows the convention for integers,
- so we use WORDS_BIG_ENDIAN here, not FLOAT_WORDS_BIG_ENDIAN. */
+ so we use WORDS_BIG_ENDIAN here, not REAL_WORDS_BIG_ENDIAN. */
static void
uditoe (di, e)
- unsigned EMUSHORT *di; /* Address of the 64-bit int. */
+ unsigned EMUSHORT *di; /* Address of the 64-bit int. */
unsigned EMUSHORT *e;
{
unsigned EMUSHORT yi[NI];
int k;
ecleaz (yi);
-#if WORDS_BIG_ENDIAN
- for (k = M; k < M + 4; k++)
- yi[k] = *di++;
-#else
- for (k = M + 3; k >= M; k--)
- yi[k] = *di++;
-#endif
+ if (WORDS_BIG_ENDIAN)
+ {
+ for (k = M; k < M + 4; k++)
+ yi[k] = *di++;
+ }
+ else
+ {
+ for (k = M + 3; k >= M; k--)
+ yi[k] = *di++;
+ }
yi[E] = EXONE + 47; /* exponent if normalize shift count were 0 */
if ((k = enormlz (yi)) > NBITS)/* normalize the significand */
ecleaz (yi); /* it was zero */
emovo (yi, e);
}
-/* Convert target computer signed 64-bit integer to e-type. */
+/* Convert target computer signed 64-bit integer to e-type. */
static void
ditoe (di, e)
- unsigned EMUSHORT *di; /* Address of the 64-bit int. */
+ unsigned EMUSHORT *di; /* Address of the 64-bit int. */
unsigned EMUSHORT *e;
{
unsigned EMULONG acc;
int k, sign;
ecleaz (yi);
-#if WORDS_BIG_ENDIAN
- for (k = M; k < M + 4; k++)
- yi[k] = *di++;
-#else
- for (k = M + 3; k >= M; k--)
- yi[k] = *di++;
-#endif
+ if (WORDS_BIG_ENDIAN)
+ {
+ for (k = M; k < M + 4; k++)
+ yi[k] = *di++;
+ }
+ else
+ {
+ for (k = M + 3; k >= M; k--)
+ yi[k] = *di++;
+ }
/* Take absolute value */
sign = 0;
if (yi[M] & 0x8000)
}
-/* Convert e-type to unsigned 64-bit int. */
+/* Convert e-type to unsigned 64-bit int. */
static void
etoudi (x, i)
if (j == 0)
j = 16;
eshift (xi, j);
-#if WORDS_BIG_ENDIAN
- *i++ = xi[M];
-#else
- i += 3;
- *i-- = xi[M];
-#endif
+ if (WORDS_BIG_ENDIAN)
+ *i++ = xi[M];
+ else
+ {
+ i += 3;
+ *i-- = xi[M];
+ }
k -= j;
do
{
eshup6 (xi);
-#if WORDS_BIG_ENDIAN
- *i++ = xi[M];
-#else
- *i-- = xi[M];
-#endif
+ if (WORDS_BIG_ENDIAN)
+ *i++ = xi[M];
+ else
+ *i-- = xi[M];
}
while ((k -= 16) > 0);
}
noshift:
-#if WORDS_BIG_ENDIAN
- i += 3;
- *i-- = xi[M];
- *i-- = 0;
- *i-- = 0;
- *i = 0;
-#else
- *i++ = xi[M];
- *i++ = 0;
- *i++ = 0;
- *i = 0;
-#endif
+ if (WORDS_BIG_ENDIAN)
+ {
+ i += 3;
+ *i-- = xi[M];
+ *i-- = 0;
+ *i-- = 0;
+ *i = 0;
+ }
+ else
+ {
+ *i++ = xi[M];
+ *i++ = 0;
+ *i++ = 0;
+ *i = 0;
+ }
}
}
-/* Convert e-type to signed 64-bit int. */
+/* Convert e-type to signed 64-bit int. */
static void
etodi (x, i)
if (j == 0)
j = 16;
eshift (xi, j);
-#if WORDS_BIG_ENDIAN
- *i++ = xi[M];
-#else
- i += 3;
- *i-- = xi[M];
-#endif
+ if (WORDS_BIG_ENDIAN)
+ *i++ = xi[M];
+ else
+ {
+ i += 3;
+ *i-- = xi[M];
+ }
k -= j;
do
{
eshup6 (xi);
-#if WORDS_BIG_ENDIAN
- *i++ = xi[M];
-#else
- *i-- = xi[M];
-#endif
+ if (WORDS_BIG_ENDIAN)
+ *i++ = xi[M];
+ else
+ *i-- = xi[M];
}
while ((k -= 16) > 0);
}
/* shift not more than 16 bits */
eshift (xi, k);
-#if WORDS_BIG_ENDIAN
- i += 3;
- *i = xi[M];
- *i-- = 0;
- *i-- = 0;
- *i = 0;
-#else
- *i++ = xi[M];
- *i++ = 0;
- *i++ = 0;
- *i = 0;
-#endif
+ if (WORDS_BIG_ENDIAN)
+ {
+ i += 3;
+ *i = xi[M];
+ *i-- = 0;
+ *i-- = 0;
+ *i = 0;
+ }
+ else
+ {
+ *i++ = xi[M];
+ *i++ = 0;
+ *i++ = 0;
+ *i = 0;
+ }
}
/* Negate if negative */
if (xi[0])
{
carry = 0;
-#if WORDS_BIG_ENDIAN
- isave += 3;
-#endif
+ if (WORDS_BIG_ENDIAN)
+ isave += 3;
for (k = 0; k < 4; k++)
{
acc = (unsigned EMULONG) (~(*isave) & 0xffff) + carry;
-#if WORDS_BIG_ENDIAN
- *isave-- = acc;
-#else
- *isave++ = acc;
-#endif
+ if (WORDS_BIG_ENDIAN)
+ *isave-- = acc;
+ else
+ *isave++ = acc;
carry = 0;
if (acc & 0x10000)
carry = 1;
}
-/* Longhand square root routine. */
+/* Longhand square root routine. */
static int esqinited = 0;
return;
}
#endif
- /* Bring in the arg and renormalize if it is denormal. */
+ /* Bring in the arg and renormalize if it is denormal. */
emovi (x, xx);
m = (EMULONG) xx[1]; /* local long word exponent */
if (m == 0)
/* bring in next word of arg */
if (j < NE)
num[NI - 1] = xx[j + 3];
- /* Do additional bit on last outer loop, for roundoff. */
+ /* Do additional bit on last outer loop, for roundoff. */
if (nlups <= 8)
n = nlups + 1;
for (i = 0; i < n; i++)
j += 1;
}
- /* Adjust for extra, roundoff loop done. */
+ /* Adjust for extra, roundoff loop done. */
exp += (NBITS - 1) - rndprc;
- /* Sticky bit = 1 if the remainder is nonzero. */
+ /* Sticky bit = 1 if the remainder is nonzero. */
k = 0;
for (i = 3; i < NI; i++)
k |= (int) num[i];
- /* Renormalize and round off. */
+ /* Renormalize and round off. */
emdnorm (sq, k, 0, exp, 64);
emovo (sq, y);
}
enum machine_mode mode;
{
-switch (mode)
+/* Don't test the modes, but their sizes, lest this
+ code won't work for BITS_PER_UNIT != 8 . */
+
+switch (GET_MODE_BITSIZE (mode))
{
- case SFmode:
+ case 32:
return 24;
- case DFmode:
+ case 64:
#if TARGET_FLOAT_FORMAT == IEEE_FLOAT_FORMAT
return 53;
#else
#endif
#endif
- case XFmode:
+ case 96:
return 64;
- case TFmode:
+ case 128:
return 113;
default:
OSDN Git Service
