1993-07-09 17:11:49 +04:00
|
|
|
|
/* real.c - implementation of REAL_ARITHMETIC, REAL_VALUE_ATOF,
|
|
|
|
|
and support for XFmode IEEE extended real floating point arithmetic.
|
|
|
|
|
Contributed by Stephen L. Moshier (moshier@world.std.com).
|
|
|
|
|
|
|
|
|
|
Copyright (C) 1993 Free Software Foundation, Inc.
|
|
|
|
|
|
|
|
|
|
This file is part of GNU CC.
|
|
|
|
|
|
|
|
|
|
GNU CC is free software; you can redistribute it and/or modify
|
|
|
|
|
it under the terms of the GNU General Public License as published by
|
|
|
|
|
the Free Software Foundation; either version 2, or (at your option)
|
|
|
|
|
any later version.
|
|
|
|
|
|
|
|
|
|
GNU CC is distributed in the hope that it will be useful,
|
|
|
|
|
but WITHOUT ANY WARRANTY; without even the implied warranty of
|
|
|
|
|
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
|
|
|
|
|
GNU General Public License for more details.
|
|
|
|
|
|
|
|
|
|
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. */
|
|
|
|
|
|
1993-08-02 21:28:42 +04:00
|
|
|
|
#ifndef lint
|
|
|
|
|
static char rcsid[] = "$Id: real.c,v 1.2 1993/08/02 17:35:25 mycroft Exp $";
|
|
|
|
|
#endif /* not lint */
|
|
|
|
|
|
1993-07-09 17:11:49 +04:00
|
|
|
|
#include <stdio.h>
|
|
|
|
|
#include <errno.h>
|
|
|
|
|
#include "config.h"
|
|
|
|
|
#include "tree.h"
|
|
|
|
|
|
|
|
|
|
#ifndef errno
|
|
|
|
|
extern int errno;
|
|
|
|
|
#endif
|
|
|
|
|
|
|
|
|
|
/* To enable support of XFmode extended real floating point, define
|
|
|
|
|
LONG_DOUBLE_TYPE_SIZE 96 in the tm.h file (m68k.h or i386.h).
|
|
|
|
|
|
|
|
|
|
To support cross compilation between IEEE and VAX floating
|
|
|
|
|
point formats, define REAL_ARITHMETIC in the tm.h file.
|
|
|
|
|
|
|
|
|
|
In either case the machine files (tm.h) must not contain any code
|
|
|
|
|
that tries to use host floating point arithmetic to convert
|
|
|
|
|
REAL_VALUE_TYPEs from `double' to `float', pass them to fprintf,
|
|
|
|
|
etc. In cross-compile situations a REAL_VALUE_TYPE may not
|
|
|
|
|
be intelligible to the host computer's native arithmetic.
|
|
|
|
|
|
|
|
|
|
The emulator defaults to the host's floating point format so that
|
|
|
|
|
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 */
|
|
|
|
|
|
|
|
|
|
/* Type of computer arithmetic.
|
|
|
|
|
* Only one of DEC, 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.
|
|
|
|
|
|
|
|
|
|
`DEC' refers specifically to the Digital Equipment Corp PDP-11
|
|
|
|
|
and VAX floating point data structure. This model currently
|
|
|
|
|
supports no type wider than DFmode.
|
|
|
|
|
|
|
|
|
|
If LONG_DOUBLE_TYPE_SIZE = 64 (the default, unless tm.h defines it)
|
|
|
|
|
then `long double' and `double' are both implemented, but they
|
|
|
|
|
both mean DFmode. In this case, the software floating-point
|
|
|
|
|
support available here is activated by writing
|
|
|
|
|
#define REAL_ARITHMETIC
|
|
|
|
|
in tm.h.
|
|
|
|
|
|
|
|
|
|
The case LONG_DOUBLE_TYPE_SIZE = 128 activates TFmode support
|
|
|
|
|
(Not Yet Implemented) and may deactivate XFmode since
|
|
|
|
|
`long double' is used to refer to both modes. */
|
|
|
|
|
|
|
|
|
|
/* The following converts gcc macros into the ones used by this file. */
|
|
|
|
|
|
|
|
|
|
/* REAL_ARITHMETIC defined means that macros in real.h are
|
|
|
|
|
defined to call emulator functions. */
|
|
|
|
|
#ifdef REAL_ARITHMETIC
|
|
|
|
|
|
|
|
|
|
#if TARGET_FLOAT_FORMAT == VAX_FLOAT_FORMAT
|
|
|
|
|
/* PDP-11, Pro350, VAX: */
|
|
|
|
|
#define DEC 1
|
|
|
|
|
#else /* it's not VAX */
|
|
|
|
|
#if TARGET_FLOAT_FORMAT == IEEE_FLOAT_FORMAT
|
|
|
|
|
#if 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 */
|
|
|
|
|
#else /* it's not IEEE either */
|
|
|
|
|
/* UNKnown arithmetic. We don't support this and can't go on. */
|
|
|
|
|
unknown arithmetic type
|
|
|
|
|
#define UNK 1
|
|
|
|
|
#endif /* not IEEE */
|
|
|
|
|
#endif /* not VAX */
|
|
|
|
|
|
|
|
|
|
#else
|
|
|
|
|
/* REAL_ARITHMETIC not defined means that the *host's* data
|
|
|
|
|
structure will be used. It may differ by endian-ness from the
|
|
|
|
|
target machine's structure and will get its ends swapped
|
|
|
|
|
accordingly (but not here). Probably only the decimal <-> binary
|
|
|
|
|
functions in this file will actually be used in this case. */
|
|
|
|
|
#if HOST_FLOAT_FORMAT == VAX_FLOAT_FORMAT
|
|
|
|
|
#define DEC 1
|
|
|
|
|
#else /* it's not VAX */
|
|
|
|
|
#if HOST_FLOAT_FORMAT == IEEE_FLOAT_FORMAT
|
|
|
|
|
#ifdef HOST_WORDS_BIG_ENDIAN
|
|
|
|
|
#define MIEEE 1
|
|
|
|
|
#else /* not big-endian */
|
|
|
|
|
#define IBMPC 1
|
|
|
|
|
#endif /* big-endian */
|
|
|
|
|
#else /* it's not IEEE either */
|
|
|
|
|
unknown arithmetic type
|
|
|
|
|
#define UNK 1
|
|
|
|
|
#endif /* not IEEE */
|
|
|
|
|
#endif /* not VAX */
|
|
|
|
|
|
|
|
|
|
#endif /* REAL_ARITHMETIC not defined */
|
|
|
|
|
|
|
|
|
|
/* Define INFINITY for support of infinity.
|
|
|
|
|
Define NANS for support of Not-a-Number's (NaN's). */
|
|
|
|
|
#ifndef DEC
|
|
|
|
|
#define INFINITY
|
|
|
|
|
#define NANS
|
|
|
|
|
#endif
|
|
|
|
|
|
|
|
|
|
/* Support of NaNs requires support of infinity. */
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
#ifndef INFINITY
|
|
|
|
|
#define INFINITY
|
|
|
|
|
#endif
|
|
|
|
|
#endif
|
|
|
|
|
|
|
|
|
|
/* ehead.h
|
|
|
|
|
*
|
|
|
|
|
* Include file for extended precision arithmetic programs.
|
|
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
/* Number of 16 bit words in external e type format */
|
|
|
|
|
#define NE 6
|
|
|
|
|
|
|
|
|
|
/* Number of 16 bit words in internal format */
|
|
|
|
|
#define NI (NE+3)
|
|
|
|
|
|
|
|
|
|
/* Array offset to exponent */
|
|
|
|
|
#define E 1
|
|
|
|
|
|
|
|
|
|
/* Array offset to high guard word */
|
|
|
|
|
#define M 2
|
|
|
|
|
|
|
|
|
|
/* Number of bits of precision */
|
|
|
|
|
#define NBITS ((NI-4)*16)
|
|
|
|
|
|
|
|
|
|
/* Maximum number of decimal digits in ASCII conversion
|
|
|
|
|
* = NBITS*log10(2)
|
|
|
|
|
*/
|
|
|
|
|
#define NDEC (NBITS*8/27)
|
|
|
|
|
|
|
|
|
|
/* The exponent of 1.0 */
|
|
|
|
|
#define EXONE (0x3fff)
|
|
|
|
|
|
|
|
|
|
/* Find a host integer type that is at least 16 bits wide,
|
|
|
|
|
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_CHAR
|
|
|
|
|
#define EMULONG_SIZE (2 * HOST_BITS_PER_CHAR)
|
|
|
|
|
#else
|
|
|
|
|
#if HOST_BITS_PER_SHORT >= 16
|
|
|
|
|
#define EMUSHORT short
|
|
|
|
|
#define EMUSHORT_SIZE HOST_BITS_PER_SHORT
|
|
|
|
|
#define EMULONG_SIZE (2 * HOST_BITS_PER_SHORT)
|
|
|
|
|
#else
|
|
|
|
|
#if HOST_BITS_PER_INT >= 16
|
|
|
|
|
#define EMUSHORT int
|
|
|
|
|
#define EMUSHORT_SIZE HOST_BITS_PER_INT
|
|
|
|
|
#define EMULONG_SIZE (2 * HOST_BITS_PER_INT)
|
|
|
|
|
#else
|
|
|
|
|
#if HOST_BITS_PER_LONG >= 16
|
|
|
|
|
#define EMUSHORT long
|
|
|
|
|
#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. */
|
|
|
|
|
#define EMU_NON_COMPILE
|
|
|
|
|
#endif
|
|
|
|
|
#endif
|
|
|
|
|
#endif
|
|
|
|
|
#endif
|
|
|
|
|
|
|
|
|
|
#if HOST_BITS_PER_SHORT >= EMULONG_SIZE
|
|
|
|
|
#define EMULONG short
|
|
|
|
|
#else
|
|
|
|
|
#if HOST_BITS_PER_INT >= EMULONG_SIZE
|
|
|
|
|
#define EMULONG int
|
|
|
|
|
#else
|
|
|
|
|
#if HOST_BITS_PER_LONG >= EMULONG_SIZE
|
|
|
|
|
#define EMULONG long
|
|
|
|
|
#else
|
|
|
|
|
#if HOST_BITS_PER_LONG_LONG >= EMULONG_SIZE
|
|
|
|
|
#define EMULONG long long int
|
|
|
|
|
#else
|
|
|
|
|
/* You will have to modify this program to have a smaller unit size. */
|
|
|
|
|
#define EMU_NON_COMPILE
|
|
|
|
|
#endif
|
|
|
|
|
#endif
|
|
|
|
|
#endif
|
|
|
|
|
#endif
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* 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. */
|
|
|
|
|
#ifndef EMU_NON_COMPILE
|
|
|
|
|
|
|
|
|
|
/* Construct macros to translate between REAL_VALUE_TYPE and e type.
|
|
|
|
|
In GET_REAL and PUT_REAL, r and e are pointers.
|
|
|
|
|
A REAL_VALUE_TYPE is guaranteed to occupy contiguous locations
|
|
|
|
|
in memory, with no holes. */
|
|
|
|
|
|
|
|
|
|
#if LONG_DOUBLE_TYPE_SIZE == 96
|
|
|
|
|
#define GET_REAL(r,e) bcopy (r, e, 2*NE)
|
|
|
|
|
#define PUT_REAL(e,r) bcopy (e, r, 2*NE)
|
|
|
|
|
#else /* no XFmode */
|
|
|
|
|
|
|
|
|
|
#ifdef REAL_ARITHMETIC
|
|
|
|
|
/* Emulator uses target format internally
|
|
|
|
|
but host stores it in host endian-ness. */
|
|
|
|
|
|
|
|
|
|
#if defined (HOST_WORDS_BIG_ENDIAN) == WORDS_BIG_ENDIAN
|
|
|
|
|
#define GET_REAL(r,e) e53toe ((r), (e))
|
|
|
|
|
#define PUT_REAL(e,r) etoe53 ((e), (r))
|
|
|
|
|
|
|
|
|
|
#else /* endian-ness differs */
|
|
|
|
|
/* emulator uses target endian-ness internally */
|
|
|
|
|
#define GET_REAL(r,e) \
|
|
|
|
|
do { 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 { 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 */
|
|
|
|
|
|
|
|
|
|
#else /* not REAL_ARITHMETIC */
|
|
|
|
|
|
|
|
|
|
/* emulator uses host format */
|
|
|
|
|
#define GET_REAL(r,e) e53toe ((r), (e))
|
|
|
|
|
#define PUT_REAL(e,r) etoe53 ((e), (r))
|
|
|
|
|
|
|
|
|
|
#endif /* not REAL_ARITHMETIC */
|
|
|
|
|
#endif /* no XFmode */
|
|
|
|
|
|
|
|
|
|
void warning ();
|
|
|
|
|
extern int extra_warnings;
|
|
|
|
|
int ecmp (), enormlz (), eshift ();
|
|
|
|
|
int eisneg (), eisinf (), eisnan (), eiisinf (), eiisnan ();
|
|
|
|
|
void eadd (), esub (), emul (), ediv ();
|
|
|
|
|
void eshup1 (), eshup8 (), eshup6 (), eshdn1 (), eshdn8 (), eshdn6 ();
|
|
|
|
|
void eabs (), eneg (), emov (), eclear (), einfin (), efloor ();
|
|
|
|
|
void eldexp (), efrexp (), eifrac (), euifrac (), ltoe (), ultoe ();
|
|
|
|
|
void eround (), ereal_to_decimal (), eiinfin (), einan ();
|
|
|
|
|
void esqrt (), elog (), eexp (), etanh (), epow ();
|
|
|
|
|
void asctoe (), asctoe24 (), asctoe53 (), asctoe64 ();
|
|
|
|
|
void etoasc (), e24toasc (), e53toasc (), e64toasc ();
|
|
|
|
|
void etoe64 (), etoe53 (), etoe24 (), e64toe (), e53toe (), e24toe ();
|
|
|
|
|
void mtherr (), make_nan ();
|
|
|
|
|
void enan ();
|
|
|
|
|
extern unsigned EMUSHORT ezero[], ehalf[], eone[], etwo[];
|
|
|
|
|
extern unsigned EMUSHORT elog2[], esqrt2[];
|
|
|
|
|
|
|
|
|
|
/* Pack output array with 32-bit numbers obtained from
|
|
|
|
|
array containing 16-bit numbers, swapping ends if required. */
|
|
|
|
|
void
|
|
|
|
|
endian (e, x, mode)
|
|
|
|
|
unsigned EMUSHORT e[];
|
|
|
|
|
long x[];
|
|
|
|
|
enum machine_mode mode;
|
|
|
|
|
{
|
|
|
|
|
unsigned long th, t;
|
|
|
|
|
|
|
|
|
|
#if WORDS_BIG_ENDIAN
|
|
|
|
|
switch (mode)
|
|
|
|
|
{
|
|
|
|
|
|
|
|
|
|
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 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;
|
|
|
|
|
|
|
|
|
|
default:
|
|
|
|
|
abort ();
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
#else
|
|
|
|
|
|
|
|
|
|
/* Pack the output array without swapping. */
|
|
|
|
|
|
|
|
|
|
switch (mode)
|
|
|
|
|
{
|
|
|
|
|
|
|
|
|
|
case XFmode:
|
|
|
|
|
|
|
|
|
|
/* Pack the third long.
|
|
|
|
|
Each element of the input REAL_VALUE_TYPE array has 16 bit 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 SFmode:
|
|
|
|
|
|
|
|
|
|
/* pack the first long */
|
|
|
|
|
th = (unsigned long) e[1] & 0xffff;
|
|
|
|
|
t = (unsigned long) e[0] & 0xffff;
|
|
|
|
|
t |= th << 16;
|
|
|
|
|
x[0] = t;
|
|
|
|
|
break;
|
|
|
|
|
|
|
|
|
|
default:
|
|
|
|
|
abort ();
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
#endif
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* This is the implementation of the REAL_ARITHMETIC macro.
|
|
|
|
|
*/
|
|
|
|
|
void
|
|
|
|
|
earith (value, icode, r1, r2)
|
|
|
|
|
REAL_VALUE_TYPE *value;
|
|
|
|
|
int icode;
|
|
|
|
|
REAL_VALUE_TYPE *r1;
|
|
|
|
|
REAL_VALUE_TYPE *r2;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT d1[NE], d2[NE], v[NE];
|
|
|
|
|
enum tree_code code;
|
|
|
|
|
|
|
|
|
|
GET_REAL (r1, d1);
|
|
|
|
|
GET_REAL (r2, d2);
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
/* Return NaN input back to the caller. */
|
|
|
|
|
if (eisnan (d1))
|
|
|
|
|
{
|
|
|
|
|
PUT_REAL (d1, value);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
if (eisnan (d2))
|
|
|
|
|
{
|
|
|
|
|
PUT_REAL (d2, value);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
code = (enum tree_code) icode;
|
|
|
|
|
switch (code)
|
|
|
|
|
{
|
|
|
|
|
case PLUS_EXPR:
|
|
|
|
|
eadd (d2, d1, v);
|
|
|
|
|
break;
|
|
|
|
|
|
|
|
|
|
case MINUS_EXPR:
|
|
|
|
|
esub (d2, d1, v); /* d1 - d2 */
|
|
|
|
|
break;
|
|
|
|
|
|
|
|
|
|
case MULT_EXPR:
|
|
|
|
|
emul (d2, d1, v);
|
|
|
|
|
break;
|
|
|
|
|
|
|
|
|
|
case RDIV_EXPR:
|
|
|
|
|
#ifndef REAL_INFINITY
|
|
|
|
|
if (ecmp (d2, ezero) == 0)
|
|
|
|
|
{
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
enan (v);
|
|
|
|
|
break;
|
|
|
|
|
#else
|
|
|
|
|
abort ();
|
|
|
|
|
#endif
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
ediv (d2, d1, v); /* d1/d2 */
|
|
|
|
|
break;
|
|
|
|
|
|
|
|
|
|
case MIN_EXPR: /* min (d1,d2) */
|
|
|
|
|
if (ecmp (d1, d2) < 0)
|
|
|
|
|
emov (d1, v);
|
|
|
|
|
else
|
|
|
|
|
emov (d2, v);
|
|
|
|
|
break;
|
|
|
|
|
|
|
|
|
|
case MAX_EXPR: /* max (d1,d2) */
|
|
|
|
|
if (ecmp (d1, d2) > 0)
|
|
|
|
|
emov (d1, v);
|
|
|
|
|
else
|
|
|
|
|
emov (d2, v);
|
|
|
|
|
break;
|
|
|
|
|
default:
|
|
|
|
|
emov (ezero, v);
|
|
|
|
|
break;
|
|
|
|
|
}
|
|
|
|
|
PUT_REAL (v, value);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* Truncate REAL_VALUE_TYPE toward zero to signed HOST_WIDE_INT
|
|
|
|
|
* implements REAL_VALUE_RNDZINT (x) (etrunci (x))
|
|
|
|
|
*/
|
|
|
|
|
REAL_VALUE_TYPE
|
|
|
|
|
etrunci (x)
|
|
|
|
|
REAL_VALUE_TYPE x;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT f[NE], g[NE];
|
|
|
|
|
REAL_VALUE_TYPE r;
|
|
|
|
|
long l;
|
|
|
|
|
|
|
|
|
|
GET_REAL (&x, g);
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
if (eisnan (g))
|
|
|
|
|
return (x);
|
|
|
|
|
#endif
|
|
|
|
|
eifrac (g, &l, f);
|
|
|
|
|
ltoe (&l, g);
|
|
|
|
|
PUT_REAL (g, &r);
|
|
|
|
|
return (r);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* Truncate REAL_VALUE_TYPE toward zero to unsigned HOST_WIDE_INT
|
|
|
|
|
* implements REAL_VALUE_UNSIGNED_RNDZINT (x) (etruncui (x))
|
|
|
|
|
*/
|
|
|
|
|
REAL_VALUE_TYPE
|
|
|
|
|
etruncui (x)
|
|
|
|
|
REAL_VALUE_TYPE x;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT f[NE], g[NE];
|
|
|
|
|
REAL_VALUE_TYPE r;
|
|
|
|
|
unsigned long l;
|
|
|
|
|
|
|
|
|
|
GET_REAL (&x, g);
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
if (eisnan (g))
|
|
|
|
|
return (x);
|
|
|
|
|
#endif
|
|
|
|
|
euifrac (g, &l, f);
|
|
|
|
|
ultoe (&l, g);
|
|
|
|
|
PUT_REAL (g, &r);
|
|
|
|
|
return (r);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* This is the REAL_VALUE_ATOF function.
|
|
|
|
|
* It converts a decimal string to binary, rounding off
|
|
|
|
|
* as indicated by the machine_mode argument. Then it
|
|
|
|
|
* promotes the rounded value to REAL_VALUE_TYPE.
|
|
|
|
|
*/
|
|
|
|
|
REAL_VALUE_TYPE
|
|
|
|
|
ereal_atof (s, t)
|
|
|
|
|
char *s;
|
|
|
|
|
enum machine_mode t;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT tem[NE], e[NE];
|
|
|
|
|
REAL_VALUE_TYPE r;
|
|
|
|
|
|
|
|
|
|
switch (t)
|
|
|
|
|
{
|
|
|
|
|
case SFmode:
|
|
|
|
|
asctoe24 (s, tem);
|
|
|
|
|
e24toe (tem, e);
|
|
|
|
|
break;
|
|
|
|
|
case DFmode:
|
|
|
|
|
asctoe53 (s, tem);
|
|
|
|
|
e53toe (tem, e);
|
|
|
|
|
break;
|
|
|
|
|
case XFmode:
|
|
|
|
|
asctoe64 (s, tem);
|
|
|
|
|
e64toe (tem, e);
|
|
|
|
|
break;
|
|
|
|
|
default:
|
|
|
|
|
asctoe (s, e);
|
|
|
|
|
}
|
|
|
|
|
PUT_REAL (e, &r);
|
|
|
|
|
return (r);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* Expansion of REAL_NEGATE.
|
|
|
|
|
*/
|
|
|
|
|
REAL_VALUE_TYPE
|
|
|
|
|
ereal_negate (x)
|
|
|
|
|
REAL_VALUE_TYPE x;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT e[NE];
|
|
|
|
|
REAL_VALUE_TYPE r;
|
|
|
|
|
|
|
|
|
|
GET_REAL (&x, e);
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
if (eisnan (e))
|
|
|
|
|
return (x);
|
|
|
|
|
#endif
|
|
|
|
|
eneg (e);
|
|
|
|
|
PUT_REAL (e, &r);
|
|
|
|
|
return (r);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* Round real to int
|
|
|
|
|
* implements REAL_VALUE_FIX (x) (eroundi (x))
|
|
|
|
|
* The type of rounding is left unspecified by real.h.
|
|
|
|
|
* It is implemented here as round to nearest (add .5 and chop).
|
|
|
|
|
*/
|
|
|
|
|
int
|
|
|
|
|
eroundi (x)
|
|
|
|
|
REAL_VALUE_TYPE x;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT f[NE], g[NE];
|
|
|
|
|
EMULONG l;
|
|
|
|
|
|
|
|
|
|
GET_REAL (&x, f);
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
if (eisnan (f))
|
|
|
|
|
{
|
|
|
|
|
warning ("conversion from NaN to int");
|
|
|
|
|
return (-1);
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
eround (f, g);
|
|
|
|
|
eifrac (g, &l, f);
|
|
|
|
|
return ((int) l);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/* Round real to nearest unsigned int
|
|
|
|
|
* implements REAL_VALUE_UNSIGNED_FIX (x) ((unsigned int) eroundi (x))
|
|
|
|
|
* Negative input returns zero.
|
|
|
|
|
* The type of rounding is left unspecified by real.h.
|
|
|
|
|
* It is implemented here as round to nearest (add .5 and chop).
|
|
|
|
|
*/
|
|
|
|
|
unsigned int
|
|
|
|
|
eroundui (x)
|
|
|
|
|
REAL_VALUE_TYPE x;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT f[NE], g[NE];
|
|
|
|
|
unsigned EMULONG l;
|
|
|
|
|
|
|
|
|
|
GET_REAL (&x, f);
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
if (eisnan (f))
|
|
|
|
|
{
|
|
|
|
|
warning ("conversion from NaN to unsigned int");
|
|
|
|
|
return (-1);
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
eround (f, g);
|
|
|
|
|
euifrac (g, &l, f);
|
|
|
|
|
return ((unsigned int)l);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* REAL_VALUE_FROM_INT macro.
|
|
|
|
|
*/
|
|
|
|
|
void
|
|
|
|
|
ereal_from_int (d, i, j)
|
|
|
|
|
REAL_VALUE_TYPE *d;
|
|
|
|
|
long i, j;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT df[NE], dg[NE];
|
|
|
|
|
long low, high;
|
|
|
|
|
int sign;
|
|
|
|
|
|
|
|
|
|
sign = 0;
|
|
|
|
|
low = i;
|
|
|
|
|
if ((high = j) < 0)
|
|
|
|
|
{
|
|
|
|
|
sign = 1;
|
|
|
|
|
/* complement and add 1 */
|
|
|
|
|
high = ~high;
|
|
|
|
|
if (low)
|
|
|
|
|
low = -low;
|
|
|
|
|
else
|
|
|
|
|
high += 1;
|
|
|
|
|
}
|
|
|
|
|
eldexp (eone, HOST_BITS_PER_LONG, df);
|
|
|
|
|
ultoe (&high, dg);
|
|
|
|
|
emul (dg, df, dg);
|
|
|
|
|
ultoe (&low, df);
|
|
|
|
|
eadd (df, dg, dg);
|
|
|
|
|
if (sign)
|
|
|
|
|
eneg (dg);
|
|
|
|
|
PUT_REAL (dg, d);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* REAL_VALUE_FROM_UNSIGNED_INT macro.
|
|
|
|
|
*/
|
|
|
|
|
void
|
|
|
|
|
ereal_from_uint (d, i, j)
|
|
|
|
|
REAL_VALUE_TYPE *d;
|
|
|
|
|
unsigned long i, j;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT df[NE], dg[NE];
|
|
|
|
|
unsigned long low, high;
|
|
|
|
|
|
|
|
|
|
low = i;
|
|
|
|
|
high = j;
|
|
|
|
|
eldexp (eone, HOST_BITS_PER_LONG, df);
|
|
|
|
|
ultoe (&high, dg);
|
|
|
|
|
emul (dg, df, dg);
|
|
|
|
|
ultoe (&low, df);
|
|
|
|
|
eadd (df, dg, dg);
|
|
|
|
|
PUT_REAL (dg, d);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* REAL_VALUE_TO_INT macro
|
|
|
|
|
*/
|
|
|
|
|
void
|
|
|
|
|
ereal_to_int (low, high, rr)
|
|
|
|
|
long *low, *high;
|
|
|
|
|
REAL_VALUE_TYPE rr;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT d[NE], df[NE], dg[NE], dh[NE];
|
|
|
|
|
int s;
|
|
|
|
|
|
|
|
|
|
GET_REAL (&rr, d);
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
if (eisnan (d))
|
|
|
|
|
{
|
|
|
|
|
warning ("conversion from NaN to int");
|
|
|
|
|
*low = -1;
|
|
|
|
|
*high = -1;
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
/* convert positive value */
|
|
|
|
|
s = 0;
|
|
|
|
|
if (eisneg (d))
|
|
|
|
|
{
|
|
|
|
|
eneg (d);
|
|
|
|
|
s = 1;
|
|
|
|
|
}
|
|
|
|
|
eldexp (eone, HOST_BITS_PER_LONG, df);
|
|
|
|
|
ediv (df, d, dg); /* dg = d / 2^32 is the high word */
|
|
|
|
|
euifrac (dg, high, dh);
|
|
|
|
|
emul (df, dh, dg); /* fractional part is the low word */
|
|
|
|
|
euifrac (dg, low, dh);
|
|
|
|
|
if (s)
|
|
|
|
|
{
|
|
|
|
|
/* complement and add 1 */
|
|
|
|
|
*high = ~(*high);
|
|
|
|
|
if (*low)
|
|
|
|
|
*low = -(*low);
|
|
|
|
|
else
|
|
|
|
|
*high += 1;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* REAL_VALUE_LDEXP macro.
|
|
|
|
|
*/
|
|
|
|
|
REAL_VALUE_TYPE
|
|
|
|
|
ereal_ldexp (x, n)
|
|
|
|
|
REAL_VALUE_TYPE x;
|
|
|
|
|
int n;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT e[NE], y[NE];
|
|
|
|
|
REAL_VALUE_TYPE r;
|
|
|
|
|
|
|
|
|
|
GET_REAL (&x, e);
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
if (eisnan (e))
|
|
|
|
|
return (x);
|
|
|
|
|
#endif
|
|
|
|
|
eldexp (e, n, y);
|
|
|
|
|
PUT_REAL (y, &r);
|
|
|
|
|
return (r);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/* These routines are conditionally compiled because functions
|
|
|
|
|
* of the same names may be defined in fold-const.c. */
|
|
|
|
|
#ifdef REAL_ARITHMETIC
|
|
|
|
|
|
|
|
|
|
/* Check for infinity in a REAL_VALUE_TYPE. */
|
|
|
|
|
int
|
|
|
|
|
target_isinf (x)
|
|
|
|
|
REAL_VALUE_TYPE x;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT e[NE];
|
|
|
|
|
|
|
|
|
|
#ifdef INFINITY
|
|
|
|
|
GET_REAL (&x, e);
|
|
|
|
|
return (eisinf (e));
|
|
|
|
|
#else
|
|
|
|
|
return 0;
|
|
|
|
|
#endif
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* Check whether a REAL_VALUE_TYPE item is a NaN. */
|
|
|
|
|
|
|
|
|
|
int
|
|
|
|
|
target_isnan (x)
|
|
|
|
|
REAL_VALUE_TYPE x;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT e[NE];
|
|
|
|
|
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
GET_REAL (&x, e);
|
|
|
|
|
return (eisnan (e));
|
|
|
|
|
#else
|
|
|
|
|
return (0);
|
|
|
|
|
#endif
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* Check for a negative REAL_VALUE_TYPE number.
|
|
|
|
|
* this means strictly less than zero, not -0.
|
|
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
int
|
|
|
|
|
target_negative (x)
|
|
|
|
|
REAL_VALUE_TYPE x;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT e[NE];
|
|
|
|
|
|
|
|
|
|
GET_REAL (&x, e);
|
|
|
|
|
if (ecmp (e, ezero) == -1)
|
|
|
|
|
return (1);
|
|
|
|
|
return (0);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/* Expansion of REAL_VALUE_TRUNCATE.
|
|
|
|
|
* The result is in floating point, rounded to nearest or even.
|
|
|
|
|
*/
|
|
|
|
|
REAL_VALUE_TYPE
|
|
|
|
|
real_value_truncate (mode, arg)
|
|
|
|
|
enum machine_mode mode;
|
|
|
|
|
REAL_VALUE_TYPE arg;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT e[NE], t[NE];
|
|
|
|
|
REAL_VALUE_TYPE r;
|
|
|
|
|
|
|
|
|
|
GET_REAL (&arg, e);
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
if (eisnan (e))
|
|
|
|
|
return (arg);
|
|
|
|
|
#endif
|
|
|
|
|
eclear (t);
|
|
|
|
|
switch (mode)
|
|
|
|
|
{
|
|
|
|
|
case XFmode:
|
|
|
|
|
etoe64 (e, t);
|
|
|
|
|
e64toe (t, t);
|
|
|
|
|
break;
|
|
|
|
|
|
|
|
|
|
case DFmode:
|
|
|
|
|
etoe53 (e, t);
|
|
|
|
|
e53toe (t, t);
|
|
|
|
|
break;
|
|
|
|
|
|
|
|
|
|
case SFmode:
|
|
|
|
|
etoe24 (e, t);
|
|
|
|
|
e24toe (t, t);
|
|
|
|
|
break;
|
|
|
|
|
|
|
|
|
|
case SImode:
|
|
|
|
|
r = etrunci (e);
|
|
|
|
|
return (r);
|
|
|
|
|
|
|
|
|
|
default:
|
|
|
|
|
abort ();
|
|
|
|
|
}
|
|
|
|
|
PUT_REAL (t, &r);
|
|
|
|
|
return (r);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
#endif /* REAL_ARITHMETIC defined */
|
|
|
|
|
|
|
|
|
|
/* Target values are arrays of host longs. A long is guaranteed
|
|
|
|
|
to be at least 32 bits wide. */
|
|
|
|
|
void
|
|
|
|
|
etarldouble (r, l)
|
|
|
|
|
REAL_VALUE_TYPE r;
|
|
|
|
|
long l[];
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT e[NE];
|
|
|
|
|
|
|
|
|
|
GET_REAL (&r, e);
|
|
|
|
|
etoe64 (e, e);
|
|
|
|
|
endian (e, l, XFmode);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
etardouble (r, l)
|
|
|
|
|
REAL_VALUE_TYPE r;
|
|
|
|
|
long l[];
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT e[NE];
|
|
|
|
|
|
|
|
|
|
GET_REAL (&r, e);
|
|
|
|
|
etoe53 (e, e);
|
|
|
|
|
endian (e, l, DFmode);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
long
|
|
|
|
|
etarsingle (r)
|
|
|
|
|
REAL_VALUE_TYPE r;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT e[NE];
|
|
|
|
|
unsigned long l;
|
|
|
|
|
|
|
|
|
|
GET_REAL (&r, e);
|
|
|
|
|
etoe24 (e, e);
|
|
|
|
|
endian (e, &l, SFmode);
|
|
|
|
|
return ((long) l);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
ereal_to_decimal (x, s)
|
|
|
|
|
REAL_VALUE_TYPE x;
|
|
|
|
|
char *s;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT e[NE];
|
|
|
|
|
|
|
|
|
|
GET_REAL (&x, e);
|
|
|
|
|
etoasc (e, s, 20);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
int
|
|
|
|
|
ereal_cmp (x, y)
|
|
|
|
|
REAL_VALUE_TYPE x, y;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT ex[NE], ey[NE];
|
|
|
|
|
|
|
|
|
|
GET_REAL (&x, ex);
|
|
|
|
|
GET_REAL (&y, ey);
|
|
|
|
|
return (ecmp (ex, ey));
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
int
|
|
|
|
|
ereal_isneg (x)
|
|
|
|
|
REAL_VALUE_TYPE x;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT ex[NE];
|
|
|
|
|
|
|
|
|
|
GET_REAL (&x, ex);
|
|
|
|
|
return (eisneg (ex));
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/* End of REAL_ARITHMETIC interface */
|
|
|
|
|
|
|
|
|
|
/* ieee.c
|
|
|
|
|
*
|
|
|
|
|
* Extended precision IEEE binary floating point arithmetic routines
|
|
|
|
|
*
|
|
|
|
|
* Numbers are stored in C language as arrays of 16-bit unsigned
|
|
|
|
|
* short integers. The arguments of the routines are pointers to
|
|
|
|
|
* the arrays.
|
|
|
|
|
*
|
|
|
|
|
*
|
|
|
|
|
* External e type data structure, simulates Intel 8087 chip
|
|
|
|
|
* temporary real format but possibly with a larger significand:
|
|
|
|
|
*
|
|
|
|
|
* NE-1 significand words (least significant word first,
|
|
|
|
|
* most significant bit is normally set)
|
|
|
|
|
* exponent (value = EXONE for 1.0,
|
|
|
|
|
* top bit is the sign)
|
|
|
|
|
*
|
|
|
|
|
*
|
|
|
|
|
* Internal 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)
|
|
|
|
|
* ei[2] high guard word (always zero after normalization)
|
|
|
|
|
* ei[3]
|
|
|
|
|
* to ei[NI-2] significand (NI-4 significand words,
|
|
|
|
|
* most significant word first,
|
|
|
|
|
* most significant bit is set)
|
|
|
|
|
* ei[NI-1] low guard word (0x8000 bit is rounding place)
|
|
|
|
|
*
|
|
|
|
|
*
|
|
|
|
|
*
|
|
|
|
|
* Routines for external format numbers
|
|
|
|
|
*
|
|
|
|
|
* asctoe (string, e) ASCII string to extended double e type
|
|
|
|
|
* asctoe64 (string, &d) ASCII string to long double
|
|
|
|
|
* asctoe53 (string, &d) ASCII string to double
|
|
|
|
|
* asctoe24 (string, &f) ASCII string to single
|
|
|
|
|
* asctoeg (string, e, prec) ASCII string to specified precision
|
|
|
|
|
* e24toe (&f, e) IEEE single precision to e type
|
|
|
|
|
* e53toe (&d, e) IEEE double precision to e type
|
|
|
|
|
* e64toe (&d, e) IEEE long double precision to e type
|
|
|
|
|
* eabs (e) absolute value
|
|
|
|
|
* eadd (a, b, c) c = b + a
|
|
|
|
|
* eclear (e) e = 0
|
|
|
|
|
* ecmp (a, b) Returns 1 if a > b, 0 if a == b,
|
|
|
|
|
* -1 if a < b, -2 if either a or b is a NaN.
|
|
|
|
|
* ediv (a, b, c) c = b / a
|
|
|
|
|
* efloor (a, b) truncate to integer, toward -infinity
|
|
|
|
|
* efrexp (a, exp, s) extract exponent and significand
|
|
|
|
|
* eifrac (e, &l, frac) e to long integer and e type fraction
|
|
|
|
|
* euifrac (e, &l, frac) e to unsigned long integer and e type fraction
|
|
|
|
|
* einfin (e) set e to infinity, leaving its sign alone
|
|
|
|
|
* eldexp (a, n, b) multiply by 2**n
|
|
|
|
|
* emov (a, b) b = a
|
|
|
|
|
* emul (a, b, c) c = b * a
|
|
|
|
|
* eneg (e) e = -e
|
|
|
|
|
* eround (a, b) b = nearest integer value to a
|
|
|
|
|
* esub (a, b, c) c = b - a
|
|
|
|
|
* e24toasc (&f, str, n) single to ASCII string, n digits after decimal
|
|
|
|
|
* e53toasc (&d, str, n) double to ASCII string, n digits after decimal
|
|
|
|
|
* e64toasc (&d, str, n) long double to ASCII string
|
|
|
|
|
* etoasc (e, str, n) e to ASCII string, n digits after decimal
|
|
|
|
|
* etoe24 (e, &f) convert e type to IEEE single precision
|
|
|
|
|
* etoe53 (e, &d) convert e type to IEEE double precision
|
|
|
|
|
* etoe64 (e, &d) convert e type to IEEE long double precision
|
|
|
|
|
* ltoe (&l, e) long (32 bit) integer to e type
|
|
|
|
|
* ultoe (&l, e) unsigned long (32 bit) integer to e type
|
|
|
|
|
* eisneg (e) 1 if sign bit of e != 0, else 0
|
|
|
|
|
* eisinf (e) 1 if e has maximum exponent (non-IEEE)
|
|
|
|
|
* or is infinite (IEEE)
|
|
|
|
|
* eisnan (e) 1 if e is a NaN
|
|
|
|
|
*
|
|
|
|
|
*
|
|
|
|
|
* Routines for internal format numbers
|
|
|
|
|
*
|
|
|
|
|
* eaddm (ai, bi) add significands, bi = bi + ai
|
|
|
|
|
* ecleaz (ei) ei = 0
|
|
|
|
|
* ecleazs (ei) set ei = 0 but leave its sign alone
|
|
|
|
|
* ecmpm (ai, bi) compare significands, return 1, 0, or -1
|
|
|
|
|
* edivm (ai, bi) divide significands, bi = bi / ai
|
|
|
|
|
* emdnorm (ai,l,s,exp) normalize and round off
|
|
|
|
|
* emovi (a, ai) convert external a to internal ai
|
|
|
|
|
* emovo (ai, a) convert internal ai to external a
|
|
|
|
|
* emovz (ai, bi) bi = ai, low guard word of bi = 0
|
|
|
|
|
* emulm (ai, bi) multiply significands, bi = bi * ai
|
|
|
|
|
* enormlz (ei) left-justify the significand
|
|
|
|
|
* eshdn1 (ai) shift significand and guards down 1 bit
|
|
|
|
|
* eshdn8 (ai) shift down 8 bits
|
|
|
|
|
* eshdn6 (ai) shift down 16 bits
|
|
|
|
|
* eshift (ai, n) shift ai n bits up (or down if n < 0)
|
|
|
|
|
* eshup1 (ai) shift significand and guards up 1 bit
|
|
|
|
|
* eshup8 (ai) shift up 8 bits
|
|
|
|
|
* eshup6 (ai) shift up 16 bits
|
|
|
|
|
* esubm (ai, bi) subtract significands, bi = bi - ai
|
|
|
|
|
* eiisinf (ai) 1 if infinite
|
|
|
|
|
* eiisnan (ai) 1 if a NaN
|
|
|
|
|
* einan (ai) set ai = NaN
|
|
|
|
|
* eiinfin (ai) set ai = infinity
|
|
|
|
|
*
|
|
|
|
|
*
|
|
|
|
|
* The result is always normalized and rounded to NI-4 word precision
|
|
|
|
|
* after each arithmetic operation.
|
|
|
|
|
*
|
|
|
|
|
* Exception flags are NOT fully supported.
|
|
|
|
|
*
|
|
|
|
|
* Signaling NaN's are NOT supported; they are treated the same
|
|
|
|
|
* as quiet NaN's.
|
|
|
|
|
*
|
|
|
|
|
* Define INFINITY for support of infinity; otherwise a
|
|
|
|
|
* saturation arithmetic is implemented.
|
|
|
|
|
*
|
|
|
|
|
* Define NANS for support of Not-a-Number items; otherwise the
|
|
|
|
|
* arithmetic will never produce a NaN output, and might be confused
|
|
|
|
|
* by a NaN input.
|
|
|
|
|
* If NaN's are supported, the output of `ecmp (a,b)' is -2 if
|
|
|
|
|
* either a or b is a NaN. This means asking `if (ecmp (a,b) < 0)'
|
|
|
|
|
* may not be legitimate. Use `if (ecmp (a,b) == -1)' for `less than'
|
|
|
|
|
* if in doubt.
|
|
|
|
|
*
|
|
|
|
|
* Denormals are always supported here where appropriate (e.g., not
|
|
|
|
|
* for conversion to DEC numbers).
|
|
|
|
|
*
|
|
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* mconf.h
|
|
|
|
|
*
|
|
|
|
|
* Common include file for math routines
|
|
|
|
|
*
|
|
|
|
|
*
|
|
|
|
|
*
|
|
|
|
|
* SYNOPSIS:
|
|
|
|
|
*
|
|
|
|
|
* #include "mconf.h"
|
|
|
|
|
*
|
|
|
|
|
*
|
|
|
|
|
*
|
|
|
|
|
* DESCRIPTION:
|
|
|
|
|
*
|
|
|
|
|
* This file contains definitions for error codes that are
|
|
|
|
|
* passed to the common error handling routine mtherr
|
|
|
|
|
* (which see).
|
|
|
|
|
*
|
|
|
|
|
* The file also includes a conditional assembly definition
|
|
|
|
|
* for the type of computer arithmetic (Intel IEEE, DEC, Motorola
|
|
|
|
|
* IEEE, or UNKnown).
|
|
|
|
|
*
|
|
|
|
|
* For Digital Equipment PDP-11 and VAX computers, certain
|
|
|
|
|
* IBM systems, and others that use numbers with a 56-bit
|
|
|
|
|
* significand, the symbol DEC should be defined. In this
|
|
|
|
|
* mode, most floating point constants are given as arrays
|
|
|
|
|
* of octal integers to eliminate decimal to binary conversion
|
|
|
|
|
* errors that might be introduced by the compiler.
|
|
|
|
|
*
|
|
|
|
|
* 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.
|
|
|
|
|
* These numbers have 53-bit significands. In this mode, constants
|
|
|
|
|
* are provided as arrays of hexadecimal 16 bit integers.
|
|
|
|
|
*
|
|
|
|
|
* To accommodate other types of computer arithmetic, all
|
|
|
|
|
* constants are also provided in a normal decimal radix
|
|
|
|
|
* which one can hope are correctly converted to a suitable
|
|
|
|
|
* format by the available C language compiler. To invoke
|
|
|
|
|
* this mode, the symbol UNK is defined.
|
|
|
|
|
*
|
|
|
|
|
* An important difference among these modes is a predefined
|
|
|
|
|
* set of machine arithmetic constants for each. The numbers
|
|
|
|
|
* MACHEP (the machine roundoff error), MAXNUM (largest number
|
|
|
|
|
* represented), and several other parameters are preset by
|
|
|
|
|
* the configuration symbol. Check the file const.c to
|
|
|
|
|
* 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.
|
|
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
/* Constant definitions for math error conditions. */
|
|
|
|
|
|
|
|
|
|
#define DOMAIN 1 /* argument domain error */
|
|
|
|
|
#define SING 2 /* argument singularity */
|
|
|
|
|
#define OVERFLOW 3 /* overflow range error */
|
|
|
|
|
#define UNDERFLOW 4 /* underflow range error */
|
|
|
|
|
#define TLOSS 5 /* total loss of precision */
|
|
|
|
|
#define PLOSS 6 /* partial loss of precision */
|
|
|
|
|
#define INVALID 7 /* NaN-producing operation */
|
|
|
|
|
|
|
|
|
|
/* e type constants used by high precision check routines */
|
|
|
|
|
|
|
|
|
|
/*include "ehead.h"*/
|
|
|
|
|
/* 0.0 */
|
|
|
|
|
unsigned EMUSHORT ezero[NE] =
|
|
|
|
|
{
|
|
|
|
|
0, 0000000, 0000000, 0000000, 0000000, 0000000,};
|
|
|
|
|
extern unsigned EMUSHORT ezero[];
|
|
|
|
|
|
|
|
|
|
/* 5.0E-1 */
|
|
|
|
|
unsigned EMUSHORT ehalf[NE] =
|
|
|
|
|
{
|
|
|
|
|
0, 0000000, 0000000, 0000000, 0100000, 0x3ffe,};
|
|
|
|
|
extern unsigned EMUSHORT ehalf[];
|
|
|
|
|
|
|
|
|
|
/* 1.0E0 */
|
|
|
|
|
unsigned EMUSHORT eone[NE] =
|
|
|
|
|
{
|
|
|
|
|
0, 0000000, 0000000, 0000000, 0100000, 0x3fff,};
|
|
|
|
|
extern unsigned EMUSHORT eone[];
|
|
|
|
|
|
|
|
|
|
/* 2.0E0 */
|
|
|
|
|
unsigned EMUSHORT etwo[NE] =
|
|
|
|
|
{
|
|
|
|
|
0, 0000000, 0000000, 0000000, 0100000, 0040000,};
|
|
|
|
|
extern unsigned EMUSHORT etwo[];
|
|
|
|
|
|
|
|
|
|
/* 3.2E1 */
|
|
|
|
|
unsigned EMUSHORT e32[NE] =
|
|
|
|
|
{
|
|
|
|
|
0, 0000000, 0000000, 0000000, 0100000, 0040004,};
|
|
|
|
|
extern unsigned EMUSHORT e32[];
|
|
|
|
|
|
|
|
|
|
/* 6.93147180559945309417232121458176568075500134360255E-1 */
|
|
|
|
|
unsigned EMUSHORT elog2[NE] =
|
|
|
|
|
{
|
|
|
|
|
0xc9e4, 0x79ab, 0150717, 0013767, 0130562, 0x3ffe,};
|
|
|
|
|
extern unsigned EMUSHORT elog2[];
|
|
|
|
|
|
|
|
|
|
/* 1.41421356237309504880168872420969807856967187537695E0 */
|
|
|
|
|
unsigned EMUSHORT esqrt2[NE] =
|
|
|
|
|
{
|
|
|
|
|
0x597e, 0x6484, 0174736, 0171463, 0132404, 0x3fff,};
|
|
|
|
|
extern unsigned EMUSHORT esqrt2[];
|
|
|
|
|
|
|
|
|
|
/* 2/sqrt (PI) =
|
|
|
|
|
* 1.12837916709551257389615890312154517168810125865800E0 */
|
|
|
|
|
unsigned EMUSHORT eoneopi[NE] =
|
|
|
|
|
{
|
|
|
|
|
0x71d5, 0x688d, 0012333, 0135202, 0110156, 0x3fff,};
|
|
|
|
|
extern unsigned EMUSHORT eoneopi[];
|
|
|
|
|
|
|
|
|
|
/* 3.14159265358979323846264338327950288419716939937511E0 */
|
|
|
|
|
unsigned EMUSHORT epi[NE] =
|
|
|
|
|
{
|
|
|
|
|
0xc4c6, 0xc234, 0020550, 0155242, 0144417, 0040000,};
|
|
|
|
|
extern unsigned EMUSHORT epi[];
|
|
|
|
|
|
|
|
|
|
/* 5.7721566490153286060651209008240243104215933593992E-1 */
|
|
|
|
|
unsigned EMUSHORT eeul[NE] =
|
|
|
|
|
{
|
|
|
|
|
0xd1be, 0xc7a4, 0076660, 0063743, 0111704, 0x3ffe,};
|
|
|
|
|
extern unsigned EMUSHORT eeul[];
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
include "ehead.h"
|
|
|
|
|
include "mconf.h"
|
|
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* Control register for rounding precision.
|
|
|
|
|
* This can be set to 80 (if NE=6), 64, 56, 53, or 24 bits.
|
|
|
|
|
*/
|
|
|
|
|
int rndprc = NBITS;
|
|
|
|
|
extern int rndprc;
|
|
|
|
|
|
|
|
|
|
void eaddm (), esubm (), emdnorm (), asctoeg ();
|
|
|
|
|
static void toe24 (), toe53 (), toe64 ();
|
|
|
|
|
void eremain (), einit (), eiremain ();
|
|
|
|
|
int ecmpm (), edivm (), emulm ();
|
|
|
|
|
void emovi (), emovo (), emovz (), ecleaz (), ecleazs (), eadd1 ();
|
|
|
|
|
void etodec (), todec (), dectoe ();
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
einit ()
|
|
|
|
|
{
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; Clear out entire external format number.
|
|
|
|
|
;
|
|
|
|
|
; unsigned EMUSHORT x[];
|
|
|
|
|
; eclear (x);
|
|
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
eclear (x)
|
|
|
|
|
register unsigned EMUSHORT *x;
|
|
|
|
|
{
|
|
|
|
|
register int i;
|
|
|
|
|
|
|
|
|
|
for (i = 0; i < NE; i++)
|
|
|
|
|
*x++ = 0;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* Move external format number from a to b.
|
|
|
|
|
*
|
|
|
|
|
* emov (a, b);
|
|
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
emov (a, b)
|
|
|
|
|
register unsigned EMUSHORT *a, *b;
|
|
|
|
|
{
|
|
|
|
|
register int i;
|
|
|
|
|
|
|
|
|
|
for (i = 0; i < NE; i++)
|
|
|
|
|
*b++ = *a++;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; Absolute value of external format number
|
|
|
|
|
;
|
|
|
|
|
; EMUSHORT x[NE];
|
|
|
|
|
; eabs (x);
|
|
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
eabs (x)
|
|
|
|
|
unsigned EMUSHORT x[]; /* x is the memory address of a short */
|
|
|
|
|
{
|
|
|
|
|
|
|
|
|
|
x[NE - 1] &= 0x7fff; /* sign is top bit of last word of external format */
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; Negate external format number
|
|
|
|
|
;
|
|
|
|
|
; unsigned EMUSHORT x[NE];
|
|
|
|
|
; eneg (x);
|
|
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
eneg (x)
|
|
|
|
|
unsigned EMUSHORT x[];
|
|
|
|
|
{
|
|
|
|
|
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
if (eisnan (x))
|
|
|
|
|
return;
|
|
|
|
|
#endif
|
|
|
|
|
x[NE - 1] ^= 0x8000; /* Toggle the sign bit */
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* Return 1 if external format number is negative,
|
|
|
|
|
* else return zero, including when it is a NaN.
|
|
|
|
|
*/
|
|
|
|
|
int
|
|
|
|
|
eisneg (x)
|
|
|
|
|
unsigned EMUSHORT x[];
|
|
|
|
|
{
|
|
|
|
|
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
if (eisnan (x))
|
|
|
|
|
return (0);
|
|
|
|
|
#endif
|
|
|
|
|
if (x[NE - 1] & 0x8000)
|
|
|
|
|
return (1);
|
|
|
|
|
else
|
|
|
|
|
return (0);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* Return 1 if external format number is infinity.
|
|
|
|
|
* else return zero.
|
|
|
|
|
*/
|
|
|
|
|
int
|
|
|
|
|
eisinf (x)
|
|
|
|
|
unsigned EMUSHORT x[];
|
|
|
|
|
{
|
|
|
|
|
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
if (eisnan (x))
|
|
|
|
|
return (0);
|
|
|
|
|
#endif
|
|
|
|
|
if ((x[NE - 1] & 0x7fff) == 0x7fff)
|
|
|
|
|
return (1);
|
|
|
|
|
else
|
|
|
|
|
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. */
|
|
|
|
|
|
|
|
|
|
int
|
|
|
|
|
eisnan (x)
|
|
|
|
|
unsigned EMUSHORT x[];
|
|
|
|
|
{
|
|
|
|
|
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
int i;
|
|
|
|
|
/* NaN has maximum exponent */
|
|
|
|
|
if ((x[NE - 1] & 0x7fff) != 0x7fff)
|
|
|
|
|
return (0);
|
|
|
|
|
/* ... and non-zero significand field. */
|
|
|
|
|
for (i = 0; i < NE - 1; i++)
|
|
|
|
|
{
|
|
|
|
|
if (*x++ != 0)
|
|
|
|
|
return (1);
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
return (0);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/* Fill external format number with infinity pattern (IEEE)
|
|
|
|
|
or largest possible number (non-IEEE).
|
|
|
|
|
Before calling einfin, you should either call eclear
|
|
|
|
|
or set up the sign bit by hand. */
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
einfin (x)
|
|
|
|
|
register unsigned EMUSHORT *x;
|
|
|
|
|
{
|
|
|
|
|
register int i;
|
|
|
|
|
|
|
|
|
|
#ifdef INFINITY
|
|
|
|
|
for (i = 0; i < NE - 1; i++)
|
|
|
|
|
*x++ = 0;
|
|
|
|
|
*x |= 32767;
|
|
|
|
|
#else
|
|
|
|
|
for (i = 0; i < NE - 1; i++)
|
|
|
|
|
*x++ = 0xffff;
|
|
|
|
|
*x |= 32766;
|
|
|
|
|
if (rndprc < NBITS)
|
|
|
|
|
{
|
|
|
|
|
if (rndprc == 64)
|
|
|
|
|
{
|
|
|
|
|
*(x - 5) = 0;
|
|
|
|
|
}
|
|
|
|
|
if (rndprc == 53)
|
|
|
|
|
{
|
|
|
|
|
*(x - 4) = 0xf800;
|
|
|
|
|
}
|
|
|
|
|
else
|
|
|
|
|
{
|
|
|
|
|
*(x - 4) = 0;
|
|
|
|
|
*(x - 3) = 0;
|
|
|
|
|
*(x - 2) = 0xff00;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
#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. */
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
enan (x)
|
|
|
|
|
register unsigned EMUSHORT *x;
|
|
|
|
|
{
|
|
|
|
|
register int i;
|
|
|
|
|
|
|
|
|
|
for (i = 0; i < NE - 2; i++)
|
|
|
|
|
*x++ = 0;
|
|
|
|
|
*x++ = 0xc000;
|
|
|
|
|
*x = 0x7fff;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* Move in external format number,
|
|
|
|
|
* converting it to internal format.
|
|
|
|
|
*/
|
|
|
|
|
void
|
|
|
|
|
emovi (a, b)
|
|
|
|
|
unsigned EMUSHORT *a, *b;
|
|
|
|
|
{
|
|
|
|
|
register unsigned EMUSHORT *p, *q;
|
|
|
|
|
int i;
|
|
|
|
|
|
|
|
|
|
q = b;
|
|
|
|
|
p = a + (NE - 1); /* point to last word of external number */
|
|
|
|
|
/* get the sign bit */
|
|
|
|
|
if (*p & 0x8000)
|
|
|
|
|
*q++ = 0xffff;
|
|
|
|
|
else
|
|
|
|
|
*q++ = 0;
|
|
|
|
|
/* get the exponent */
|
|
|
|
|
*q = *p--;
|
|
|
|
|
*q++ &= 0x7fff; /* delete the sign bit */
|
|
|
|
|
#ifdef INFINITY
|
|
|
|
|
if ((*(q - 1) & 0x7fff) == 0x7fff)
|
|
|
|
|
{
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
if (eisnan (a))
|
|
|
|
|
{
|
|
|
|
|
*q++ = 0;
|
|
|
|
|
for (i = 3; i < NI; i++)
|
|
|
|
|
*q++ = *p--;
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
for (i = 2; i < NI; i++)
|
|
|
|
|
*q++ = 0;
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
/* clear high guard word */
|
|
|
|
|
*q++ = 0;
|
|
|
|
|
/* move in the significand */
|
|
|
|
|
for (i = 0; i < NE - 1; i++)
|
|
|
|
|
*q++ = *p--;
|
|
|
|
|
/* clear low guard word */
|
|
|
|
|
*q = 0;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* Move internal format number out,
|
|
|
|
|
* converting it to external format.
|
|
|
|
|
*/
|
|
|
|
|
void
|
|
|
|
|
emovo (a, b)
|
|
|
|
|
unsigned EMUSHORT *a, *b;
|
|
|
|
|
{
|
|
|
|
|
register unsigned EMUSHORT *p, *q;
|
|
|
|
|
unsigned EMUSHORT i;
|
|
|
|
|
|
|
|
|
|
p = a;
|
|
|
|
|
q = b + (NE - 1); /* point to output exponent */
|
|
|
|
|
/* combine sign and exponent */
|
|
|
|
|
i = *p++;
|
|
|
|
|
if (i)
|
|
|
|
|
*q-- = *p++ | 0x8000;
|
|
|
|
|
else
|
|
|
|
|
*q-- = *p++;
|
|
|
|
|
#ifdef INFINITY
|
|
|
|
|
if (*(p - 1) == 0x7fff)
|
|
|
|
|
{
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
if (eiisnan (a))
|
|
|
|
|
{
|
|
|
|
|
enan (b);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
einfin (b);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
/* skip over guard word */
|
|
|
|
|
++p;
|
|
|
|
|
/* move the significand */
|
|
|
|
|
for (i = 0; i < NE - 1; i++)
|
|
|
|
|
*q-- = *p++;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* Clear out internal format number.
|
|
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
ecleaz (xi)
|
|
|
|
|
register unsigned EMUSHORT *xi;
|
|
|
|
|
{
|
|
|
|
|
register int i;
|
|
|
|
|
|
|
|
|
|
for (i = 0; i < NI; i++)
|
|
|
|
|
*xi++ = 0;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* same, but don't touch the sign. */
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
ecleazs (xi)
|
|
|
|
|
register unsigned EMUSHORT *xi;
|
|
|
|
|
{
|
|
|
|
|
register int i;
|
|
|
|
|
|
|
|
|
|
++xi;
|
|
|
|
|
for (i = 0; i < NI - 1; i++)
|
|
|
|
|
*xi++ = 0;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* Move internal format number from a to b.
|
|
|
|
|
*/
|
|
|
|
|
void
|
|
|
|
|
emovz (a, b)
|
|
|
|
|
register unsigned EMUSHORT *a, *b;
|
|
|
|
|
{
|
|
|
|
|
register int i;
|
|
|
|
|
|
|
|
|
|
for (i = 0; i < NI - 1; i++)
|
|
|
|
|
*b++ = *a++;
|
|
|
|
|
/* clear low guard word */
|
|
|
|
|
*b = 0;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/* Generate internal format NaN.
|
|
|
|
|
The explicit pattern for this is maximum exponent and
|
|
|
|
|
top two significand bits set. */
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
einan (x)
|
|
|
|
|
unsigned EMUSHORT x[];
|
|
|
|
|
{
|
|
|
|
|
|
|
|
|
|
ecleaz (x);
|
|
|
|
|
x[E] = 0x7fff;
|
|
|
|
|
x[M + 1] = 0xc000;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/* Return nonzero if internal format number is a NaN. */
|
|
|
|
|
|
|
|
|
|
int
|
|
|
|
|
eiisnan (x)
|
|
|
|
|
unsigned EMUSHORT x[];
|
|
|
|
|
{
|
|
|
|
|
int i;
|
|
|
|
|
|
|
|
|
|
if ((x[E] & 0x7fff) == 0x7fff)
|
|
|
|
|
{
|
|
|
|
|
for (i = M + 1; i < NI; i++)
|
|
|
|
|
{
|
|
|
|
|
if (x[i] != 0)
|
|
|
|
|
return (1);
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
return (0);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/* Fill internal format number with infinity pattern.
|
|
|
|
|
This has maximum exponent and significand all zeros. */
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
eiinfin (x)
|
|
|
|
|
unsigned EMUSHORT x[];
|
|
|
|
|
{
|
|
|
|
|
|
|
|
|
|
ecleaz (x);
|
|
|
|
|
x[E] = 0x7fff;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/* Return nonzero if internal format number is infinite. */
|
|
|
|
|
|
|
|
|
|
int
|
|
|
|
|
eiisinf (x)
|
|
|
|
|
unsigned EMUSHORT x[];
|
|
|
|
|
{
|
|
|
|
|
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
if (eiisnan (x))
|
|
|
|
|
return (0);
|
|
|
|
|
#endif
|
|
|
|
|
if ((x[E] & 0x7fff) == 0x7fff)
|
|
|
|
|
return (1);
|
|
|
|
|
return (0);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; Compare significands of numbers in internal format.
|
|
|
|
|
; Guard words are included in the comparison.
|
|
|
|
|
;
|
|
|
|
|
; unsigned EMUSHORT a[NI], b[NI];
|
|
|
|
|
; cmpm (a, b);
|
|
|
|
|
;
|
|
|
|
|
; for the significands:
|
|
|
|
|
; returns +1 if a > b
|
|
|
|
|
; 0 if a == b
|
|
|
|
|
; -1 if a < b
|
|
|
|
|
*/
|
|
|
|
|
int
|
|
|
|
|
ecmpm (a, b)
|
|
|
|
|
register unsigned EMUSHORT *a, *b;
|
|
|
|
|
{
|
|
|
|
|
int i;
|
|
|
|
|
|
|
|
|
|
a += M; /* skip up to significand area */
|
|
|
|
|
b += M;
|
|
|
|
|
for (i = M; i < NI; i++)
|
|
|
|
|
{
|
|
|
|
|
if (*a++ != *b++)
|
|
|
|
|
goto difrnt;
|
|
|
|
|
}
|
|
|
|
|
return (0);
|
|
|
|
|
|
|
|
|
|
difrnt:
|
|
|
|
|
if (*(--a) > *(--b))
|
|
|
|
|
return (1);
|
|
|
|
|
else
|
|
|
|
|
return (-1);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; Shift significand down by 1 bit
|
|
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
eshdn1 (x)
|
|
|
|
|
register unsigned EMUSHORT *x;
|
|
|
|
|
{
|
|
|
|
|
register unsigned EMUSHORT bits;
|
|
|
|
|
int i;
|
|
|
|
|
|
|
|
|
|
x += M; /* point to significand area */
|
|
|
|
|
|
|
|
|
|
bits = 0;
|
|
|
|
|
for (i = M; i < NI; i++)
|
|
|
|
|
{
|
|
|
|
|
if (*x & 1)
|
|
|
|
|
bits |= 1;
|
|
|
|
|
*x >>= 1;
|
|
|
|
|
if (bits & 2)
|
|
|
|
|
*x |= 0x8000;
|
|
|
|
|
bits <<= 1;
|
|
|
|
|
++x;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; Shift significand up by 1 bit
|
|
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
eshup1 (x)
|
|
|
|
|
register unsigned EMUSHORT *x;
|
|
|
|
|
{
|
|
|
|
|
register unsigned EMUSHORT bits;
|
|
|
|
|
int i;
|
|
|
|
|
|
|
|
|
|
x += NI - 1;
|
|
|
|
|
bits = 0;
|
|
|
|
|
|
|
|
|
|
for (i = M; i < NI; i++)
|
|
|
|
|
{
|
|
|
|
|
if (*x & 0x8000)
|
|
|
|
|
bits |= 1;
|
|
|
|
|
*x <<= 1;
|
|
|
|
|
if (bits & 2)
|
|
|
|
|
*x |= 1;
|
|
|
|
|
bits <<= 1;
|
|
|
|
|
--x;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; Shift significand down by 8 bits
|
|
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
eshdn8 (x)
|
|
|
|
|
register unsigned EMUSHORT *x;
|
|
|
|
|
{
|
|
|
|
|
register unsigned EMUSHORT newbyt, oldbyt;
|
|
|
|
|
int i;
|
|
|
|
|
|
|
|
|
|
x += M;
|
|
|
|
|
oldbyt = 0;
|
|
|
|
|
for (i = M; i < NI; i++)
|
|
|
|
|
{
|
|
|
|
|
newbyt = *x << 8;
|
|
|
|
|
*x >>= 8;
|
|
|
|
|
*x |= oldbyt;
|
|
|
|
|
oldbyt = newbyt;
|
|
|
|
|
++x;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; Shift significand up by 8 bits
|
|
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
eshup8 (x)
|
|
|
|
|
register unsigned EMUSHORT *x;
|
|
|
|
|
{
|
|
|
|
|
int i;
|
|
|
|
|
register unsigned EMUSHORT newbyt, oldbyt;
|
|
|
|
|
|
|
|
|
|
x += NI - 1;
|
|
|
|
|
oldbyt = 0;
|
|
|
|
|
|
|
|
|
|
for (i = M; i < NI; i++)
|
|
|
|
|
{
|
|
|
|
|
newbyt = *x >> 8;
|
|
|
|
|
*x <<= 8;
|
|
|
|
|
*x |= oldbyt;
|
|
|
|
|
oldbyt = newbyt;
|
|
|
|
|
--x;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; Shift significand up by 16 bits
|
|
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
eshup6 (x)
|
|
|
|
|
register unsigned EMUSHORT *x;
|
|
|
|
|
{
|
|
|
|
|
int i;
|
|
|
|
|
register unsigned EMUSHORT *p;
|
|
|
|
|
|
|
|
|
|
p = x + M;
|
|
|
|
|
x += M + 1;
|
|
|
|
|
|
|
|
|
|
for (i = M; i < NI - 1; i++)
|
|
|
|
|
*p++ = *x++;
|
|
|
|
|
|
|
|
|
|
*p = 0;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; Shift significand down by 16 bits
|
|
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
eshdn6 (x)
|
|
|
|
|
register unsigned EMUSHORT *x;
|
|
|
|
|
{
|
|
|
|
|
int i;
|
|
|
|
|
register unsigned EMUSHORT *p;
|
|
|
|
|
|
|
|
|
|
x += NI - 1;
|
|
|
|
|
p = x + 1;
|
|
|
|
|
|
|
|
|
|
for (i = M; i < NI - 1; i++)
|
|
|
|
|
*(--p) = *(--x);
|
|
|
|
|
|
|
|
|
|
*(--p) = 0;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; Add significands
|
|
|
|
|
; x + y replaces y
|
|
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
eaddm (x, y)
|
|
|
|
|
unsigned EMUSHORT *x, *y;
|
|
|
|
|
{
|
|
|
|
|
register unsigned EMULONG a;
|
|
|
|
|
int i;
|
|
|
|
|
unsigned int carry;
|
|
|
|
|
|
|
|
|
|
x += NI - 1;
|
|
|
|
|
y += NI - 1;
|
|
|
|
|
carry = 0;
|
|
|
|
|
for (i = M; i < NI; i++)
|
|
|
|
|
{
|
|
|
|
|
a = (unsigned EMULONG) (*x) + (unsigned EMULONG) (*y) + carry;
|
|
|
|
|
if (a & 0x10000)
|
|
|
|
|
carry = 1;
|
|
|
|
|
else
|
|
|
|
|
carry = 0;
|
|
|
|
|
*y = (unsigned EMUSHORT) a;
|
|
|
|
|
--x;
|
|
|
|
|
--y;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; Subtract significands
|
|
|
|
|
; y - x replaces y
|
|
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
esubm (x, y)
|
|
|
|
|
unsigned EMUSHORT *x, *y;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMULONG a;
|
|
|
|
|
int i;
|
|
|
|
|
unsigned int carry;
|
|
|
|
|
|
|
|
|
|
x += NI - 1;
|
|
|
|
|
y += NI - 1;
|
|
|
|
|
carry = 0;
|
|
|
|
|
for (i = M; i < NI; i++)
|
|
|
|
|
{
|
|
|
|
|
a = (unsigned EMULONG) (*y) - (unsigned EMULONG) (*x) - carry;
|
|
|
|
|
if (a & 0x10000)
|
|
|
|
|
carry = 1;
|
|
|
|
|
else
|
|
|
|
|
carry = 0;
|
|
|
|
|
*y = (unsigned EMUSHORT) a;
|
|
|
|
|
--x;
|
|
|
|
|
--y;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* Divide significands */
|
|
|
|
|
|
|
|
|
|
static unsigned EMUSHORT equot[NI];
|
|
|
|
|
|
|
|
|
|
int
|
|
|
|
|
edivm (den, num)
|
|
|
|
|
unsigned EMUSHORT den[], num[];
|
|
|
|
|
{
|
|
|
|
|
int i;
|
|
|
|
|
register unsigned EMUSHORT *p, *q;
|
|
|
|
|
unsigned EMUSHORT j;
|
|
|
|
|
|
|
|
|
|
p = &equot[0];
|
|
|
|
|
*p++ = num[0];
|
|
|
|
|
*p++ = num[1];
|
|
|
|
|
|
|
|
|
|
for (i = M; i < NI; i++)
|
|
|
|
|
{
|
|
|
|
|
*p++ = 0;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/* Use faster compare and subtraction if denominator
|
|
|
|
|
* has only 15 bits of significance.
|
|
|
|
|
*/
|
|
|
|
|
p = &den[M + 2];
|
|
|
|
|
if (*p++ == 0)
|
|
|
|
|
{
|
|
|
|
|
for (i = M + 3; i < NI; i++)
|
|
|
|
|
{
|
|
|
|
|
if (*p++ != 0)
|
|
|
|
|
goto fulldiv;
|
|
|
|
|
}
|
|
|
|
|
if ((den[M + 1] & 1) != 0)
|
|
|
|
|
goto fulldiv;
|
|
|
|
|
eshdn1 (num);
|
|
|
|
|
eshdn1 (den);
|
|
|
|
|
|
|
|
|
|
p = &den[M + 1];
|
|
|
|
|
q = &num[M + 1];
|
|
|
|
|
|
|
|
|
|
for (i = 0; i < NBITS + 2; i++)
|
|
|
|
|
{
|
|
|
|
|
if (*p <= *q)
|
|
|
|
|
{
|
|
|
|
|
*q -= *p;
|
|
|
|
|
j = 1;
|
|
|
|
|
}
|
|
|
|
|
else
|
|
|
|
|
{
|
|
|
|
|
j = 0;
|
|
|
|
|
}
|
|
|
|
|
eshup1 (equot);
|
|
|
|
|
equot[NI - 2] |= j;
|
|
|
|
|
eshup1 (num);
|
|
|
|
|
}
|
|
|
|
|
goto divdon;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/* The number of quotient bits to calculate is
|
|
|
|
|
* NBITS + 1 scaling guard bit + 1 roundoff bit.
|
|
|
|
|
*/
|
|
|
|
|
fulldiv:
|
|
|
|
|
|
|
|
|
|
p = &equot[NI - 2];
|
|
|
|
|
for (i = 0; i < NBITS + 2; i++)
|
|
|
|
|
{
|
|
|
|
|
if (ecmpm (den, num) <= 0)
|
|
|
|
|
{
|
|
|
|
|
esubm (den, num);
|
|
|
|
|
j = 1; /* quotient bit = 1 */
|
|
|
|
|
}
|
|
|
|
|
else
|
|
|
|
|
j = 0;
|
|
|
|
|
eshup1 (equot);
|
|
|
|
|
*p |= j;
|
|
|
|
|
eshup1 (num);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
divdon:
|
|
|
|
|
|
|
|
|
|
eshdn1 (equot);
|
|
|
|
|
eshdn1 (equot);
|
|
|
|
|
|
|
|
|
|
/* test for nonzero remainder after roundoff bit */
|
|
|
|
|
p = &num[M];
|
|
|
|
|
j = 0;
|
|
|
|
|
for (i = M; i < NI; i++)
|
|
|
|
|
{
|
|
|
|
|
j |= *p++;
|
|
|
|
|
}
|
|
|
|
|
if (j)
|
|
|
|
|
j = 1;
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
for (i = 0; i < NI; i++)
|
|
|
|
|
num[i] = equot[i];
|
|
|
|
|
return ((int) j);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* Multiply significands */
|
|
|
|
|
int
|
|
|
|
|
emulm (a, b)
|
|
|
|
|
unsigned EMUSHORT a[], b[];
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT *p, *q;
|
|
|
|
|
int i, j, k;
|
|
|
|
|
|
|
|
|
|
equot[0] = b[0];
|
|
|
|
|
equot[1] = b[1];
|
|
|
|
|
for (i = M; i < NI; i++)
|
|
|
|
|
equot[i] = 0;
|
|
|
|
|
|
|
|
|
|
p = &a[NI - 2];
|
|
|
|
|
k = NBITS;
|
|
|
|
|
while (*p == 0) /* significand is not supposed to be all zero */
|
|
|
|
|
{
|
|
|
|
|
eshdn6 (a);
|
|
|
|
|
k -= 16;
|
|
|
|
|
}
|
|
|
|
|
if ((*p & 0xff) == 0)
|
|
|
|
|
{
|
|
|
|
|
eshdn8 (a);
|
|
|
|
|
k -= 8;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
q = &equot[NI - 1];
|
|
|
|
|
j = 0;
|
|
|
|
|
for (i = 0; i < k; i++)
|
|
|
|
|
{
|
|
|
|
|
if (*p & 1)
|
|
|
|
|
eaddm (b, equot);
|
|
|
|
|
/* remember if there were any nonzero bits shifted out */
|
|
|
|
|
if (*q & 1)
|
|
|
|
|
j |= 1;
|
|
|
|
|
eshdn1 (a);
|
|
|
|
|
eshdn1 (equot);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
for (i = 0; i < NI; i++)
|
|
|
|
|
b[i] = equot[i];
|
|
|
|
|
|
|
|
|
|
/* return flag for lost nonzero bits */
|
|
|
|
|
return (j);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
* 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.
|
|
|
|
|
*
|
|
|
|
|
* 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 "rcntrl" is the rounding control.
|
|
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
static int rlast = -1;
|
|
|
|
|
static int rw = 0;
|
|
|
|
|
static unsigned EMUSHORT rmsk = 0;
|
|
|
|
|
static unsigned EMUSHORT rmbit = 0;
|
|
|
|
|
static unsigned EMUSHORT rebit = 0;
|
|
|
|
|
static int re = 0;
|
|
|
|
|
static unsigned EMUSHORT rbit[NI];
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
emdnorm (s, lost, subflg, exp, rcntrl)
|
|
|
|
|
unsigned EMUSHORT s[];
|
|
|
|
|
int lost;
|
|
|
|
|
int subflg;
|
|
|
|
|
EMULONG exp;
|
|
|
|
|
int rcntrl;
|
|
|
|
|
{
|
|
|
|
|
int i, j;
|
|
|
|
|
unsigned EMUSHORT r;
|
|
|
|
|
|
|
|
|
|
/* Normalize */
|
|
|
|
|
j = enormlz (s);
|
|
|
|
|
|
|
|
|
|
/* a blank significand could mean either zero or infinity. */
|
|
|
|
|
#ifndef INFINITY
|
|
|
|
|
if (j > NBITS)
|
|
|
|
|
{
|
|
|
|
|
ecleazs (s);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
exp -= j;
|
|
|
|
|
#ifndef INFINITY
|
|
|
|
|
if (exp >= 32767L)
|
|
|
|
|
goto overf;
|
|
|
|
|
#else
|
|
|
|
|
if ((j > NBITS) && (exp < 32767))
|
|
|
|
|
{
|
|
|
|
|
ecleazs (s);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
if (exp < 0L)
|
|
|
|
|
{
|
|
|
|
|
if (exp > (EMULONG) (-NBITS - 1))
|
|
|
|
|
{
|
|
|
|
|
j = (int) exp;
|
|
|
|
|
i = eshift (s, j);
|
|
|
|
|
if (i)
|
|
|
|
|
lost = 1;
|
|
|
|
|
}
|
|
|
|
|
else
|
|
|
|
|
{
|
|
|
|
|
ecleazs (s);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
/* Round off, unless told not to by rcntrl. */
|
|
|
|
|
if (rcntrl == 0)
|
|
|
|
|
goto mdfin;
|
|
|
|
|
/* Set up rounding parameters if the control register changed. */
|
|
|
|
|
if (rndprc != rlast)
|
|
|
|
|
{
|
|
|
|
|
ecleaz (rbit);
|
|
|
|
|
switch (rndprc)
|
|
|
|
|
{
|
|
|
|
|
default:
|
|
|
|
|
case NBITS:
|
|
|
|
|
rw = NI - 1; /* low guard word */
|
|
|
|
|
rmsk = 0xffff;
|
|
|
|
|
rmbit = 0x8000;
|
|
|
|
|
rbit[rw - 1] = 1;
|
|
|
|
|
re = NI - 2;
|
|
|
|
|
rebit = 1;
|
|
|
|
|
break;
|
|
|
|
|
case 64:
|
|
|
|
|
rw = 7;
|
|
|
|
|
rmsk = 0xffff;
|
|
|
|
|
rmbit = 0x8000;
|
|
|
|
|
rbit[rw - 1] = 1;
|
|
|
|
|
re = rw - 1;
|
|
|
|
|
rebit = 1;
|
|
|
|
|
break;
|
|
|
|
|
/* For DEC arithmetic */
|
|
|
|
|
case 56:
|
|
|
|
|
rw = 6;
|
|
|
|
|
rmsk = 0xff;
|
|
|
|
|
rmbit = 0x80;
|
|
|
|
|
rbit[rw] = 0x100;
|
|
|
|
|
re = rw;
|
|
|
|
|
rebit = 0x100;
|
|
|
|
|
break;
|
|
|
|
|
case 53:
|
|
|
|
|
rw = 6;
|
|
|
|
|
rmsk = 0x7ff;
|
|
|
|
|
rmbit = 0x0400;
|
|
|
|
|
rbit[rw] = 0x800;
|
|
|
|
|
re = rw;
|
|
|
|
|
rebit = 0x800;
|
|
|
|
|
break;
|
|
|
|
|
case 24:
|
|
|
|
|
rw = 4;
|
|
|
|
|
rmsk = 0xff;
|
|
|
|
|
rmbit = 0x80;
|
|
|
|
|
rbit[rw] = 0x100;
|
|
|
|
|
re = rw;
|
|
|
|
|
rebit = 0x100;
|
|
|
|
|
break;
|
|
|
|
|
}
|
|
|
|
|
rlast = rndprc;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
if (rndprc >= 64)
|
|
|
|
|
{
|
|
|
|
|
r = s[rw] & rmsk;
|
|
|
|
|
if (rndprc == 64)
|
|
|
|
|
{
|
|
|
|
|
i = rw + 1;
|
|
|
|
|
while (i < NI)
|
|
|
|
|
{
|
|
|
|
|
if (s[i])
|
|
|
|
|
r |= 1;
|
|
|
|
|
s[i] = 0;
|
|
|
|
|
++i;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
else
|
|
|
|
|
{
|
|
|
|
|
if (exp <= 0)
|
|
|
|
|
eshdn1 (s);
|
|
|
|
|
r = s[rw] & rmsk;
|
|
|
|
|
/* These tests assume NI = 8 */
|
|
|
|
|
i = rw + 1;
|
|
|
|
|
while (i < NI)
|
|
|
|
|
{
|
|
|
|
|
if (s[i])
|
|
|
|
|
r |= 1;
|
|
|
|
|
s[i] = 0;
|
|
|
|
|
++i;
|
|
|
|
|
}
|
|
|
|
|
/*
|
|
|
|
|
if (rndprc == 24)
|
|
|
|
|
{
|
|
|
|
|
if (s[5] || s[6])
|
|
|
|
|
r |= 1;
|
|
|
|
|
s[5] = 0;
|
|
|
|
|
s[6] = 0;
|
|
|
|
|
}
|
|
|
|
|
*/
|
|
|
|
|
}
|
|
|
|
|
s[rw] &= ~rmsk;
|
|
|
|
|
if ((r & rmbit) != 0)
|
|
|
|
|
{
|
|
|
|
|
if (r == rmbit)
|
|
|
|
|
{
|
|
|
|
|
if (lost == 0)
|
|
|
|
|
{ /* round to even */
|
|
|
|
|
if ((s[re] & rebit) == 0)
|
|
|
|
|
goto mddone;
|
|
|
|
|
}
|
|
|
|
|
else
|
|
|
|
|
{
|
|
|
|
|
if (subflg != 0)
|
|
|
|
|
goto mddone;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
eaddm (rbit, s);
|
|
|
|
|
}
|
|
|
|
|
mddone:
|
|
|
|
|
if ((rndprc < 64) && (exp <= 0))
|
|
|
|
|
{
|
|
|
|
|
eshup1 (s);
|
|
|
|
|
}
|
|
|
|
|
if (s[2] != 0)
|
|
|
|
|
{ /* overflow on roundoff */
|
|
|
|
|
eshdn1 (s);
|
|
|
|
|
exp += 1;
|
|
|
|
|
}
|
|
|
|
|
mdfin:
|
|
|
|
|
s[NI - 1] = 0;
|
|
|
|
|
if (exp >= 32767L)
|
|
|
|
|
{
|
|
|
|
|
#ifndef INFINITY
|
|
|
|
|
overf:
|
|
|
|
|
#endif
|
|
|
|
|
#ifdef INFINITY
|
|
|
|
|
s[1] = 32767;
|
|
|
|
|
for (i = 2; i < NI - 1; i++)
|
|
|
|
|
s[i] = 0;
|
|
|
|
|
if (extra_warnings)
|
|
|
|
|
warning ("floating point overflow");
|
|
|
|
|
#else
|
|
|
|
|
s[1] = 32766;
|
|
|
|
|
s[2] = 0;
|
|
|
|
|
for (i = M + 1; i < NI - 1; i++)
|
|
|
|
|
s[i] = 0xffff;
|
|
|
|
|
s[NI - 1] = 0;
|
|
|
|
|
if (rndprc < 64)
|
|
|
|
|
{
|
|
|
|
|
s[rw] &= ~rmsk;
|
|
|
|
|
if (rndprc == 24)
|
|
|
|
|
{
|
|
|
|
|
s[5] = 0;
|
|
|
|
|
s[6] = 0;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
if (exp < 0)
|
|
|
|
|
s[1] = 0;
|
|
|
|
|
else
|
|
|
|
|
s[1] = (unsigned EMUSHORT) exp;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; Subtract external format numbers.
|
|
|
|
|
;
|
|
|
|
|
; unsigned EMUSHORT a[NE], b[NE], c[NE];
|
|
|
|
|
; esub (a, b, c); c = b - a
|
|
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
static int subflg = 0;
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
esub (a, b, c)
|
|
|
|
|
unsigned EMUSHORT *a, *b, *c;
|
|
|
|
|
{
|
|
|
|
|
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
if (eisnan (a))
|
|
|
|
|
{
|
|
|
|
|
emov (a, c);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
if (eisnan (b))
|
|
|
|
|
{
|
|
|
|
|
emov (b, c);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
/* Infinity minus infinity is a NaN.
|
|
|
|
|
Test for subtracting infinities of the same sign. */
|
|
|
|
|
if (eisinf (a) && eisinf (b)
|
|
|
|
|
&& ((eisneg (a) ^ eisneg (b)) == 0))
|
|
|
|
|
{
|
|
|
|
|
mtherr ("esub", INVALID);
|
|
|
|
|
enan (c);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
subflg = 1;
|
|
|
|
|
eadd1 (a, b, c);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; Add.
|
|
|
|
|
;
|
|
|
|
|
; unsigned EMUSHORT a[NE], b[NE], c[NE];
|
|
|
|
|
; eadd (a, b, c); c = b + a
|
|
|
|
|
*/
|
|
|
|
|
void
|
|
|
|
|
eadd (a, b, c)
|
|
|
|
|
unsigned EMUSHORT *a, *b, *c;
|
|
|
|
|
{
|
|
|
|
|
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
/* NaN plus anything is a NaN. */
|
|
|
|
|
if (eisnan (a))
|
|
|
|
|
{
|
|
|
|
|
emov (a, c);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
if (eisnan (b))
|
|
|
|
|
{
|
|
|
|
|
emov (b, c);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
/* Infinity minus infinity is a NaN.
|
|
|
|
|
Test for adding infinities of opposite signs. */
|
|
|
|
|
if (eisinf (a) && eisinf (b)
|
|
|
|
|
&& ((eisneg (a) ^ eisneg (b)) != 0))
|
|
|
|
|
{
|
|
|
|
|
mtherr ("esub", INVALID);
|
|
|
|
|
enan (c);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
subflg = 0;
|
|
|
|
|
eadd1 (a, b, c);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
eadd1 (a, b, c)
|
|
|
|
|
unsigned EMUSHORT *a, *b, *c;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT ai[NI], bi[NI], ci[NI];
|
|
|
|
|
int i, lost, j, k;
|
|
|
|
|
EMULONG lt, lta, ltb;
|
|
|
|
|
|
|
|
|
|
#ifdef INFINITY
|
|
|
|
|
if (eisinf (a))
|
|
|
|
|
{
|
|
|
|
|
emov (a, c);
|
|
|
|
|
if (subflg)
|
|
|
|
|
eneg (c);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
if (eisinf (b))
|
|
|
|
|
{
|
|
|
|
|
emov (b, c);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
emovi (a, ai);
|
|
|
|
|
emovi (b, bi);
|
|
|
|
|
if (subflg)
|
|
|
|
|
ai[0] = ~ai[0];
|
|
|
|
|
|
|
|
|
|
/* compare exponents */
|
|
|
|
|
lta = ai[E];
|
|
|
|
|
ltb = bi[E];
|
|
|
|
|
lt = lta - ltb;
|
|
|
|
|
if (lt > 0L)
|
|
|
|
|
{ /* put the larger number in bi */
|
|
|
|
|
emovz (bi, ci);
|
|
|
|
|
emovz (ai, bi);
|
|
|
|
|
emovz (ci, ai);
|
|
|
|
|
ltb = bi[E];
|
|
|
|
|
lt = -lt;
|
|
|
|
|
}
|
|
|
|
|
lost = 0;
|
|
|
|
|
if (lt != 0L)
|
|
|
|
|
{
|
|
|
|
|
if (lt < (EMULONG) (-NBITS - 1))
|
|
|
|
|
goto done; /* answer same as larger addend */
|
|
|
|
|
k = (int) lt;
|
|
|
|
|
lost = eshift (ai, k); /* shift the smaller number down */
|
|
|
|
|
}
|
|
|
|
|
else
|
|
|
|
|
{
|
|
|
|
|
/* exponents were the same, so must compare significands */
|
|
|
|
|
i = ecmpm (ai, bi);
|
|
|
|
|
if (i == 0)
|
|
|
|
|
{ /* the numbers are identical in magnitude */
|
|
|
|
|
/* if different signs, result is zero */
|
|
|
|
|
if (ai[0] != bi[0])
|
|
|
|
|
{
|
|
|
|
|
eclear (c);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
/* if same sign, result is double */
|
|
|
|
|
/* double denomalized tiny number */
|
|
|
|
|
if ((bi[E] == 0) && ((bi[3] & 0x8000) == 0))
|
|
|
|
|
{
|
|
|
|
|
eshup1 (bi);
|
|
|
|
|
goto done;
|
|
|
|
|
}
|
|
|
|
|
/* add 1 to exponent unless both are zero! */
|
|
|
|
|
for (j = 1; j < NI - 1; j++)
|
|
|
|
|
{
|
|
|
|
|
if (bi[j] != 0)
|
|
|
|
|
{
|
|
|
|
|
/* This could overflow, but let emovo take care of that. */
|
|
|
|
|
ltb += 1;
|
|
|
|
|
break;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
bi[E] = (unsigned EMUSHORT) ltb;
|
|
|
|
|
goto done;
|
|
|
|
|
}
|
|
|
|
|
if (i > 0)
|
|
|
|
|
{ /* put the larger number in bi */
|
|
|
|
|
emovz (bi, ci);
|
|
|
|
|
emovz (ai, bi);
|
|
|
|
|
emovz (ci, ai);
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
if (ai[0] == bi[0])
|
|
|
|
|
{
|
|
|
|
|
eaddm (ai, bi);
|
|
|
|
|
subflg = 0;
|
|
|
|
|
}
|
|
|
|
|
else
|
|
|
|
|
{
|
|
|
|
|
esubm (ai, bi);
|
|
|
|
|
subflg = 1;
|
|
|
|
|
}
|
|
|
|
|
emdnorm (bi, lost, subflg, ltb, 64);
|
|
|
|
|
|
|
|
|
|
done:
|
|
|
|
|
emovo (bi, c);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; Divide.
|
|
|
|
|
;
|
|
|
|
|
; unsigned EMUSHORT a[NE], b[NE], c[NE];
|
|
|
|
|
; ediv (a, b, c); c = b / a
|
|
|
|
|
*/
|
|
|
|
|
void
|
|
|
|
|
ediv (a, b, c)
|
|
|
|
|
unsigned EMUSHORT *a, *b, *c;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT ai[NI], bi[NI];
|
|
|
|
|
int i;
|
|
|
|
|
EMULONG lt, lta, ltb;
|
|
|
|
|
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
/* Return any NaN input. */
|
|
|
|
|
if (eisnan (a))
|
|
|
|
|
{
|
|
|
|
|
emov (a, c);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
if (eisnan (b))
|
|
|
|
|
{
|
|
|
|
|
emov (b, c);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
/* 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);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
/* 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;
|
|
|
|
|
}
|
|
|
|
|
/* Anything else over infinity is zero. */
|
|
|
|
|
if (eisinf (a))
|
|
|
|
|
{
|
|
|
|
|
eclear (c);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
emovi (a, ai);
|
|
|
|
|
emovi (b, bi);
|
|
|
|
|
lta = ai[E];
|
|
|
|
|
ltb = bi[E];
|
|
|
|
|
if (bi[E] == 0)
|
|
|
|
|
{ /* See if numerator is zero. */
|
|
|
|
|
for (i = 1; i < NI - 1; i++)
|
|
|
|
|
{
|
|
|
|
|
if (bi[i] != 0)
|
|
|
|
|
{
|
|
|
|
|
ltb -= enormlz (bi);
|
|
|
|
|
goto dnzro1;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
eclear (c);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
dnzro1:
|
|
|
|
|
|
|
|
|
|
if (ai[E] == 0)
|
|
|
|
|
{ /* possible divide by zero */
|
|
|
|
|
for (i = 1; i < NI - 1; i++)
|
|
|
|
|
{
|
|
|
|
|
if (ai[i] != 0)
|
|
|
|
|
{
|
|
|
|
|
lta -= enormlz (ai);
|
|
|
|
|
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;
|
|
|
|
|
}
|
|
|
|
|
dnzro2:
|
|
|
|
|
|
|
|
|
|
i = edivm (ai, bi);
|
|
|
|
|
/* 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);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; Multiply.
|
|
|
|
|
;
|
|
|
|
|
; unsigned EMUSHORT a[NE], b[NE], c[NE];
|
|
|
|
|
; emul (a, b, c); c = b * a
|
|
|
|
|
*/
|
|
|
|
|
void
|
|
|
|
|
emul (a, b, c)
|
|
|
|
|
unsigned EMUSHORT *a, *b, *c;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT ai[NI], bi[NI];
|
|
|
|
|
int i, j;
|
|
|
|
|
EMULONG lt, lta, ltb;
|
|
|
|
|
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
/* NaN times anything is the same NaN. */
|
|
|
|
|
if (eisnan (a))
|
|
|
|
|
{
|
|
|
|
|
emov (a, c);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
if (eisnan (b))
|
|
|
|
|
{
|
|
|
|
|
emov (b, c);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
/* Zero times infinity is a NaN. */
|
|
|
|
|
if ((eisinf (a) && (ecmp (b, ezero) == 0))
|
|
|
|
|
|| (eisinf (b) && (ecmp (a, ezero) == 0)))
|
|
|
|
|
{
|
|
|
|
|
mtherr ("emul", INVALID);
|
|
|
|
|
enan (c);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
/* 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;
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
emovi (a, ai);
|
|
|
|
|
emovi (b, bi);
|
|
|
|
|
lta = ai[E];
|
|
|
|
|
ltb = bi[E];
|
|
|
|
|
if (ai[E] == 0)
|
|
|
|
|
{
|
|
|
|
|
for (i = 1; i < NI - 1; i++)
|
|
|
|
|
{
|
|
|
|
|
if (ai[i] != 0)
|
|
|
|
|
{
|
|
|
|
|
lta -= enormlz (ai);
|
|
|
|
|
goto mnzer1;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
eclear (c);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
mnzer1:
|
|
|
|
|
|
|
|
|
|
if (bi[E] == 0)
|
|
|
|
|
{
|
|
|
|
|
for (i = 1; i < NI - 1; i++)
|
|
|
|
|
{
|
|
|
|
|
if (bi[i] != 0)
|
|
|
|
|
{
|
|
|
|
|
ltb -= enormlz (bi);
|
|
|
|
|
goto mnzer2;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
eclear (c);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
mnzer2:
|
|
|
|
|
|
|
|
|
|
/* Multiply significands */
|
|
|
|
|
j = emulm (ai, bi);
|
|
|
|
|
/* 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);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; Convert IEEE double precision to e type
|
|
|
|
|
; double d;
|
|
|
|
|
; unsigned EMUSHORT x[N+2];
|
|
|
|
|
; e53toe (&d, x);
|
|
|
|
|
*/
|
|
|
|
|
void
|
|
|
|
|
e53toe (pe, y)
|
|
|
|
|
unsigned EMUSHORT *pe, *y;
|
|
|
|
|
{
|
|
|
|
|
#ifdef DEC
|
|
|
|
|
|
|
|
|
|
dectoe (pe, y); /* see etodec.c */
|
|
|
|
|
|
|
|
|
|
#else
|
|
|
|
|
|
|
|
|
|
register unsigned EMUSHORT r;
|
|
|
|
|
register unsigned EMUSHORT *e, *p;
|
|
|
|
|
unsigned EMUSHORT yy[NI];
|
|
|
|
|
int denorm, k;
|
|
|
|
|
|
|
|
|
|
e = pe;
|
|
|
|
|
denorm = 0; /* flag if denormalized number */
|
|
|
|
|
ecleaz (yy);
|
|
|
|
|
#ifdef IBMPC
|
|
|
|
|
e += 3;
|
|
|
|
|
#endif
|
|
|
|
|
r = *e;
|
|
|
|
|
yy[0] = 0;
|
|
|
|
|
if (r & 0x8000)
|
|
|
|
|
yy[0] = 0xffff;
|
|
|
|
|
yy[M] = (r & 0x0f) | 0x10;
|
|
|
|
|
r &= ~0x800f; /* strip sign and 4 significand bits */
|
|
|
|
|
#ifdef INFINITY
|
|
|
|
|
if (r == 0x7ff0)
|
|
|
|
|
{
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
#ifdef IBMPC
|
|
|
|
|
if (((pe[3] & 0xf) != 0) || (pe[2] != 0)
|
|
|
|
|
|| (pe[1] != 0) || (pe[0] != 0))
|
|
|
|
|
{
|
|
|
|
|
enan (y);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#else
|
|
|
|
|
if (((pe[0] & 0xf) != 0) || (pe[1] != 0)
|
|
|
|
|
|| (pe[2] != 0) || (pe[3] != 0))
|
|
|
|
|
{
|
|
|
|
|
enan (y);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
#endif /* NANS */
|
|
|
|
|
eclear (y);
|
|
|
|
|
einfin (y);
|
|
|
|
|
if (yy[0])
|
|
|
|
|
eneg (y);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#endif /* INFINITY */
|
|
|
|
|
r >>= 4;
|
|
|
|
|
/* If zero exponent, then the significand is denormalized.
|
|
|
|
|
* So, take back the understood high significand bit. */
|
|
|
|
|
if (r == 0)
|
|
|
|
|
{
|
|
|
|
|
denorm = 1;
|
|
|
|
|
yy[M] &= ~0x10;
|
|
|
|
|
}
|
|
|
|
|
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++;
|
|
|
|
|
#endif
|
|
|
|
|
eshift (yy, -5);
|
|
|
|
|
if (denorm)
|
|
|
|
|
{ /* if zero exponent, then normalize the significand */
|
|
|
|
|
if ((k = enormlz (yy)) > NBITS)
|
|
|
|
|
ecleazs (yy);
|
|
|
|
|
else
|
|
|
|
|
yy[E] -= (unsigned EMUSHORT) (k - 1);
|
|
|
|
|
}
|
|
|
|
|
emovo (yy, y);
|
|
|
|
|
#endif /* not DEC */
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
e64toe (pe, y)
|
|
|
|
|
unsigned EMUSHORT *pe, *y;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT yy[NI];
|
|
|
|
|
unsigned EMUSHORT *e, *p, *q;
|
|
|
|
|
int i;
|
|
|
|
|
|
|
|
|
|
e = pe;
|
|
|
|
|
p = yy;
|
|
|
|
|
for (i = 0; i < NE - 5; i++)
|
|
|
|
|
*p++ = 0;
|
|
|
|
|
#ifdef IBMPC
|
|
|
|
|
for (i = 0; i < 5; i++)
|
|
|
|
|
*p++ = *e++;
|
|
|
|
|
#endif
|
|
|
|
|
#ifdef DEC
|
|
|
|
|
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++;
|
|
|
|
|
#endif
|
|
|
|
|
p = yy;
|
|
|
|
|
q = y;
|
|
|
|
|
#ifdef INFINITY
|
|
|
|
|
if (*p == 0x7fff)
|
|
|
|
|
{
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
#ifdef IBMPC
|
|
|
|
|
for (i = 0; i < 4; i++)
|
|
|
|
|
{
|
|
|
|
|
if (pe[i] != 0)
|
|
|
|
|
{
|
|
|
|
|
enan (y);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
#else
|
|
|
|
|
for (i = 1; i <= 4; i++)
|
|
|
|
|
{
|
|
|
|
|
if (pe[i] != 0)
|
|
|
|
|
{
|
|
|
|
|
enan (y);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
#endif /* NANS */
|
|
|
|
|
eclear (y);
|
|
|
|
|
einfin (y);
|
|
|
|
|
if (*p & 0x8000)
|
|
|
|
|
eneg (y);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#endif /* INFINITY */
|
|
|
|
|
for (i = 0; i < NE; i++)
|
|
|
|
|
*q++ = *p++;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; Convert IEEE single precision to e type
|
|
|
|
|
; float d;
|
|
|
|
|
; unsigned EMUSHORT x[N+2];
|
|
|
|
|
; dtox (&d, x);
|
|
|
|
|
*/
|
|
|
|
|
void
|
|
|
|
|
e24toe (pe, y)
|
|
|
|
|
unsigned EMUSHORT *pe, *y;
|
|
|
|
|
{
|
|
|
|
|
register unsigned EMUSHORT r;
|
|
|
|
|
register unsigned EMUSHORT *e, *p;
|
|
|
|
|
unsigned EMUSHORT yy[NI];
|
|
|
|
|
int denorm, k;
|
|
|
|
|
|
|
|
|
|
e = pe;
|
|
|
|
|
denorm = 0; /* flag if denormalized number */
|
|
|
|
|
ecleaz (yy);
|
|
|
|
|
#ifdef IBMPC
|
|
|
|
|
e += 1;
|
|
|
|
|
#endif
|
|
|
|
|
#ifdef DEC
|
|
|
|
|
e += 1;
|
|
|
|
|
#endif
|
|
|
|
|
r = *e;
|
|
|
|
|
yy[0] = 0;
|
|
|
|
|
if (r & 0x8000)
|
|
|
|
|
yy[0] = 0xffff;
|
|
|
|
|
yy[M] = (r & 0x7f) | 0200;
|
|
|
|
|
r &= ~0x807f; /* strip sign and 7 significand bits */
|
|
|
|
|
#ifdef INFINITY
|
|
|
|
|
if (r == 0x7f80)
|
|
|
|
|
{
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
#ifdef MIEEE
|
|
|
|
|
if (((pe[0] & 0x7f) != 0) || (pe[1] != 0))
|
|
|
|
|
{
|
|
|
|
|
enan (y);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#else
|
|
|
|
|
if (((pe[1] & 0x7f) != 0) || (pe[0] != 0))
|
|
|
|
|
{
|
|
|
|
|
enan (y);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
#endif /* NANS */
|
|
|
|
|
eclear (y);
|
|
|
|
|
einfin (y);
|
|
|
|
|
if (yy[0])
|
|
|
|
|
eneg (y);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#endif /* INFINITY */
|
|
|
|
|
r >>= 7;
|
|
|
|
|
/* If zero exponent, then the significand is denormalized.
|
|
|
|
|
* So, take back the understood high significand bit. */
|
|
|
|
|
if (r == 0)
|
|
|
|
|
{
|
|
|
|
|
denorm = 1;
|
|
|
|
|
yy[M] &= ~0200;
|
|
|
|
|
}
|
|
|
|
|
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++;
|
|
|
|
|
#endif
|
|
|
|
|
eshift (yy, -8);
|
|
|
|
|
if (denorm)
|
|
|
|
|
{ /* if zero exponent, then normalize the significand */
|
|
|
|
|
if ((k = enormlz (yy)) > NBITS)
|
|
|
|
|
ecleazs (yy);
|
|
|
|
|
else
|
|
|
|
|
yy[E] -= (unsigned EMUSHORT) (k - 1);
|
|
|
|
|
}
|
|
|
|
|
emovo (yy, y);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
etoe64 (x, e)
|
|
|
|
|
unsigned EMUSHORT *x, *e;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT xi[NI];
|
|
|
|
|
EMULONG exp;
|
|
|
|
|
int rndsav;
|
|
|
|
|
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
if (eisnan (x))
|
|
|
|
|
{
|
|
|
|
|
make_nan (e, XFmode);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
emovi (x, xi);
|
|
|
|
|
/* adjust exponent for offset */
|
|
|
|
|
exp = (EMULONG) xi[E];
|
|
|
|
|
#ifdef INFINITY
|
|
|
|
|
if (eisinf (x))
|
|
|
|
|
goto nonorm;
|
|
|
|
|
#endif
|
|
|
|
|
/* round off to nearest or even */
|
|
|
|
|
rndsav = rndprc;
|
|
|
|
|
rndprc = 64;
|
|
|
|
|
emdnorm (xi, 0, 0, exp, 64);
|
|
|
|
|
rndprc = rndsav;
|
|
|
|
|
nonorm:
|
|
|
|
|
toe64 (xi, e);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/* move out internal format to ieee long double */
|
|
|
|
|
static void
|
|
|
|
|
toe64 (a, b)
|
|
|
|
|
unsigned EMUSHORT *a, *b;
|
|
|
|
|
{
|
|
|
|
|
register unsigned EMUSHORT *p, *q;
|
|
|
|
|
unsigned EMUSHORT i;
|
|
|
|
|
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
if (eiisnan (a))
|
|
|
|
|
{
|
|
|
|
|
make_nan (b, XFmode);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
p = a;
|
|
|
|
|
#ifdef MIEEE
|
|
|
|
|
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;
|
|
|
|
|
#endif
|
|
|
|
|
#endif
|
|
|
|
|
|
|
|
|
|
/* combine sign and exponent */
|
|
|
|
|
i = *p++;
|
|
|
|
|
#ifdef MIEEE
|
|
|
|
|
if (i)
|
|
|
|
|
*q++ = *p++ | 0x8000;
|
|
|
|
|
else
|
|
|
|
|
*q++ = *p++;
|
|
|
|
|
*q++ = 0;
|
|
|
|
|
#else
|
|
|
|
|
if (i)
|
|
|
|
|
*q-- = *p++ | 0x8000;
|
|
|
|
|
else
|
|
|
|
|
*q-- = *p++;
|
|
|
|
|
#endif
|
|
|
|
|
/* skip over guard word */
|
|
|
|
|
++p;
|
|
|
|
|
/* move the significand */
|
|
|
|
|
#ifdef MIEEE
|
|
|
|
|
for (i = 0; i < 4; i++)
|
|
|
|
|
*q++ = *p++;
|
|
|
|
|
#else
|
|
|
|
|
for (i = 0; i < 4; i++)
|
|
|
|
|
*q-- = *p++;
|
|
|
|
|
#endif
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; e type to IEEE double precision
|
|
|
|
|
; double d;
|
|
|
|
|
; unsigned EMUSHORT x[NE];
|
|
|
|
|
; etoe53 (x, &d);
|
|
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
#ifdef DEC
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
etoe53 (x, e)
|
|
|
|
|
unsigned EMUSHORT *x, *e;
|
|
|
|
|
{
|
|
|
|
|
etodec (x, e); /* see etodec.c */
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
static void
|
|
|
|
|
toe53 (x, y)
|
|
|
|
|
unsigned EMUSHORT *x, *y;
|
|
|
|
|
{
|
|
|
|
|
todec (x, y);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
#else
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
etoe53 (x, e)
|
|
|
|
|
unsigned EMUSHORT *x, *e;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT xi[NI];
|
|
|
|
|
EMULONG exp;
|
|
|
|
|
int rndsav;
|
|
|
|
|
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
if (eisnan (x))
|
|
|
|
|
{
|
|
|
|
|
make_nan (e, DFmode);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
emovi (x, xi);
|
|
|
|
|
/* adjust exponent for offsets */
|
|
|
|
|
exp = (EMULONG) xi[E] - (EXONE - 0x3ff);
|
|
|
|
|
#ifdef INFINITY
|
|
|
|
|
if (eisinf (x))
|
|
|
|
|
goto nonorm;
|
|
|
|
|
#endif
|
|
|
|
|
/* round off to nearest or even */
|
|
|
|
|
rndsav = rndprc;
|
|
|
|
|
rndprc = 53;
|
|
|
|
|
emdnorm (xi, 0, 0, exp, 64);
|
|
|
|
|
rndprc = rndsav;
|
|
|
|
|
nonorm:
|
|
|
|
|
toe53 (xi, e);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
static void
|
|
|
|
|
toe53 (x, y)
|
|
|
|
|
unsigned EMUSHORT *x, *y;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT i;
|
|
|
|
|
unsigned EMUSHORT *p;
|
|
|
|
|
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
if (eiisnan (x))
|
|
|
|
|
{
|
|
|
|
|
make_nan (y, DFmode);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
p = &x[0];
|
|
|
|
|
#ifdef IBMPC
|
|
|
|
|
y += 3;
|
|
|
|
|
#endif
|
|
|
|
|
*y = 0; /* output high order */
|
|
|
|
|
if (*p++)
|
|
|
|
|
*y = 0x8000; /* output sign bit */
|
|
|
|
|
|
|
|
|
|
i = *p++;
|
|
|
|
|
if (i >= (unsigned int) 2047)
|
|
|
|
|
{ /* 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
|
|
|
|
|
#else
|
|
|
|
|
*y |= (unsigned EMUSHORT) 0x7fef;
|
|
|
|
|
#ifdef IBMPC
|
|
|
|
|
*(--y) = 0xffff;
|
|
|
|
|
*(--y) = 0xffff;
|
|
|
|
|
*(--y) = 0xffff;
|
|
|
|
|
#endif
|
|
|
|
|
#ifdef MIEEE
|
|
|
|
|
++y;
|
|
|
|
|
*y++ = 0xffff;
|
|
|
|
|
*y++ = 0xffff;
|
|
|
|
|
*y++ = 0xffff;
|
|
|
|
|
#endif
|
|
|
|
|
#endif
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
if (i == 0)
|
|
|
|
|
{
|
|
|
|
|
eshift (x, 4);
|
|
|
|
|
}
|
|
|
|
|
else
|
|
|
|
|
{
|
|
|
|
|
i <<= 4;
|
|
|
|
|
eshift (x, 5);
|
|
|
|
|
}
|
|
|
|
|
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
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
#endif /* not DEC */
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; e type to IEEE single precision
|
|
|
|
|
; float d;
|
|
|
|
|
; unsigned EMUSHORT x[N+2];
|
|
|
|
|
; xtod (x, &d);
|
|
|
|
|
*/
|
|
|
|
|
void
|
|
|
|
|
etoe24 (x, e)
|
|
|
|
|
unsigned EMUSHORT *x, *e;
|
|
|
|
|
{
|
|
|
|
|
EMULONG exp;
|
|
|
|
|
unsigned EMUSHORT xi[NI];
|
|
|
|
|
int rndsav;
|
|
|
|
|
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
if (eisnan (x))
|
|
|
|
|
{
|
|
|
|
|
make_nan (e, SFmode);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
emovi (x, xi);
|
|
|
|
|
/* adjust exponent for offsets */
|
|
|
|
|
exp = (EMULONG) xi[E] - (EXONE - 0177);
|
|
|
|
|
#ifdef INFINITY
|
|
|
|
|
if (eisinf (x))
|
|
|
|
|
goto nonorm;
|
|
|
|
|
#endif
|
|
|
|
|
/* round off to nearest or even */
|
|
|
|
|
rndsav = rndprc;
|
|
|
|
|
rndprc = 24;
|
|
|
|
|
emdnorm (xi, 0, 0, exp, 64);
|
|
|
|
|
rndprc = rndsav;
|
|
|
|
|
nonorm:
|
|
|
|
|
toe24 (xi, e);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
static void
|
|
|
|
|
toe24 (x, y)
|
|
|
|
|
unsigned EMUSHORT *x, *y;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT i;
|
|
|
|
|
unsigned EMUSHORT *p;
|
|
|
|
|
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
if (eiisnan (x))
|
|
|
|
|
{
|
|
|
|
|
make_nan (y, SFmode);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
p = &x[0];
|
|
|
|
|
#ifdef IBMPC
|
|
|
|
|
y += 1;
|
|
|
|
|
#endif
|
|
|
|
|
#ifdef DEC
|
|
|
|
|
y += 1;
|
|
|
|
|
#endif
|
|
|
|
|
*y = 0; /* output high order */
|
|
|
|
|
if (*p++)
|
|
|
|
|
*y = 0x8000; /* output sign bit */
|
|
|
|
|
|
|
|
|
|
i = *p++;
|
|
|
|
|
/* 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;
|
|
|
|
|
#endif
|
|
|
|
|
#else /* no INFINITY */
|
|
|
|
|
*y |= (unsigned EMUSHORT) 0x7f7f;
|
|
|
|
|
#ifdef IBMPC
|
|
|
|
|
*(--y) = 0xffff;
|
|
|
|
|
#endif
|
|
|
|
|
#ifdef DEC
|
|
|
|
|
*(--y) = 0xffff;
|
|
|
|
|
#endif
|
|
|
|
|
#ifdef MIEEE
|
|
|
|
|
++y;
|
|
|
|
|
*y = 0xffff;
|
|
|
|
|
#endif
|
|
|
|
|
#ifdef ERANGE
|
|
|
|
|
errno = ERANGE;
|
|
|
|
|
#endif
|
|
|
|
|
#endif /* no INFINITY */
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
if (i == 0)
|
|
|
|
|
{
|
|
|
|
|
eshift (x, 7);
|
|
|
|
|
}
|
|
|
|
|
else
|
|
|
|
|
{
|
|
|
|
|
i <<= 7;
|
|
|
|
|
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
|
|
|
|
|
#ifdef DEC
|
|
|
|
|
*(--y) = *p;
|
|
|
|
|
#endif
|
|
|
|
|
#ifdef MIEEE
|
|
|
|
|
++y;
|
|
|
|
|
*y = *p;
|
|
|
|
|
#endif
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* Compare two e type numbers.
|
|
|
|
|
*
|
|
|
|
|
* unsigned EMUSHORT a[NE], b[NE];
|
|
|
|
|
* ecmp (a, b);
|
|
|
|
|
*
|
|
|
|
|
* returns +1 if a > b
|
|
|
|
|
* 0 if a == b
|
|
|
|
|
* -1 if a < b
|
|
|
|
|
* -2 if either a or b is a NaN.
|
|
|
|
|
*/
|
|
|
|
|
int
|
|
|
|
|
ecmp (a, b)
|
|
|
|
|
unsigned EMUSHORT *a, *b;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT ai[NI], bi[NI];
|
|
|
|
|
register unsigned EMUSHORT *p, *q;
|
|
|
|
|
register int i;
|
|
|
|
|
int msign;
|
|
|
|
|
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
if (eisnan (a) || eisnan (b))
|
|
|
|
|
return (-2);
|
|
|
|
|
#endif
|
|
|
|
|
emovi (a, ai);
|
|
|
|
|
p = ai;
|
|
|
|
|
emovi (b, bi);
|
|
|
|
|
q = bi;
|
|
|
|
|
|
|
|
|
|
if (*p != *q)
|
|
|
|
|
{ /* the signs are different */
|
|
|
|
|
/* -0 equals + 0 */
|
|
|
|
|
for (i = 1; i < NI - 1; i++)
|
|
|
|
|
{
|
|
|
|
|
if (ai[i] != 0)
|
|
|
|
|
goto nzro;
|
|
|
|
|
if (bi[i] != 0)
|
|
|
|
|
goto nzro;
|
|
|
|
|
}
|
|
|
|
|
return (0);
|
|
|
|
|
nzro:
|
|
|
|
|
if (*p == 0)
|
|
|
|
|
return (1);
|
|
|
|
|
else
|
|
|
|
|
return (-1);
|
|
|
|
|
}
|
|
|
|
|
/* both are the same sign */
|
|
|
|
|
if (*p == 0)
|
|
|
|
|
msign = 1;
|
|
|
|
|
else
|
|
|
|
|
msign = -1;
|
|
|
|
|
i = NI - 1;
|
|
|
|
|
do
|
|
|
|
|
{
|
|
|
|
|
if (*p++ != *q++)
|
|
|
|
|
{
|
|
|
|
|
goto diff;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
while (--i > 0);
|
|
|
|
|
|
|
|
|
|
return (0); /* equality */
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
diff:
|
|
|
|
|
|
|
|
|
|
if (*(--p) > *(--q))
|
|
|
|
|
return (msign); /* p is bigger */
|
|
|
|
|
else
|
|
|
|
|
return (-msign); /* p is littler */
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* Find nearest integer to x = floor (x + 0.5)
|
|
|
|
|
*
|
|
|
|
|
* unsigned EMUSHORT x[NE], y[NE]
|
|
|
|
|
* eround (x, y);
|
|
|
|
|
*/
|
|
|
|
|
void
|
|
|
|
|
eround (x, y)
|
|
|
|
|
unsigned EMUSHORT *x, *y;
|
|
|
|
|
{
|
|
|
|
|
eadd (ehalf, x, y);
|
|
|
|
|
efloor (y, y);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; convert long integer to e type
|
|
|
|
|
;
|
|
|
|
|
; long l;
|
|
|
|
|
; unsigned EMUSHORT x[NE];
|
|
|
|
|
; ltoe (&l, x);
|
|
|
|
|
; note &l is the memory address of l
|
|
|
|
|
*/
|
|
|
|
|
void
|
|
|
|
|
ltoe (lp, y)
|
|
|
|
|
long *lp; /* lp is the memory address of a long integer */
|
|
|
|
|
unsigned EMUSHORT *y; /* y is the address of a short */
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT yi[NI];
|
|
|
|
|
unsigned long ll;
|
|
|
|
|
int k;
|
|
|
|
|
|
|
|
|
|
ecleaz (yi);
|
|
|
|
|
if (*lp < 0)
|
|
|
|
|
{
|
|
|
|
|
/* make it positive */
|
|
|
|
|
ll = (unsigned long) (-(*lp));
|
|
|
|
|
yi[0] = 0xffff; /* put correct sign in the e type number */
|
|
|
|
|
}
|
|
|
|
|
else
|
|
|
|
|
{
|
|
|
|
|
ll = (unsigned long) (*lp);
|
|
|
|
|
}
|
|
|
|
|
/* move the long integer to yi significand area */
|
|
|
|
|
#if HOST_BITS_PER_LONG == 64
|
|
|
|
|
yi[M] = (unsigned EMUSHORT) (ll >> 48);
|
|
|
|
|
yi[M + 1] = (unsigned EMUSHORT) (ll >> 32);
|
|
|
|
|
yi[M + 2] = (unsigned EMUSHORT) (ll >> 16);
|
|
|
|
|
yi[M + 3] = (unsigned EMUSHORT) ll;
|
|
|
|
|
yi[E] = EXONE + 47; /* exponent if normalize shift count were 0 */
|
|
|
|
|
#else
|
|
|
|
|
yi[M] = (unsigned EMUSHORT) (ll >> 16);
|
|
|
|
|
yi[M + 1] = (unsigned EMUSHORT) ll;
|
|
|
|
|
yi[E] = EXONE + 15; /* exponent if normalize shift count were 0 */
|
|
|
|
|
#endif
|
|
|
|
|
|
|
|
|
|
if ((k = enormlz (yi)) > NBITS)/* normalize the significand */
|
|
|
|
|
ecleaz (yi); /* it was zero */
|
|
|
|
|
else
|
|
|
|
|
yi[E] -= (unsigned EMUSHORT) k;/* subtract shift count from exponent */
|
|
|
|
|
emovo (yi, y); /* output the answer */
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; convert unsigned long integer to e type
|
|
|
|
|
;
|
|
|
|
|
; unsigned long l;
|
|
|
|
|
; unsigned EMUSHORT x[NE];
|
|
|
|
|
; ltox (&l, x);
|
|
|
|
|
; note &l is the memory address of l
|
|
|
|
|
*/
|
|
|
|
|
void
|
|
|
|
|
ultoe (lp, y)
|
|
|
|
|
unsigned long *lp; /* lp is the memory address of a long integer */
|
|
|
|
|
unsigned EMUSHORT *y; /* y is the address of a short */
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT yi[NI];
|
|
|
|
|
unsigned long ll;
|
|
|
|
|
int k;
|
|
|
|
|
|
|
|
|
|
ecleaz (yi);
|
|
|
|
|
ll = *lp;
|
|
|
|
|
|
|
|
|
|
/* move the long integer to ayi significand area */
|
|
|
|
|
#if HOST_BITS_PER_LONG == 64
|
|
|
|
|
yi[M] = (unsigned EMUSHORT) (ll >> 48);
|
|
|
|
|
yi[M + 1] = (unsigned EMUSHORT) (ll >> 32);
|
|
|
|
|
yi[M + 2] = (unsigned EMUSHORT) (ll >> 16);
|
|
|
|
|
yi[M + 3] = (unsigned EMUSHORT) ll;
|
|
|
|
|
yi[E] = EXONE + 47; /* exponent if normalize shift count were 0 */
|
|
|
|
|
#else
|
|
|
|
|
yi[M] = (unsigned EMUSHORT) (ll >> 16);
|
|
|
|
|
yi[M + 1] = (unsigned EMUSHORT) ll;
|
|
|
|
|
yi[E] = EXONE + 15; /* exponent if normalize shift count were 0 */
|
|
|
|
|
#endif
|
|
|
|
|
|
|
|
|
|
if ((k = enormlz (yi)) > NBITS)/* normalize the significand */
|
|
|
|
|
ecleaz (yi); /* it was zero */
|
|
|
|
|
else
|
|
|
|
|
yi[E] -= (unsigned EMUSHORT) k; /* subtract shift count from exponent */
|
|
|
|
|
emovo (yi, y); /* output the answer */
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; Find long integer and fractional parts
|
|
|
|
|
|
|
|
|
|
; long i;
|
|
|
|
|
; unsigned EMUSHORT x[NE], frac[NE];
|
|
|
|
|
; xifrac (x, &i, frac);
|
|
|
|
|
|
|
|
|
|
The integer output has the sign of the input. The fraction is
|
|
|
|
|
the positive fractional part of abs (x).
|
|
|
|
|
*/
|
|
|
|
|
void
|
|
|
|
|
eifrac (x, i, frac)
|
|
|
|
|
unsigned EMUSHORT *x;
|
|
|
|
|
long *i;
|
|
|
|
|
unsigned EMUSHORT *frac;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT xi[NI];
|
|
|
|
|
int j, k;
|
|
|
|
|
unsigned long ll;
|
|
|
|
|
|
|
|
|
|
emovi (x, xi);
|
|
|
|
|
k = (int) xi[E] - (EXONE - 1);
|
|
|
|
|
if (k <= 0)
|
|
|
|
|
{
|
|
|
|
|
/* if exponent <= 0, integer = 0 and real output is fraction */
|
|
|
|
|
*i = 0L;
|
|
|
|
|
emovo (xi, frac);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
if (k > (HOST_BITS_PER_LONG - 1))
|
|
|
|
|
{
|
|
|
|
|
/* long integer overflow: output large integer
|
|
|
|
|
and correct fraction */
|
|
|
|
|
if (xi[0])
|
|
|
|
|
*i = ((unsigned long) 1) << (HOST_BITS_PER_LONG - 1);
|
|
|
|
|
else
|
|
|
|
|
*i = (((unsigned long) 1) << (HOST_BITS_PER_LONG - 1)) - 1;
|
|
|
|
|
eshift (xi, k);
|
|
|
|
|
if (extra_warnings)
|
|
|
|
|
warning ("overflow on truncation to integer");
|
|
|
|
|
}
|
|
|
|
|
else if (k > 16)
|
|
|
|
|
{
|
|
|
|
|
/* Shift more than 16 bits: first shift up k-16 mod 16,
|
|
|
|
|
then shift up by 16's. */
|
|
|
|
|
j = k - ((k >> 4) << 4);
|
|
|
|
|
eshift (xi, j);
|
|
|
|
|
ll = xi[M];
|
|
|
|
|
k -= j;
|
|
|
|
|
do
|
|
|
|
|
{
|
|
|
|
|
eshup6 (xi);
|
|
|
|
|
ll = (ll << 16) | xi[M];
|
|
|
|
|
}
|
|
|
|
|
while ((k -= 16) > 0);
|
|
|
|
|
*i = ll;
|
|
|
|
|
if (xi[0])
|
|
|
|
|
*i = -(*i);
|
|
|
|
|
}
|
|
|
|
|
else
|
|
|
|
|
{
|
|
|
|
|
/* shift not more than 16 bits */
|
|
|
|
|
eshift (xi, k);
|
|
|
|
|
*i = (long) xi[M] & 0xffff;
|
|
|
|
|
if (xi[0])
|
|
|
|
|
*i = -(*i);
|
|
|
|
|
}
|
|
|
|
|
xi[0] = 0;
|
|
|
|
|
xi[E] = EXONE - 1;
|
|
|
|
|
xi[M] = 0;
|
|
|
|
|
if ((k = enormlz (xi)) > NBITS)
|
|
|
|
|
ecleaz (xi);
|
|
|
|
|
else
|
|
|
|
|
xi[E] -= (unsigned EMUSHORT) k;
|
|
|
|
|
|
|
|
|
|
emovo (xi, frac);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* Find unsigned long integer and fractional parts.
|
|
|
|
|
A negative e type input yields integer output = 0
|
|
|
|
|
but correct fraction. */
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
euifrac (x, i, frac)
|
|
|
|
|
unsigned EMUSHORT *x;
|
|
|
|
|
unsigned long *i;
|
|
|
|
|
unsigned EMUSHORT *frac;
|
|
|
|
|
{
|
|
|
|
|
unsigned long ll;
|
|
|
|
|
unsigned EMUSHORT xi[NI];
|
|
|
|
|
int j, k;
|
|
|
|
|
|
|
|
|
|
emovi (x, xi);
|
|
|
|
|
k = (int) xi[E] - (EXONE - 1);
|
|
|
|
|
if (k <= 0)
|
|
|
|
|
{
|
|
|
|
|
/* if exponent <= 0, integer = 0 and argument is fraction */
|
|
|
|
|
*i = 0L;
|
|
|
|
|
emovo (xi, frac);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
if (k > HOST_BITS_PER_LONG)
|
|
|
|
|
{
|
|
|
|
|
/* Long integer overflow: output large integer
|
|
|
|
|
and correct fraction.
|
|
|
|
|
Note, the BSD microvax compiler says that ~(0UL)
|
|
|
|
|
is a syntax error. */
|
|
|
|
|
*i = ~(0L);
|
|
|
|
|
eshift (xi, k);
|
|
|
|
|
if (extra_warnings)
|
|
|
|
|
warning ("overflow on truncation to unsigned integer");
|
|
|
|
|
}
|
|
|
|
|
else if (k > 16)
|
|
|
|
|
{
|
|
|
|
|
/* Shift more than 16 bits: first shift up k-16 mod 16,
|
|
|
|
|
then shift up by 16's. */
|
|
|
|
|
j = k - ((k >> 4) << 4);
|
|
|
|
|
eshift (xi, j);
|
|
|
|
|
ll = xi[M];
|
|
|
|
|
k -= j;
|
|
|
|
|
do
|
|
|
|
|
{
|
|
|
|
|
eshup6 (xi);
|
|
|
|
|
ll = (ll << 16) | xi[M];
|
|
|
|
|
}
|
|
|
|
|
while ((k -= 16) > 0);
|
|
|
|
|
*i = ll;
|
|
|
|
|
}
|
|
|
|
|
else
|
|
|
|
|
{
|
|
|
|
|
/* shift not more than 16 bits */
|
|
|
|
|
eshift (xi, k);
|
|
|
|
|
*i = (long) xi[M] & 0xffff;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
if (xi[0]) /* A negative value yields unsigned integer 0. */
|
|
|
|
|
*i = 0L;
|
|
|
|
|
xi[0] = 0;
|
|
|
|
|
xi[E] = EXONE - 1;
|
|
|
|
|
xi[M] = 0;
|
|
|
|
|
if ((k = enormlz (xi)) > NBITS)
|
|
|
|
|
ecleaz (xi);
|
|
|
|
|
else
|
|
|
|
|
xi[E] -= (unsigned EMUSHORT) k;
|
|
|
|
|
|
|
|
|
|
emovo (xi, frac);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; Shift significand
|
|
|
|
|
;
|
|
|
|
|
; Shifts significand area up or down by the number of bits
|
|
|
|
|
; given by the variable sc.
|
|
|
|
|
*/
|
|
|
|
|
int
|
|
|
|
|
eshift (x, sc)
|
|
|
|
|
unsigned EMUSHORT *x;
|
|
|
|
|
int sc;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT lost;
|
|
|
|
|
unsigned EMUSHORT *p;
|
|
|
|
|
|
|
|
|
|
if (sc == 0)
|
|
|
|
|
return (0);
|
|
|
|
|
|
|
|
|
|
lost = 0;
|
|
|
|
|
p = x + NI - 1;
|
|
|
|
|
|
|
|
|
|
if (sc < 0)
|
|
|
|
|
{
|
|
|
|
|
sc = -sc;
|
|
|
|
|
while (sc >= 16)
|
|
|
|
|
{
|
|
|
|
|
lost |= *p; /* remember lost bits */
|
|
|
|
|
eshdn6 (x);
|
|
|
|
|
sc -= 16;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
while (sc >= 8)
|
|
|
|
|
{
|
|
|
|
|
lost |= *p & 0xff;
|
|
|
|
|
eshdn8 (x);
|
|
|
|
|
sc -= 8;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
while (sc > 0)
|
|
|
|
|
{
|
|
|
|
|
lost |= *p & 1;
|
|
|
|
|
eshdn1 (x);
|
|
|
|
|
sc -= 1;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
else
|
|
|
|
|
{
|
|
|
|
|
while (sc >= 16)
|
|
|
|
|
{
|
|
|
|
|
eshup6 (x);
|
|
|
|
|
sc -= 16;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
while (sc >= 8)
|
|
|
|
|
{
|
|
|
|
|
eshup8 (x);
|
|
|
|
|
sc -= 8;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
while (sc > 0)
|
|
|
|
|
{
|
|
|
|
|
eshup1 (x);
|
|
|
|
|
sc -= 1;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
if (lost)
|
|
|
|
|
lost = 1;
|
|
|
|
|
return ((int) lost);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; normalize
|
|
|
|
|
;
|
|
|
|
|
; Shift normalizes the significand area pointed to by argument
|
|
|
|
|
; shift count (up = positive) is returned.
|
|
|
|
|
*/
|
|
|
|
|
int
|
|
|
|
|
enormlz (x)
|
|
|
|
|
unsigned EMUSHORT x[];
|
|
|
|
|
{
|
|
|
|
|
register unsigned EMUSHORT *p;
|
|
|
|
|
int sc;
|
|
|
|
|
|
|
|
|
|
sc = 0;
|
|
|
|
|
p = &x[M];
|
|
|
|
|
if (*p != 0)
|
|
|
|
|
goto normdn;
|
|
|
|
|
++p;
|
|
|
|
|
if (*p & 0x8000)
|
|
|
|
|
return (0); /* already normalized */
|
|
|
|
|
while (*p == 0)
|
|
|
|
|
{
|
|
|
|
|
eshup6 (x);
|
|
|
|
|
sc += 16;
|
|
|
|
|
/* With guard word, there are NBITS+16 bits available.
|
|
|
|
|
* return true if all are zero.
|
|
|
|
|
*/
|
|
|
|
|
if (sc > NBITS)
|
|
|
|
|
return (sc);
|
|
|
|
|
}
|
|
|
|
|
/* see if high byte is zero */
|
|
|
|
|
while ((*p & 0xff00) == 0)
|
|
|
|
|
{
|
|
|
|
|
eshup8 (x);
|
|
|
|
|
sc += 8;
|
|
|
|
|
}
|
|
|
|
|
/* now shift 1 bit at a time */
|
|
|
|
|
while ((*p & 0x8000) == 0)
|
|
|
|
|
{
|
|
|
|
|
eshup1 (x);
|
|
|
|
|
sc += 1;
|
|
|
|
|
if (sc > NBITS)
|
|
|
|
|
{
|
|
|
|
|
mtherr ("enormlz", UNDERFLOW);
|
|
|
|
|
return (sc);
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
return (sc);
|
|
|
|
|
|
|
|
|
|
/* Normalize by shifting down out of the high guard word
|
|
|
|
|
of the significand */
|
|
|
|
|
normdn:
|
|
|
|
|
|
|
|
|
|
if (*p & 0xff00)
|
|
|
|
|
{
|
|
|
|
|
eshdn8 (x);
|
|
|
|
|
sc -= 8;
|
|
|
|
|
}
|
|
|
|
|
while (*p != 0)
|
|
|
|
|
{
|
|
|
|
|
eshdn1 (x);
|
|
|
|
|
sc -= 1;
|
|
|
|
|
|
|
|
|
|
if (sc < -NBITS)
|
|
|
|
|
{
|
|
|
|
|
mtherr ("enormlz", OVERFLOW);
|
|
|
|
|
return (sc);
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
return (sc);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* Convert e type number to decimal format ASCII string.
|
|
|
|
|
* The constants are for 64 bit precision.
|
|
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
#define NTEN 12
|
|
|
|
|
#define MAXP 4096
|
|
|
|
|
|
|
|
|
|
static unsigned EMUSHORT etens[NTEN + 1][NE] =
|
|
|
|
|
{
|
|
|
|
|
{0xc94c, 0x979a, 0x8a20, 0x5202, 0xc460, 0x7525,}, /* 10**4096 */
|
|
|
|
|
{0xa74d, 0x5de4, 0xc53d, 0x3b5d, 0x9e8b, 0x5a92,}, /* 10**2048 */
|
|
|
|
|
{0x650d, 0x0c17, 0x8175, 0x7586, 0xc976, 0x4d48,},
|
|
|
|
|
{0xcc65, 0x91c6, 0xa60e, 0xa0ae, 0xe319, 0x46a3,},
|
|
|
|
|
{0xddbc, 0xde8d, 0x9df9, 0xebfb, 0xaa7e, 0x4351,},
|
|
|
|
|
{0xc66f, 0x8cdf, 0x80e9, 0x47c9, 0x93ba, 0x41a8,},
|
|
|
|
|
{0x3cbf, 0xa6d5, 0xffcf, 0x1f49, 0xc278, 0x40d3,},
|
|
|
|
|
{0xf020, 0xb59d, 0x2b70, 0xada8, 0x9dc5, 0x4069,},
|
|
|
|
|
{0x0000, 0x0000, 0x0400, 0xc9bf, 0x8e1b, 0x4034,},
|
|
|
|
|
{0x0000, 0x0000, 0x0000, 0x2000, 0xbebc, 0x4019,},
|
|
|
|
|
{0x0000, 0x0000, 0x0000, 0x0000, 0x9c40, 0x400c,},
|
|
|
|
|
{0x0000, 0x0000, 0x0000, 0x0000, 0xc800, 0x4005,},
|
|
|
|
|
{0x0000, 0x0000, 0x0000, 0x0000, 0xa000, 0x4002,}, /* 10**1 */
|
|
|
|
|
};
|
|
|
|
|
|
|
|
|
|
static unsigned EMUSHORT emtens[NTEN + 1][NE] =
|
|
|
|
|
{
|
|
|
|
|
{0x2de4, 0x9fde, 0xd2ce, 0x04c8, 0xa6dd, 0x0ad8,}, /* 10**-4096 */
|
|
|
|
|
{0x4925, 0x2de4, 0x3436, 0x534f, 0xceae, 0x256b,}, /* 10**-2048 */
|
|
|
|
|
{0x87a6, 0xc0bd, 0xda57, 0x82a5, 0xa2a6, 0x32b5,},
|
|
|
|
|
{0x7133, 0xd21c, 0xdb23, 0xee32, 0x9049, 0x395a,},
|
|
|
|
|
{0xfa91, 0x1939, 0x637a, 0x4325, 0xc031, 0x3cac,},
|
|
|
|
|
{0xac7d, 0xe4a0, 0x64bc, 0x467c, 0xddd0, 0x3e55,},
|
|
|
|
|
{0x3f24, 0xe9a5, 0xa539, 0xea27, 0xa87f, 0x3f2a,},
|
|
|
|
|
{0x67de, 0x94ba, 0x4539, 0x1ead, 0xcfb1, 0x3f94,},
|
|
|
|
|
{0x4c2f, 0xe15b, 0xc44d, 0x94be, 0xe695, 0x3fc9,},
|
|
|
|
|
{0xfdc2, 0xcefc, 0x8461, 0x7711, 0xabcc, 0x3fe4,},
|
|
|
|
|
{0xd3c3, 0x652b, 0xe219, 0x1758, 0xd1b7, 0x3ff1,},
|
|
|
|
|
{0x3d71, 0xd70a, 0x70a3, 0x0a3d, 0xa3d7, 0x3ff8,},
|
|
|
|
|
{0xcccd, 0xcccc, 0xcccc, 0xcccc, 0xcccc, 0x3ffb,}, /* 10**-1 */
|
|
|
|
|
};
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
e24toasc (x, string, ndigs)
|
|
|
|
|
unsigned EMUSHORT x[];
|
|
|
|
|
char *string;
|
|
|
|
|
int ndigs;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT w[NI];
|
|
|
|
|
|
|
|
|
|
e24toe (x, w);
|
|
|
|
|
etoasc (w, string, ndigs);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
e53toasc (x, string, ndigs)
|
|
|
|
|
unsigned EMUSHORT x[];
|
|
|
|
|
char *string;
|
|
|
|
|
int ndigs;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT w[NI];
|
|
|
|
|
|
|
|
|
|
e53toe (x, w);
|
|
|
|
|
etoasc (w, string, ndigs);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
e64toasc (x, string, ndigs)
|
|
|
|
|
unsigned EMUSHORT x[];
|
|
|
|
|
char *string;
|
|
|
|
|
int ndigs;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT w[NI];
|
|
|
|
|
|
|
|
|
|
e64toe (x, w);
|
|
|
|
|
etoasc (w, string, ndigs);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
static char wstring[80]; /* working storage for ASCII output */
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
etoasc (x, string, ndigs)
|
|
|
|
|
unsigned EMUSHORT x[];
|
|
|
|
|
char *string;
|
|
|
|
|
int ndigs;
|
|
|
|
|
{
|
|
|
|
|
EMUSHORT digit;
|
|
|
|
|
unsigned EMUSHORT y[NI], t[NI], u[NI], w[NI];
|
|
|
|
|
unsigned EMUSHORT *p, *r, *ten;
|
|
|
|
|
unsigned EMUSHORT sign;
|
|
|
|
|
int i, j, k, expon, rndsav;
|
|
|
|
|
char *s, *ss;
|
|
|
|
|
unsigned EMUSHORT m;
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
rndsav = rndprc;
|
|
|
|
|
ss = string;
|
|
|
|
|
s = wstring;
|
|
|
|
|
*ss = '\0';
|
|
|
|
|
*s = '\0';
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
if (eisnan (x))
|
|
|
|
|
{
|
|
|
|
|
sprintf (wstring, " NaN ");
|
|
|
|
|
goto bxit;
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
rndprc = NBITS; /* set to full precision */
|
|
|
|
|
emov (x, y); /* retain external format */
|
|
|
|
|
if (y[NE - 1] & 0x8000)
|
|
|
|
|
{
|
|
|
|
|
sign = 0xffff;
|
|
|
|
|
y[NE - 1] &= 0x7fff;
|
|
|
|
|
}
|
|
|
|
|
else
|
|
|
|
|
{
|
|
|
|
|
sign = 0;
|
|
|
|
|
}
|
|
|
|
|
expon = 0;
|
|
|
|
|
ten = &etens[NTEN][0];
|
|
|
|
|
emov (eone, t);
|
|
|
|
|
/* Test for zero exponent */
|
|
|
|
|
if (y[NE - 1] == 0)
|
|
|
|
|
{
|
|
|
|
|
for (k = 0; k < NE - 1; k++)
|
|
|
|
|
{
|
|
|
|
|
if (y[k] != 0)
|
|
|
|
|
goto tnzro; /* denormalized number */
|
|
|
|
|
}
|
|
|
|
|
goto isone; /* legal all zeros */
|
|
|
|
|
}
|
|
|
|
|
tnzro:
|
|
|
|
|
|
|
|
|
|
/* Test for infinity. */
|
|
|
|
|
if (y[NE - 1] == 0x7fff)
|
|
|
|
|
{
|
|
|
|
|
if (sign)
|
|
|
|
|
sprintf (wstring, " -Infinity ");
|
|
|
|
|
else
|
|
|
|
|
sprintf (wstring, " Infinity ");
|
|
|
|
|
goto bxit;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/* Test for exponent nonzero but significand denormalized.
|
|
|
|
|
* This is an error condition.
|
|
|
|
|
*/
|
|
|
|
|
if ((y[NE - 1] != 0) && ((y[NE - 2] & 0x8000) == 0))
|
|
|
|
|
{
|
|
|
|
|
mtherr ("etoasc", DOMAIN);
|
|
|
|
|
sprintf (wstring, "NaN");
|
|
|
|
|
goto bxit;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/* Compare to 1.0 */
|
|
|
|
|
i = ecmp (eone, y);
|
|
|
|
|
if (i == 0)
|
|
|
|
|
goto isone;
|
|
|
|
|
|
|
|
|
|
if (i == -2)
|
|
|
|
|
abort ();
|
|
|
|
|
|
|
|
|
|
if (i < 0)
|
|
|
|
|
{ /* Number is greater than 1 */
|
|
|
|
|
/* Convert significand to an integer and strip trailing decimal zeros. */
|
|
|
|
|
emov (y, u);
|
|
|
|
|
u[NE - 1] = EXONE + NBITS - 1;
|
|
|
|
|
|
|
|
|
|
p = &etens[NTEN - 4][0];
|
|
|
|
|
m = 16;
|
|
|
|
|
do
|
|
|
|
|
{
|
|
|
|
|
ediv (p, u, t);
|
|
|
|
|
efloor (t, w);
|
|
|
|
|
for (j = 0; j < NE - 1; j++)
|
|
|
|
|
{
|
|
|
|
|
if (t[j] != w[j])
|
|
|
|
|
goto noint;
|
|
|
|
|
}
|
|
|
|
|
emov (t, u);
|
|
|
|
|
expon += (int) m;
|
|
|
|
|
noint:
|
|
|
|
|
p += NE;
|
|
|
|
|
m >>= 1;
|
|
|
|
|
}
|
|
|
|
|
while (m != 0);
|
|
|
|
|
|
|
|
|
|
/* Rescale from integer significand */
|
|
|
|
|
u[NE - 1] += y[NE - 1] - (unsigned int) (EXONE + NBITS - 1);
|
|
|
|
|
emov (u, y);
|
|
|
|
|
/* Find power of 10 */
|
|
|
|
|
emov (eone, t);
|
|
|
|
|
m = MAXP;
|
|
|
|
|
p = &etens[0][0];
|
|
|
|
|
/* An unordered compare result shouldn't happen here. */
|
|
|
|
|
while (ecmp (ten, u) <= 0)
|
|
|
|
|
{
|
|
|
|
|
if (ecmp (p, u) <= 0)
|
|
|
|
|
{
|
|
|
|
|
ediv (p, u, u);
|
|
|
|
|
emul (p, t, t);
|
|
|
|
|
expon += (int) m;
|
|
|
|
|
}
|
|
|
|
|
m >>= 1;
|
|
|
|
|
if (m == 0)
|
|
|
|
|
break;
|
|
|
|
|
p += NE;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
else
|
|
|
|
|
{ /* Number is less than 1.0 */
|
|
|
|
|
/* Pad significand with trailing decimal zeros. */
|
|
|
|
|
if (y[NE - 1] == 0)
|
|
|
|
|
{
|
|
|
|
|
while ((y[NE - 2] & 0x8000) == 0)
|
|
|
|
|
{
|
|
|
|
|
emul (ten, y, y);
|
|
|
|
|
expon -= 1;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
else
|
|
|
|
|
{
|
|
|
|
|
emovi (y, w);
|
|
|
|
|
for (i = 0; i < NDEC + 1; i++)
|
|
|
|
|
{
|
|
|
|
|
if ((w[NI - 1] & 0x7) != 0)
|
|
|
|
|
break;
|
|
|
|
|
/* multiply by 10 */
|
|
|
|
|
emovz (w, u);
|
|
|
|
|
eshdn1 (u);
|
|
|
|
|
eshdn1 (u);
|
|
|
|
|
eaddm (w, u);
|
|
|
|
|
u[1] += 3;
|
|
|
|
|
while (u[2] != 0)
|
|
|
|
|
{
|
|
|
|
|
eshdn1 (u);
|
|
|
|
|
u[1] += 1;
|
|
|
|
|
}
|
|
|
|
|
if (u[NI - 1] != 0)
|
|
|
|
|
break;
|
|
|
|
|
if (eone[NE - 1] <= u[1])
|
|
|
|
|
break;
|
|
|
|
|
emovz (u, w);
|
|
|
|
|
expon -= 1;
|
|
|
|
|
}
|
|
|
|
|
emovo (w, y);
|
|
|
|
|
}
|
|
|
|
|
k = -MAXP;
|
|
|
|
|
p = &emtens[0][0];
|
|
|
|
|
r = &etens[0][0];
|
|
|
|
|
emov (y, w);
|
|
|
|
|
emov (eone, t);
|
|
|
|
|
while (ecmp (eone, w) > 0)
|
|
|
|
|
{
|
|
|
|
|
if (ecmp (p, w) >= 0)
|
|
|
|
|
{
|
|
|
|
|
emul (r, w, w);
|
|
|
|
|
emul (r, t, t);
|
|
|
|
|
expon += k;
|
|
|
|
|
}
|
|
|
|
|
k /= 2;
|
|
|
|
|
if (k == 0)
|
|
|
|
|
break;
|
|
|
|
|
p += NE;
|
|
|
|
|
r += NE;
|
|
|
|
|
}
|
|
|
|
|
ediv (t, eone, t);
|
|
|
|
|
}
|
|
|
|
|
isone:
|
|
|
|
|
/* Find the first (leading) digit. */
|
|
|
|
|
emovi (t, w);
|
|
|
|
|
emovz (w, t);
|
|
|
|
|
emovi (y, w);
|
|
|
|
|
emovz (w, y);
|
|
|
|
|
eiremain (t, y);
|
|
|
|
|
digit = equot[NI - 1];
|
|
|
|
|
while ((digit == 0) && (ecmp (y, ezero) != 0))
|
|
|
|
|
{
|
|
|
|
|
eshup1 (y);
|
|
|
|
|
emovz (y, u);
|
|
|
|
|
eshup1 (u);
|
|
|
|
|
eshup1 (u);
|
|
|
|
|
eaddm (u, y);
|
|
|
|
|
eiremain (t, y);
|
|
|
|
|
digit = equot[NI - 1];
|
|
|
|
|
expon -= 1;
|
|
|
|
|
}
|
|
|
|
|
s = wstring;
|
|
|
|
|
if (sign)
|
|
|
|
|
*s++ = '-';
|
|
|
|
|
else
|
|
|
|
|
*s++ = ' ';
|
|
|
|
|
/* Examine number of digits requested by caller. */
|
|
|
|
|
if (ndigs < 0)
|
|
|
|
|
ndigs = 0;
|
|
|
|
|
if (ndigs > NDEC)
|
|
|
|
|
ndigs = NDEC;
|
|
|
|
|
if (digit == 10)
|
|
|
|
|
{
|
|
|
|
|
*s++ = '1';
|
|
|
|
|
*s++ = '.';
|
|
|
|
|
if (ndigs > 0)
|
|
|
|
|
{
|
|
|
|
|
*s++ = '0';
|
|
|
|
|
ndigs -= 1;
|
|
|
|
|
}
|
|
|
|
|
expon += 1;
|
|
|
|
|
}
|
|
|
|
|
else
|
|
|
|
|
{
|
|
|
|
|
*s++ = (char )digit + '0';
|
|
|
|
|
*s++ = '.';
|
|
|
|
|
}
|
|
|
|
|
/* Generate digits after the decimal point. */
|
|
|
|
|
for (k = 0; k <= ndigs; k++)
|
|
|
|
|
{
|
|
|
|
|
/* multiply current number by 10, without normalizing */
|
|
|
|
|
eshup1 (y);
|
|
|
|
|
emovz (y, u);
|
|
|
|
|
eshup1 (u);
|
|
|
|
|
eshup1 (u);
|
|
|
|
|
eaddm (u, y);
|
|
|
|
|
eiremain (t, y);
|
|
|
|
|
*s++ = (char) equot[NI - 1] + '0';
|
|
|
|
|
}
|
|
|
|
|
digit = equot[NI - 1];
|
|
|
|
|
--s;
|
|
|
|
|
ss = s;
|
|
|
|
|
/* round off the ASCII string */
|
|
|
|
|
if (digit > 4)
|
|
|
|
|
{
|
|
|
|
|
/* Test for critical rounding case in ASCII output. */
|
|
|
|
|
if (digit == 5)
|
|
|
|
|
{
|
|
|
|
|
emovo (y, t);
|
|
|
|
|
if (ecmp (t, ezero) != 0)
|
|
|
|
|
goto roun; /* round to nearest */
|
|
|
|
|
if ((*(s - 1) & 1) == 0)
|
|
|
|
|
goto doexp; /* round to even */
|
|
|
|
|
}
|
|
|
|
|
/* Round up and propagate carry-outs */
|
|
|
|
|
roun:
|
|
|
|
|
--s;
|
|
|
|
|
k = *s & 0x7f;
|
|
|
|
|
/* Carry out to most significant digit? */
|
|
|
|
|
if (k == '.')
|
|
|
|
|
{
|
|
|
|
|
--s;
|
|
|
|
|
k = *s;
|
|
|
|
|
k += 1;
|
|
|
|
|
*s = (char) k;
|
|
|
|
|
/* Most significant digit carries to 10? */
|
|
|
|
|
if (k > '9')
|
|
|
|
|
{
|
|
|
|
|
expon += 1;
|
|
|
|
|
*s = '1';
|
|
|
|
|
}
|
|
|
|
|
goto doexp;
|
|
|
|
|
}
|
|
|
|
|
/* Round up and carry out from less significant digits */
|
|
|
|
|
k += 1;
|
|
|
|
|
*s = (char) k;
|
|
|
|
|
if (k > '9')
|
|
|
|
|
{
|
|
|
|
|
*s = '0';
|
|
|
|
|
goto roun;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
doexp:
|
|
|
|
|
/*
|
|
|
|
|
if (expon >= 0)
|
|
|
|
|
sprintf (ss, "e+%d", expon);
|
|
|
|
|
else
|
|
|
|
|
sprintf (ss, "e%d", expon);
|
|
|
|
|
*/
|
|
|
|
|
sprintf (ss, "e%d", expon);
|
|
|
|
|
bxit:
|
|
|
|
|
rndprc = rndsav;
|
|
|
|
|
/* copy out the working string */
|
|
|
|
|
s = string;
|
|
|
|
|
ss = wstring;
|
|
|
|
|
while (*ss == ' ') /* strip possible leading space */
|
|
|
|
|
++ss;
|
|
|
|
|
while ((*s++ = *ss++) != '\0')
|
|
|
|
|
;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; ASCTOQ
|
|
|
|
|
; ASCTOQ.MAC LATEST REV: 11 JAN 84
|
|
|
|
|
; SLM, 3 JAN 78
|
|
|
|
|
;
|
|
|
|
|
; Convert ASCII string to quadruple precision 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).
|
|
|
|
|
;
|
|
|
|
|
; Usage:
|
|
|
|
|
; asctoq (string, q);
|
|
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
/* ASCII to single */
|
|
|
|
|
void
|
|
|
|
|
asctoe24 (s, y)
|
|
|
|
|
char *s;
|
|
|
|
|
unsigned EMUSHORT *y;
|
|
|
|
|
{
|
|
|
|
|
asctoeg (s, y, 24);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* ASCII to double */
|
|
|
|
|
void
|
|
|
|
|
asctoe53 (s, y)
|
|
|
|
|
char *s;
|
|
|
|
|
unsigned EMUSHORT *y;
|
|
|
|
|
{
|
|
|
|
|
#ifdef DEC
|
|
|
|
|
asctoeg (s, y, 56);
|
|
|
|
|
#else
|
|
|
|
|
asctoeg (s, y, 53);
|
|
|
|
|
#endif
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* ASCII to long double */
|
|
|
|
|
void
|
|
|
|
|
asctoe64 (s, y)
|
|
|
|
|
char *s;
|
|
|
|
|
unsigned EMUSHORT *y;
|
|
|
|
|
{
|
|
|
|
|
asctoeg (s, y, 64);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/* ASCII to super double */
|
|
|
|
|
void
|
|
|
|
|
asctoe (s, y)
|
|
|
|
|
char *s;
|
|
|
|
|
unsigned EMUSHORT *y;
|
|
|
|
|
{
|
|
|
|
|
asctoeg (s, y, NBITS);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/* Space to make a copy of the input string: */
|
|
|
|
|
static char lstr[82];
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
asctoeg (ss, y, oprec)
|
|
|
|
|
char *ss;
|
|
|
|
|
unsigned EMUSHORT *y;
|
|
|
|
|
int oprec;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT yy[NI], xt[NI], tt[NI];
|
|
|
|
|
int esign, decflg, sgnflg, nexp, exp, prec, lost;
|
|
|
|
|
int k, trail, c, rndsav;
|
|
|
|
|
EMULONG lexp;
|
|
|
|
|
unsigned EMUSHORT nsign, *p;
|
|
|
|
|
char *sp, *s;
|
|
|
|
|
|
|
|
|
|
/* Copy the input string. */
|
|
|
|
|
s = ss;
|
|
|
|
|
while (*s == ' ') /* skip leading spaces */
|
|
|
|
|
++s;
|
|
|
|
|
sp = lstr;
|
|
|
|
|
for (k = 0; k < 79; k++)
|
|
|
|
|
{
|
|
|
|
|
if ((*sp++ = *s++) == '\0')
|
|
|
|
|
break;
|
|
|
|
|
}
|
|
|
|
|
*sp = '\0';
|
|
|
|
|
s = lstr;
|
|
|
|
|
|
|
|
|
|
rndsav = rndprc;
|
|
|
|
|
rndprc = NBITS; /* Set to full precision */
|
|
|
|
|
lost = 0;
|
|
|
|
|
nsign = 0;
|
|
|
|
|
decflg = 0;
|
|
|
|
|
sgnflg = 0;
|
|
|
|
|
nexp = 0;
|
|
|
|
|
exp = 0;
|
|
|
|
|
prec = 0;
|
|
|
|
|
ecleaz (yy);
|
|
|
|
|
trail = 0;
|
|
|
|
|
|
|
|
|
|
nxtcom:
|
|
|
|
|
k = *s - '0';
|
|
|
|
|
if ((k >= 0) && (k <= 9))
|
|
|
|
|
{
|
|
|
|
|
/* Ignore leading zeros */
|
|
|
|
|
if ((prec == 0) && (decflg == 0) && (k == 0))
|
|
|
|
|
goto donchr;
|
|
|
|
|
/* Identify and strip trailing zeros after the decimal point. */
|
|
|
|
|
if ((trail == 0) && (decflg != 0))
|
|
|
|
|
{
|
|
|
|
|
sp = s;
|
|
|
|
|
while ((*sp >= '0') && (*sp <= '9'))
|
|
|
|
|
++sp;
|
|
|
|
|
/* Check for syntax error */
|
|
|
|
|
c = *sp & 0x7f;
|
|
|
|
|
if ((c != 'e') && (c != 'E') && (c != '\0')
|
|
|
|
|
&& (c != '\n') && (c != '\r') && (c != ' ')
|
|
|
|
|
&& (c != ','))
|
|
|
|
|
goto error;
|
|
|
|
|
--sp;
|
|
|
|
|
while (*sp == '0')
|
|
|
|
|
*sp-- = 'z';
|
|
|
|
|
trail = 1;
|
|
|
|
|
if (*s == 'z')
|
|
|
|
|
goto donchr;
|
|
|
|
|
}
|
|
|
|
|
/* If enough digits were given to more than fill up the yy register,
|
|
|
|
|
* continuing until overflow into the high guard word yy[2]
|
|
|
|
|
* guarantees that there will be a roundoff bit at the top
|
|
|
|
|
* of the low guard word after normalization.
|
|
|
|
|
*/
|
|
|
|
|
if (yy[2] == 0)
|
|
|
|
|
{
|
|
|
|
|
if (decflg)
|
|
|
|
|
nexp += 1; /* count digits after decimal point */
|
|
|
|
|
eshup1 (yy); /* multiply current number by 10 */
|
|
|
|
|
emovz (yy, xt);
|
|
|
|
|
eshup1 (xt);
|
|
|
|
|
eshup1 (xt);
|
|
|
|
|
eaddm (xt, yy);
|
|
|
|
|
ecleaz (xt);
|
|
|
|
|
xt[NI - 2] = (unsigned EMUSHORT) k;
|
|
|
|
|
eaddm (xt, yy);
|
|
|
|
|
}
|
|
|
|
|
else
|
|
|
|
|
{
|
|
|
|
|
lost |= k;
|
|
|
|
|
}
|
|
|
|
|
prec += 1;
|
|
|
|
|
goto donchr;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
switch (*s)
|
|
|
|
|
{
|
|
|
|
|
case 'z':
|
|
|
|
|
break;
|
|
|
|
|
case 'E':
|
|
|
|
|
case 'e':
|
|
|
|
|
goto expnt;
|
|
|
|
|
case '.': /* decimal point */
|
|
|
|
|
if (decflg)
|
|
|
|
|
goto error;
|
|
|
|
|
++decflg;
|
|
|
|
|
break;
|
|
|
|
|
case '-':
|
|
|
|
|
nsign = 0xffff;
|
|
|
|
|
if (sgnflg)
|
|
|
|
|
goto error;
|
|
|
|
|
++sgnflg;
|
|
|
|
|
break;
|
|
|
|
|
case '+':
|
|
|
|
|
if (sgnflg)
|
|
|
|
|
goto error;
|
|
|
|
|
++sgnflg;
|
|
|
|
|
break;
|
|
|
|
|
case ',':
|
|
|
|
|
case ' ':
|
|
|
|
|
case '\0':
|
|
|
|
|
case '\n':
|
|
|
|
|
case '\r':
|
|
|
|
|
goto daldone;
|
|
|
|
|
case 'i':
|
|
|
|
|
case 'I':
|
|
|
|
|
goto infinite;
|
|
|
|
|
default:
|
|
|
|
|
error:
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
einan (yy);
|
|
|
|
|
#else
|
|
|
|
|
mtherr ("asctoe", DOMAIN);
|
|
|
|
|
eclear (yy);
|
|
|
|
|
#endif
|
|
|
|
|
goto aexit;
|
|
|
|
|
}
|
|
|
|
|
donchr:
|
|
|
|
|
++s;
|
|
|
|
|
goto nxtcom;
|
|
|
|
|
|
|
|
|
|
/* Exponent interpretation */
|
|
|
|
|
expnt:
|
|
|
|
|
|
|
|
|
|
esign = 1;
|
|
|
|
|
exp = 0;
|
|
|
|
|
++s;
|
|
|
|
|
/* check for + or - */
|
|
|
|
|
if (*s == '-')
|
|
|
|
|
{
|
|
|
|
|
esign = -1;
|
|
|
|
|
++s;
|
|
|
|
|
}
|
|
|
|
|
if (*s == '+')
|
|
|
|
|
++s;
|
|
|
|
|
while ((*s >= '0') && (*s <= '9'))
|
|
|
|
|
{
|
|
|
|
|
exp *= 10;
|
|
|
|
|
exp += *s++ - '0';
|
|
|
|
|
if (exp > 4956)
|
|
|
|
|
{
|
|
|
|
|
if (esign < 0)
|
|
|
|
|
goto zero;
|
|
|
|
|
else
|
|
|
|
|
goto infinite;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
if (esign < 0)
|
|
|
|
|
exp = -exp;
|
|
|
|
|
if (exp > 4932)
|
|
|
|
|
{
|
|
|
|
|
infinite:
|
|
|
|
|
ecleaz (yy);
|
|
|
|
|
yy[E] = 0x7fff; /* infinity */
|
|
|
|
|
goto aexit;
|
|
|
|
|
}
|
|
|
|
|
if (exp < -4956)
|
|
|
|
|
{
|
|
|
|
|
zero:
|
|
|
|
|
ecleaz (yy);
|
|
|
|
|
goto aexit;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
daldone:
|
|
|
|
|
nexp = exp - nexp;
|
|
|
|
|
/* Pad trailing zeros to minimize power of 10, per IEEE spec. */
|
|
|
|
|
while ((nexp > 0) && (yy[2] == 0))
|
|
|
|
|
{
|
|
|
|
|
emovz (yy, xt);
|
|
|
|
|
eshup1 (xt);
|
|
|
|
|
eshup1 (xt);
|
|
|
|
|
eaddm (yy, xt);
|
|
|
|
|
eshup1 (xt);
|
|
|
|
|
if (xt[2] != 0)
|
|
|
|
|
break;
|
|
|
|
|
nexp -= 1;
|
|
|
|
|
emovz (xt, yy);
|
|
|
|
|
}
|
|
|
|
|
if ((k = enormlz (yy)) > NBITS)
|
|
|
|
|
{
|
|
|
|
|
ecleaz (yy);
|
|
|
|
|
goto aexit;
|
|
|
|
|
}
|
|
|
|
|
lexp = (EXONE - 1 + NBITS) - k;
|
|
|
|
|
emdnorm (yy, lost, 0, lexp, 64);
|
|
|
|
|
/* convert to external format */
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* Multiply by 10**nexp. If precision is 64 bits,
|
|
|
|
|
* the maximum relative error incurred in forming 10**n
|
|
|
|
|
* for 0 <= n <= 324 is 8.2e-20, at 10**180.
|
|
|
|
|
* For 0 <= n <= 999, the peak relative error is 1.4e-19 at 10**947.
|
|
|
|
|
* For 0 >= n >= -999, it is -1.55e-19 at 10**-435.
|
|
|
|
|
*/
|
|
|
|
|
lexp = yy[E];
|
|
|
|
|
if (nexp == 0)
|
|
|
|
|
{
|
|
|
|
|
k = 0;
|
|
|
|
|
goto expdon;
|
|
|
|
|
}
|
|
|
|
|
esign = 1;
|
|
|
|
|
if (nexp < 0)
|
|
|
|
|
{
|
|
|
|
|
nexp = -nexp;
|
|
|
|
|
esign = -1;
|
|
|
|
|
if (nexp > 4096)
|
|
|
|
|
{ /* Punt. Can't handle this without 2 divides. */
|
|
|
|
|
emovi (etens[0], tt);
|
|
|
|
|
lexp -= tt[E];
|
|
|
|
|
k = edivm (tt, yy);
|
|
|
|
|
lexp += EXONE;
|
|
|
|
|
nexp -= 4096;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
p = &etens[NTEN][0];
|
|
|
|
|
emov (eone, xt);
|
|
|
|
|
exp = 1;
|
|
|
|
|
do
|
|
|
|
|
{
|
|
|
|
|
if (exp & nexp)
|
|
|
|
|
emul (p, xt, xt);
|
|
|
|
|
p -= NE;
|
|
|
|
|
exp = exp + exp;
|
|
|
|
|
}
|
|
|
|
|
while (exp <= MAXP);
|
|
|
|
|
|
|
|
|
|
emovi (xt, tt);
|
|
|
|
|
if (esign < 0)
|
|
|
|
|
{
|
|
|
|
|
lexp -= tt[E];
|
|
|
|
|
k = edivm (tt, yy);
|
|
|
|
|
lexp += EXONE;
|
|
|
|
|
}
|
|
|
|
|
else
|
|
|
|
|
{
|
|
|
|
|
lexp += tt[E];
|
|
|
|
|
k = emulm (tt, yy);
|
|
|
|
|
lexp -= EXONE - 1;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
expdon:
|
|
|
|
|
|
|
|
|
|
/* Round and convert directly to the destination type */
|
|
|
|
|
if (oprec == 53)
|
|
|
|
|
lexp -= EXONE - 0x3ff;
|
|
|
|
|
else if (oprec == 24)
|
|
|
|
|
lexp -= EXONE - 0177;
|
|
|
|
|
#ifdef DEC
|
|
|
|
|
else if (oprec == 56)
|
|
|
|
|
lexp -= EXONE - 0201;
|
|
|
|
|
#endif
|
|
|
|
|
rndprc = oprec;
|
|
|
|
|
emdnorm (yy, k, 0, lexp, 64);
|
|
|
|
|
|
|
|
|
|
aexit:
|
|
|
|
|
|
|
|
|
|
rndprc = rndsav;
|
|
|
|
|
yy[0] = nsign;
|
|
|
|
|
switch (oprec)
|
|
|
|
|
{
|
|
|
|
|
#ifdef DEC
|
|
|
|
|
case 56:
|
|
|
|
|
todec (yy, y); /* see etodec.c */
|
|
|
|
|
break;
|
|
|
|
|
#endif
|
|
|
|
|
case 53:
|
|
|
|
|
toe53 (yy, y);
|
|
|
|
|
break;
|
|
|
|
|
case 24:
|
|
|
|
|
toe24 (yy, y);
|
|
|
|
|
break;
|
|
|
|
|
case 64:
|
|
|
|
|
toe64 (yy, y);
|
|
|
|
|
break;
|
|
|
|
|
case NBITS:
|
|
|
|
|
emovo (yy, y);
|
|
|
|
|
break;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* y = largest integer not greater than x
|
|
|
|
|
* (truncated toward minus infinity)
|
|
|
|
|
*
|
|
|
|
|
* unsigned EMUSHORT x[NE], y[NE]
|
|
|
|
|
*
|
|
|
|
|
* efloor (x, y);
|
|
|
|
|
*/
|
|
|
|
|
static unsigned EMUSHORT bmask[] =
|
|
|
|
|
{
|
|
|
|
|
0xffff,
|
|
|
|
|
0xfffe,
|
|
|
|
|
0xfffc,
|
|
|
|
|
0xfff8,
|
|
|
|
|
0xfff0,
|
|
|
|
|
0xffe0,
|
|
|
|
|
0xffc0,
|
|
|
|
|
0xff80,
|
|
|
|
|
0xff00,
|
|
|
|
|
0xfe00,
|
|
|
|
|
0xfc00,
|
|
|
|
|
0xf800,
|
|
|
|
|
0xf000,
|
|
|
|
|
0xe000,
|
|
|
|
|
0xc000,
|
|
|
|
|
0x8000,
|
|
|
|
|
0x0000,
|
|
|
|
|
};
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
efloor (x, y)
|
|
|
|
|
unsigned EMUSHORT x[], y[];
|
|
|
|
|
{
|
|
|
|
|
register unsigned EMUSHORT *p;
|
|
|
|
|
int e, expon, i;
|
|
|
|
|
unsigned EMUSHORT f[NE];
|
|
|
|
|
|
|
|
|
|
emov (x, f); /* leave in external format */
|
|
|
|
|
expon = (int) f[NE - 1];
|
|
|
|
|
e = (expon & 0x7fff) - (EXONE - 1);
|
|
|
|
|
if (e <= 0)
|
|
|
|
|
{
|
|
|
|
|
eclear (y);
|
|
|
|
|
goto isitneg;
|
|
|
|
|
}
|
|
|
|
|
/* number of bits to clear out */
|
|
|
|
|
e = NBITS - e;
|
|
|
|
|
emov (f, y);
|
|
|
|
|
if (e <= 0)
|
|
|
|
|
return;
|
|
|
|
|
|
|
|
|
|
p = &y[0];
|
|
|
|
|
while (e >= 16)
|
|
|
|
|
{
|
|
|
|
|
*p++ = 0;
|
|
|
|
|
e -= 16;
|
|
|
|
|
}
|
|
|
|
|
/* clear the remaining bits */
|
|
|
|
|
*p &= bmask[e];
|
|
|
|
|
/* truncate negatives toward minus infinity */
|
|
|
|
|
isitneg:
|
|
|
|
|
|
|
|
|
|
if ((unsigned EMUSHORT) expon & (unsigned EMUSHORT) 0x8000)
|
|
|
|
|
{
|
|
|
|
|
for (i = 0; i < NE - 1; i++)
|
|
|
|
|
{
|
|
|
|
|
if (f[i] != y[i])
|
|
|
|
|
{
|
|
|
|
|
esub (eone, y, y);
|
|
|
|
|
break;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* unsigned EMUSHORT x[], s[];
|
|
|
|
|
* int *exp;
|
|
|
|
|
*
|
|
|
|
|
* efrexp (x, exp, s);
|
|
|
|
|
*
|
|
|
|
|
* 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.
|
|
|
|
|
*/
|
|
|
|
|
void
|
|
|
|
|
efrexp (x, exp, s)
|
|
|
|
|
unsigned EMUSHORT x[];
|
|
|
|
|
int *exp;
|
|
|
|
|
unsigned EMUSHORT s[];
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT xi[NI];
|
|
|
|
|
EMULONG li;
|
|
|
|
|
|
|
|
|
|
emovi (x, xi);
|
|
|
|
|
li = (EMULONG) ((EMUSHORT) xi[1]);
|
|
|
|
|
|
|
|
|
|
if (li == 0)
|
|
|
|
|
{
|
|
|
|
|
li -= enormlz (xi);
|
|
|
|
|
}
|
|
|
|
|
xi[1] = 0x3ffe;
|
|
|
|
|
emovo (xi, s);
|
|
|
|
|
*exp = (int) (li - 0x3ffe);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* unsigned EMUSHORT x[], y[];
|
|
|
|
|
* long pwr2;
|
|
|
|
|
*
|
|
|
|
|
* eldexp (x, pwr2, y);
|
|
|
|
|
*
|
|
|
|
|
* Returns y = x * 2**pwr2.
|
|
|
|
|
*/
|
|
|
|
|
void
|
|
|
|
|
eldexp (x, pwr2, y)
|
|
|
|
|
unsigned EMUSHORT x[];
|
|
|
|
|
int pwr2;
|
|
|
|
|
unsigned EMUSHORT y[];
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT xi[NI];
|
|
|
|
|
EMULONG li;
|
|
|
|
|
int i;
|
|
|
|
|
|
|
|
|
|
emovi (x, xi);
|
|
|
|
|
li = xi[1];
|
|
|
|
|
li += pwr2;
|
|
|
|
|
i = 0;
|
|
|
|
|
emdnorm (xi, i, i, li, 64);
|
|
|
|
|
emovo (xi, y);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* c = remainder after dividing b by a
|
|
|
|
|
* Least significant integer quotient bits left in equot[].
|
|
|
|
|
*/
|
|
|
|
|
void
|
|
|
|
|
eremain (a, b, c)
|
|
|
|
|
unsigned EMUSHORT a[], b[], c[];
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT den[NI], num[NI];
|
|
|
|
|
|
|
|
|
|
#ifdef NANS
|
|
|
|
|
if ( eisinf (b)
|
|
|
|
|
|| (ecmp (a, ezero) == 0)
|
|
|
|
|
|| eisnan (a)
|
|
|
|
|
|| eisnan (b))
|
|
|
|
|
{
|
|
|
|
|
enan (c);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
#endif
|
|
|
|
|
if (ecmp (a, ezero) == 0)
|
|
|
|
|
{
|
|
|
|
|
mtherr ("eremain", SING);
|
|
|
|
|
eclear (c);
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
emovi (a, den);
|
|
|
|
|
emovi (b, num);
|
|
|
|
|
eiremain (den, num);
|
|
|
|
|
/* Sign of remainder = sign of quotient */
|
|
|
|
|
if (a[0] == b[0])
|
|
|
|
|
num[0] = 0;
|
|
|
|
|
else
|
|
|
|
|
num[0] = 0xffff;
|
|
|
|
|
emovo (num, c);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
eiremain (den, num)
|
|
|
|
|
unsigned EMUSHORT den[], num[];
|
|
|
|
|
{
|
|
|
|
|
EMULONG ld, ln;
|
|
|
|
|
unsigned EMUSHORT j;
|
|
|
|
|
|
|
|
|
|
ld = den[E];
|
|
|
|
|
ld -= enormlz (den);
|
|
|
|
|
ln = num[E];
|
|
|
|
|
ln -= enormlz (num);
|
|
|
|
|
ecleaz (equot);
|
|
|
|
|
while (ln >= ld)
|
|
|
|
|
{
|
|
|
|
|
if (ecmpm (den, num) <= 0)
|
|
|
|
|
{
|
|
|
|
|
esubm (den, num);
|
|
|
|
|
j = 1;
|
|
|
|
|
}
|
|
|
|
|
else
|
|
|
|
|
{
|
|
|
|
|
j = 0;
|
|
|
|
|
}
|
|
|
|
|
eshup1 (equot);
|
|
|
|
|
equot[NI - 1] |= j;
|
|
|
|
|
eshup1 (num);
|
|
|
|
|
ln -= 1;
|
|
|
|
|
}
|
|
|
|
|
emdnorm (num, 0, 0, ln, 0);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/* mtherr.c
|
|
|
|
|
*
|
|
|
|
|
* Library common error handling routine
|
|
|
|
|
*
|
|
|
|
|
*
|
|
|
|
|
*
|
|
|
|
|
* SYNOPSIS:
|
|
|
|
|
*
|
|
|
|
|
* char *fctnam;
|
|
|
|
|
* int code;
|
|
|
|
|
* void mtherr ();
|
|
|
|
|
*
|
|
|
|
|
* mtherr (fctnam, code);
|
|
|
|
|
*
|
|
|
|
|
*
|
|
|
|
|
*
|
|
|
|
|
* DESCRIPTION:
|
|
|
|
|
*
|
|
|
|
|
* This routine may be called to report one of the following
|
|
|
|
|
* error conditions (in the include file mconf.h).
|
|
|
|
|
*
|
|
|
|
|
* Mnemonic Value Significance
|
|
|
|
|
*
|
|
|
|
|
* DOMAIN 1 argument domain error
|
|
|
|
|
* SING 2 function singularity
|
|
|
|
|
* OVERFLOW 3 overflow range error
|
|
|
|
|
* UNDERFLOW 4 underflow range error
|
|
|
|
|
* TLOSS 5 total loss of precision
|
|
|
|
|
* PLOSS 6 partial loss of precision
|
|
|
|
|
* INVALID 7 NaN - producing operation
|
|
|
|
|
* 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.
|
|
|
|
|
*
|
|
|
|
|
* SEE ALSO:
|
|
|
|
|
*
|
|
|
|
|
* mconf.h
|
|
|
|
|
*
|
|
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
Cephes Math Library Release 2.0: April, 1987
|
|
|
|
|
Copyright 1984, 1987 by Stephen L. Moshier
|
|
|
|
|
Direct inquiries to 30 Frost Street, Cambridge, MA 02140
|
|
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
/* include "mconf.h" */
|
|
|
|
|
|
|
|
|
|
/* Notice: the order of appearance of the following
|
|
|
|
|
* messages is bound to the error codes defined
|
|
|
|
|
* in mconf.h.
|
|
|
|
|
*/
|
|
|
|
|
#define NMSGS 8
|
|
|
|
|
static char *ermsg[NMSGS] =
|
|
|
|
|
{
|
|
|
|
|
"unknown", /* error code 0 */
|
|
|
|
|
"domain", /* error code 1 */
|
|
|
|
|
"singularity", /* et seq. */
|
|
|
|
|
"overflow",
|
|
|
|
|
"underflow",
|
|
|
|
|
"total loss of precision",
|
|
|
|
|
"partial loss of precision",
|
|
|
|
|
"invalid operation"
|
|
|
|
|
};
|
|
|
|
|
|
|
|
|
|
int merror = 0;
|
|
|
|
|
extern int merror;
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
mtherr (name, code)
|
|
|
|
|
char *name;
|
|
|
|
|
int code;
|
|
|
|
|
{
|
|
|
|
|
char errstr[80];
|
|
|
|
|
|
|
|
|
|
/* Display string passed by calling program,
|
|
|
|
|
* which is supposed to be the name of the
|
|
|
|
|
* function in which the error occurred.
|
|
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
/* Display error message defined
|
|
|
|
|
* by the code argument.
|
|
|
|
|
*/
|
|
|
|
|
if ((code <= 0) || (code >= NMSGS))
|
|
|
|
|
code = 0;
|
|
|
|
|
sprintf (errstr, " %s %s error", name, ermsg[code]);
|
|
|
|
|
if (extra_warnings)
|
|
|
|
|
warning (errstr);
|
|
|
|
|
/* Set global error message word */
|
|
|
|
|
merror = code + 1;
|
|
|
|
|
|
|
|
|
|
/* Return to calling
|
|
|
|
|
* program
|
|
|
|
|
*/
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/* Here is etodec.c .
|
|
|
|
|
*
|
|
|
|
|
*/
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; convert DEC double precision to e type
|
|
|
|
|
; double d;
|
|
|
|
|
; EMUSHORT e[NE];
|
|
|
|
|
; dectoe (&d, e);
|
|
|
|
|
*/
|
|
|
|
|
void
|
|
|
|
|
dectoe (d, e)
|
|
|
|
|
unsigned EMUSHORT *d;
|
|
|
|
|
unsigned EMUSHORT *e;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT y[NI];
|
|
|
|
|
register unsigned EMUSHORT r, *p;
|
|
|
|
|
|
|
|
|
|
ecleaz (y); /* start with a zero */
|
|
|
|
|
p = y; /* point to our number */
|
|
|
|
|
r = *d; /* get DEC exponent word */
|
|
|
|
|
if (*d & (unsigned int) 0x8000)
|
|
|
|
|
*p = 0xffff; /* fill in our sign */
|
|
|
|
|
++p; /* bump pointer to our exponent word */
|
|
|
|
|
r &= 0x7fff; /* strip the sign bit */
|
|
|
|
|
if (r == 0) /* answer = 0 if high order DEC word = 0 */
|
|
|
|
|
goto done;
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
r >>= 7; /* shift exponent word down 7 bits */
|
|
|
|
|
r += EXONE - 0201; /* subtract DEC exponent offset */
|
|
|
|
|
/* add our e type exponent offset */
|
|
|
|
|
*p++ = r; /* to form our exponent */
|
|
|
|
|
|
|
|
|
|
r = *d++; /* now do the high order mantissa */
|
|
|
|
|
r &= 0177; /* strip off the DEC exponent and sign bits */
|
|
|
|
|
r |= 0200; /* the DEC understood high order mantissa bit */
|
|
|
|
|
*p++ = r; /* put result in our high guard word */
|
|
|
|
|
|
|
|
|
|
*p++ = *d++; /* fill in the rest of our mantissa */
|
|
|
|
|
*p++ = *d++;
|
|
|
|
|
*p = *d;
|
|
|
|
|
|
|
|
|
|
eshdn8 (y); /* shift our mantissa down 8 bits */
|
|
|
|
|
done:
|
|
|
|
|
emovo (y, e);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/*
|
|
|
|
|
; convert e type to DEC double precision
|
|
|
|
|
; double d;
|
|
|
|
|
; EMUSHORT e[NE];
|
|
|
|
|
; etodec (e, &d);
|
|
|
|
|
*/
|
|
|
|
|
#if 0
|
|
|
|
|
static unsigned EMUSHORT decbit[NI] = {0, 0, 0, 0, 0, 0, 0200, 0};
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
etodec (x, d)
|
|
|
|
|
unsigned EMUSHORT *x, *d;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT xi[NI];
|
|
|
|
|
register unsigned EMUSHORT r;
|
|
|
|
|
int i, j;
|
|
|
|
|
|
|
|
|
|
emovi (x, xi);
|
|
|
|
|
*d = 0;
|
|
|
|
|
if (xi[0] != 0)
|
|
|
|
|
*d = 0100000;
|
|
|
|
|
r = xi[E];
|
|
|
|
|
if (r < (EXONE - 128))
|
|
|
|
|
goto zout;
|
|
|
|
|
i = xi[M + 4];
|
|
|
|
|
if ((i & 0200) != 0)
|
|
|
|
|
{
|
|
|
|
|
if ((i & 0377) == 0200)
|
|
|
|
|
{
|
|
|
|
|
if ((i & 0400) != 0)
|
|
|
|
|
{
|
|
|
|
|
/* check all less significant bits */
|
|
|
|
|
for (j = M + 5; j < NI; j++)
|
|
|
|
|
{
|
|
|
|
|
if (xi[j] != 0)
|
|
|
|
|
goto yesrnd;
|
|
|
|
|
}
|
|
|
|
|
}
|
|
|
|
|
goto nornd;
|
|
|
|
|
}
|
|
|
|
|
yesrnd:
|
|
|
|
|
eaddm (decbit, xi);
|
|
|
|
|
r -= enormlz (xi);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
nornd:
|
|
|
|
|
|
|
|
|
|
r -= EXONE;
|
|
|
|
|
r += 0201;
|
|
|
|
|
if (r < 0)
|
|
|
|
|
{
|
|
|
|
|
zout:
|
|
|
|
|
*d++ = 0;
|
|
|
|
|
*d++ = 0;
|
|
|
|
|
*d++ = 0;
|
|
|
|
|
*d++ = 0;
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
if (r >= 0377)
|
|
|
|
|
{
|
|
|
|
|
*d++ = 077777;
|
|
|
|
|
*d++ = -1;
|
|
|
|
|
*d++ = -1;
|
|
|
|
|
*d++ = -1;
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
r &= 0377;
|
|
|
|
|
r <<= 7;
|
|
|
|
|
eshup8 (xi);
|
|
|
|
|
xi[M] &= 0177;
|
|
|
|
|
r |= xi[M];
|
|
|
|
|
*d++ |= r;
|
|
|
|
|
*d++ = xi[M + 1];
|
|
|
|
|
*d++ = xi[M + 2];
|
|
|
|
|
*d++ = xi[M + 3];
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
#else
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
etodec (x, d)
|
|
|
|
|
unsigned EMUSHORT *x, *d;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT xi[NI];
|
|
|
|
|
EMULONG exp;
|
|
|
|
|
int rndsav;
|
|
|
|
|
|
|
|
|
|
emovi (x, xi);
|
|
|
|
|
exp = (EMULONG) xi[E] - (EXONE - 0201); /* adjust exponent for offsets */
|
|
|
|
|
/* round off to nearest or even */
|
|
|
|
|
rndsav = rndprc;
|
|
|
|
|
rndprc = 56;
|
|
|
|
|
emdnorm (xi, 0, 0, exp, 64);
|
|
|
|
|
rndprc = rndsav;
|
|
|
|
|
todec (xi, d);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
todec (x, y)
|
|
|
|
|
unsigned EMUSHORT *x, *y;
|
|
|
|
|
{
|
|
|
|
|
unsigned EMUSHORT i;
|
|
|
|
|
unsigned EMUSHORT *p;
|
|
|
|
|
|
|
|
|
|
p = x;
|
|
|
|
|
*y = 0;
|
|
|
|
|
if (*p++)
|
|
|
|
|
*y = 0100000;
|
|
|
|
|
i = *p++;
|
|
|
|
|
if (i == 0)
|
|
|
|
|
{
|
|
|
|
|
*y++ = 0;
|
|
|
|
|
*y++ = 0;
|
|
|
|
|
*y++ = 0;
|
|
|
|
|
*y++ = 0;
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
if (i > 0377)
|
|
|
|
|
{
|
|
|
|
|
*y++ |= 077777;
|
|
|
|
|
*y++ = 0xffff;
|
|
|
|
|
*y++ = 0xffff;
|
|
|
|
|
*y++ = 0xffff;
|
|
|
|
|
#ifdef ERANGE
|
|
|
|
|
errno = ERANGE;
|
|
|
|
|
#endif
|
|
|
|
|
return;
|
|
|
|
|
}
|
|
|
|
|
i &= 0377;
|
|
|
|
|
i <<= 7;
|
|
|
|
|
eshup8 (x);
|
|
|
|
|
x[M] &= 0177;
|
|
|
|
|
i |= x[M];
|
|
|
|
|
*y++ |= i;
|
|
|
|
|
*y++ = x[M + 1];
|
|
|
|
|
*y++ = x[M + 2];
|
|
|
|
|
*y++ = x[M + 3];
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
#endif /* not 0 */
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/* Output a binary NaN bit pattern in the target machine's format. */
|
|
|
|
|
|
|
|
|
|
/* 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. */
|
|
|
|
|
#ifdef TFMODE_NAN
|
|
|
|
|
TFMODE_NAN;
|
|
|
|
|
#else
|
|
|
|
|
#ifdef MIEEE
|
|
|
|
|
unsigned EMUSHORT TFnan[8] =
|
|
|
|
|
{0x7fff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff};
|
|
|
|
|
#endif
|
|
|
|
|
#ifdef IBMPC
|
|
|
|
|
unsigned EMUSHORT TFnan[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};
|
|
|
|
|
#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};
|
|
|
|
|
#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};
|
|
|
|
|
#endif
|
|
|
|
|
#endif
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
void
|
|
|
|
|
make_nan (nan, mode)
|
|
|
|
|
unsigned EMUSHORT *nan;
|
|
|
|
|
enum machine_mode mode;
|
|
|
|
|
{
|
|
|
|
|
int i, n;
|
|
|
|
|
unsigned EMUSHORT *p;
|
|
|
|
|
|
|
|
|
|
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. */
|
|
|
|
|
#ifndef DEC
|
|
|
|
|
case TFmode:
|
|
|
|
|
n = 8;
|
|
|
|
|
p = TFnan;
|
|
|
|
|
break;
|
|
|
|
|
case XFmode:
|
|
|
|
|
n = 6;
|
|
|
|
|
p = XFnan;
|
|
|
|
|
break;
|
|
|
|
|
case DFmode:
|
|
|
|
|
n = 4;
|
|
|
|
|
p = DFnan;
|
|
|
|
|
break;
|
|
|
|
|
case SFmode:
|
|
|
|
|
n = 2;
|
|
|
|
|
p = SFnan;
|
|
|
|
|
break;
|
|
|
|
|
#endif
|
|
|
|
|
default:
|
|
|
|
|
abort ();
|
|
|
|
|
}
|
|
|
|
|
for (i=0; i < n; i++)
|
|
|
|
|
*nan++ = *p++;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/* Convert an SFmode target `float' value to a REAL_VALUE_TYPE.
|
|
|
|
|
This is the inverse of the function `etarsingle' invoked by
|
|
|
|
|
REAL_VALUE_TO_TARGET_SINGLE. */
|
|
|
|
|
|
|
|
|
|
REAL_VALUE_TYPE
|
|
|
|
|
ereal_from_float (f)
|
|
|
|
|
unsigned long 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 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
|
|
|
|
|
/* 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.
|
|
|
|
|
|
|
|
|
|
The DFmode is stored as an array of longs (i.e., HOST_WIDE_INTs)
|
|
|
|
|
with 32 bits of the value per each long. The first element
|
|
|
|
|
of the input array holds the bits that would come first in the
|
|
|
|
|
target computer's memory. */
|
|
|
|
|
|
|
|
|
|
REAL_VALUE_TYPE
|
|
|
|
|
ereal_from_double (d)
|
|
|
|
|
unsigned long d[];
|
|
|
|
|
{
|
|
|
|
|
REAL_VALUE_TYPE r;
|
|
|
|
|
unsigned EMUSHORT s[4];
|
|
|
|
|
unsigned EMUSHORT e[NE];
|
|
|
|
|
|
|
|
|
|
/* Convert array of 32 bit pieces to equivalent array of 16 bit pieces.
|
|
|
|
|
This is the inverse of `endian'. */
|
|
|
|
|
#if 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
|
|
|
|
|
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);
|
|
|
|
|
#endif
|
|
|
|
|
/* Convert target double to E-type. */
|
|
|
|
|
e53toe (s, e);
|
|
|
|
|
/* Output E-type to REAL_VALUE_TYPE. */
|
|
|
|
|
PUT_REAL (e, &r);
|
|
|
|
|
return r;
|
|
|
|
|
}
|
|
|
|
|
#endif /* EMU_NON_COMPILE not defined */
|