| /* real.c - implementation of REAL_ARITHMETIC, REAL_VALUE_ATOF, |
| and support for XFmode IEEE extended real floating point arithmetic. |
| Copyright (C) 1993, 1994, 1995, 1996 Free Software Foundation, Inc. |
| Contributed by Stephen L. Moshier (moshier@world.std.com). |
| |
| 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, 59 Temple Place - Suite 330, |
| Boston, MA 02111-1307, USA. */ |
| |
| #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, VAX and IBM 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 a floating point |
| arithmetic suite that was not written with gcc in mind. Avoid |
| changing the low-level arithmetic routines unless you have suitable |
| test programs available. A special version of the PARANOIA floating |
| point arithmetic tester, modified for this purpose, can be found on |
| usc.edu: /pub/C-numanal/ieeetest.zoo. Other tests, and libraries of |
| XFmode and TFmode transcendental functions, can be obtained by ftp from |
| netlib.att.com: netlib/cephes. */ |
| |
| /* Type of computer arithmetic. |
| Only one of DEC, IBM, IEEE, or UNK should get defined. |
| |
| `IEEE', when REAL_WORDS_BIG_ENDIAN is non-zero, refers generically |
| to big-endian IEEE floating-point data structure. This definition |
| should work in SFmode `float' type and DFmode `double' type on |
| virtually all big-endian IEEE machines. If LONG_DOUBLE_TYPE_SIZE |
| has been defined to be 96, then IEEE also invokes the particular |
| XFmode (`long double' type) data structure used by the Motorola |
| 680x0 series processors. |
| |
| `IEEE', when REAL_WORDS_BIG_ENDIAN is zero, refers generally to |
| little-endian IEEE machines. In this case, if LONG_DOUBLE_TYPE_SIZE |
| has been defined to be 96, then IEEE also invokes the particular |
| XFmode `long double' data structure used by the Intel 80x86 series |
| processors. |
| |
| `DEC' refers specifically to the Digital Equipment Corp PDP-11 |
| and VAX floating point data structure. This model currently |
| supports no type wider than DFmode. |
| |
| `IBM' refers specifically to the IBM System/370 and compatible |
| floating point data structure. This model currently supports |
| no type wider than DFmode. The IBM conversions were contributed by |
| frank@atom.ansto.gov.au (Frank Crawford). |
| |
| 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 |
| and may deactivate XFmode since `long double' is used to refer |
| to both modes. |
| |
| The macros FLOAT_WORDS_BIG_ENDIAN, HOST_FLOAT_WORDS_BIG_ENDIAN, |
| contributed by Richard Earnshaw <Richard.Earnshaw@cl.cam.ac.uk>, |
| separate the floating point unit's endian-ness from that of |
| the integer addressing. This permits one to define a big-endian |
| FPU on a little-endian machine (e.g., ARM). An extension to |
| BYTES_BIG_ENDIAN may be required for some machines in the future. |
| These optional macros may be defined in tm.h. In real.h, they |
| default to WORDS_BIG_ENDIAN, etc., so there is no need to define |
| them for any normal host or target machine on which the floats |
| and the integers have the same endian-ness. */ |
| |
| |
| /* 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 == IBM_FLOAT_FORMAT |
| /* IBM System/370 style */ |
| #define IBM 1 |
| #else /* it's also not an IBM */ |
| #if TARGET_FLOAT_FORMAT == IEEE_FLOAT_FORMAT |
| #define IEEE |
| #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 IBM */ |
| #endif /* not VAX */ |
| |
| #define REAL_WORDS_BIG_ENDIAN FLOAT_WORDS_BIG_ENDIAN |
| |
| #else |
| /* REAL_ARITHMETIC not defined means that the *host's* data |
| structure will be used. It may differ by endian-ness from the |
| 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 == IBM_FLOAT_FORMAT |
| /* IBM System/370 style */ |
| #define IBM 1 |
| #else /* it's also not an IBM */ |
| #if HOST_FLOAT_FORMAT == IEEE_FLOAT_FORMAT |
| #define IEEE |
| #else /* it's not IEEE either */ |
| unknown arithmetic type |
| #define UNK 1 |
| #endif /* not IEEE */ |
| #endif /* not IBM */ |
| #endif /* not VAX */ |
| |
| #define REAL_WORDS_BIG_ENDIAN HOST_FLOAT_WORDS_BIG_ENDIAN |
| |
| #endif /* REAL_ARITHMETIC not defined */ |
| |
| /* Define INFINITY for support of infinity. |
| Define NANS for support of Not-a-Number's (NaN's). */ |
| #if !defined(DEC) && !defined(IBM) |
| #define INFINITY |
| #define NANS |
| #endif |
| |
| /* Support of NaNs requires support of infinity. */ |
| #ifdef NANS |
| #ifndef INFINITY |
| #define INFINITY |
| #endif |
| #endif |
| |
| /* 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 |
| /* Number of 16 bit words in external e type format */ |
| #define NE 6 |
| #define MAXDECEXP 4932 |
| #define MINDECEXP -4956 |
| #define GET_REAL(r,e) bcopy ((char *) r, (char *) e, 2*NE) |
| #define PUT_REAL(e,r) bcopy ((char *) e, (char *) r, 2*NE) |
| #else /* no XFmode */ |
| #if LONG_DOUBLE_TYPE_SIZE == 128 |
| #define NE 10 |
| #define MAXDECEXP 4932 |
| #define MINDECEXP -4977 |
| #define GET_REAL(r,e) bcopy ((char *) r, (char *) e, 2*NE) |
| #define PUT_REAL(e,r) bcopy ((char *) e, (char *) r, 2*NE) |
| #else |
| #define NE 6 |
| #define MAXDECEXP 4932 |
| #define MINDECEXP -4956 |
| #ifdef REAL_ARITHMETIC |
| /* Emulator uses target format internally |
| but host stores it in host endian-ness. */ |
| |
| #define GET_REAL(r,e) \ |
| do { \ |
| if (HOST_FLOAT_WORDS_BIG_ENDIAN == REAL_WORDS_BIG_ENDIAN) \ |
| e53toe ((unsigned EMUSHORT *) (r), (e)); \ |
| else \ |
| { \ |
| unsigned EMUSHORT w[4]; \ |
| w[3] = ((EMUSHORT *) r)[0]; \ |
| w[2] = ((EMUSHORT *) r)[1]; \ |
| w[1] = ((EMUSHORT *) r)[2]; \ |
| w[0] = ((EMUSHORT *) r)[3]; \ |
| e53toe (w, (e)); \ |
| } \ |
| } while (0) |
| |
| #define PUT_REAL(e,r) \ |
| do { \ |
| if (HOST_FLOAT_WORDS_BIG_ENDIAN == REAL_WORDS_BIG_ENDIAN) \ |
| etoe53 ((e), (unsigned EMUSHORT *) (r)); \ |
| else \ |
| { \ |
| unsigned EMUSHORT w[4]; \ |
| etoe53 ((e), w); \ |
| *((EMUSHORT *) r) = w[3]; \ |
| *((EMUSHORT *) r + 1) = w[2]; \ |
| *((EMUSHORT *) r + 2) = w[1]; \ |
| *((EMUSHORT *) r + 3) = w[0]; \ |
| } \ |
| } while (0) |
| |
| #else /* not REAL_ARITHMETIC */ |
| |
| /* emulator uses host format */ |
| #define GET_REAL(r,e) e53toe ((unsigned EMUSHORT *) (r), (e)) |
| #define PUT_REAL(e,r) etoe53 ((e), (unsigned EMUSHORT *) (r)) |
| |
| #endif /* not REAL_ARITHMETIC */ |
| #endif /* not TFmode */ |
| #endif /* no XFmode */ |
| |
| |
| /* 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) |
| |
| extern int extra_warnings; |
| extern unsigned EMUSHORT ezero[], ehalf[], eone[], etwo[]; |
| extern unsigned EMUSHORT elog2[], esqrt2[]; |
| |
| static void endian PROTO((unsigned EMUSHORT *, long *, |
| enum machine_mode)); |
| static void eclear PROTO((unsigned EMUSHORT *)); |
| static void emov PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void eabs PROTO((unsigned EMUSHORT *)); |
| static void eneg PROTO((unsigned EMUSHORT *)); |
| static int eisneg PROTO((unsigned EMUSHORT *)); |
| static int eisinf PROTO((unsigned EMUSHORT *)); |
| static int eisnan PROTO((unsigned EMUSHORT *)); |
| static void einfin PROTO((unsigned EMUSHORT *)); |
| static void enan PROTO((unsigned EMUSHORT *, int)); |
| static void emovi PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void emovo PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void ecleaz PROTO((unsigned EMUSHORT *)); |
| static void ecleazs PROTO((unsigned EMUSHORT *)); |
| static void emovz PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void einan PROTO((unsigned EMUSHORT *)); |
| static int eiisnan PROTO((unsigned EMUSHORT *)); |
| static int eiisneg PROTO((unsigned EMUSHORT *)); |
| static void eiinfin PROTO((unsigned EMUSHORT *)); |
| static int eiisinf PROTO((unsigned EMUSHORT *)); |
| static int ecmpm PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void eshdn1 PROTO((unsigned EMUSHORT *)); |
| static void eshup1 PROTO((unsigned EMUSHORT *)); |
| static void eshdn8 PROTO((unsigned EMUSHORT *)); |
| static void eshup8 PROTO((unsigned EMUSHORT *)); |
| static void eshup6 PROTO((unsigned EMUSHORT *)); |
| static void eshdn6 PROTO((unsigned EMUSHORT *)); |
| static void eaddm PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void esubm PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void m16m PROTO((unsigned int, unsigned short *, |
| unsigned short *)); |
| static int edivm PROTO((unsigned short *, unsigned short *)); |
| static int emulm PROTO((unsigned short *, unsigned short *)); |
| static void emdnorm PROTO((unsigned EMUSHORT *, int, int, EMULONG, int)); |
| static void esub PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *, |
| unsigned EMUSHORT *)); |
| static void eadd PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *, |
| unsigned EMUSHORT *)); |
| static void eadd1 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *, |
| unsigned EMUSHORT *)); |
| static void ediv PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *, |
| unsigned EMUSHORT *)); |
| static void emul PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *, |
| unsigned EMUSHORT *)); |
| static void e53toe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void e64toe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void e113toe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void e24toe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void etoe113 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void toe113 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void etoe64 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void toe64 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void etoe53 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void toe53 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void etoe24 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void toe24 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static int ecmp PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void eround PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void ltoe PROTO((HOST_WIDE_INT *, unsigned EMUSHORT *)); |
| static void ultoe PROTO((unsigned HOST_WIDE_INT *, unsigned EMUSHORT *)); |
| static void eifrac PROTO((unsigned EMUSHORT *, HOST_WIDE_INT *, |
| unsigned EMUSHORT *)); |
| static void euifrac PROTO((unsigned EMUSHORT *, unsigned HOST_WIDE_INT *, |
| unsigned EMUSHORT *)); |
| static int eshift PROTO((unsigned EMUSHORT *, int)); |
| static int enormlz PROTO((unsigned EMUSHORT *)); |
| static void e24toasc PROTO((unsigned EMUSHORT *, char *, int)); |
| static void e53toasc PROTO((unsigned EMUSHORT *, char *, int)); |
| static void e64toasc PROTO((unsigned EMUSHORT *, char *, int)); |
| static void e113toasc PROTO((unsigned EMUSHORT *, char *, int)); |
| static void etoasc PROTO((unsigned EMUSHORT *, char *, int)); |
| static void asctoe24 PROTO((char *, unsigned EMUSHORT *)); |
| static void asctoe53 PROTO((char *, unsigned EMUSHORT *)); |
| static void asctoe64 PROTO((char *, unsigned EMUSHORT *)); |
| static void asctoe113 PROTO((char *, unsigned EMUSHORT *)); |
| static void asctoe PROTO((char *, unsigned EMUSHORT *)); |
| static void asctoeg PROTO((char *, unsigned EMUSHORT *, int)); |
| static void efloor PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void efrexp PROTO((unsigned EMUSHORT *, int *, |
| unsigned EMUSHORT *)); |
| static void eldexp PROTO((unsigned EMUSHORT *, int, unsigned EMUSHORT *)); |
| static void eremain PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *, |
| unsigned EMUSHORT *)); |
| static void eiremain PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void mtherr PROTO((char *, int)); |
| static void dectoe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void etodec PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void todec PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void ibmtoe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *, |
| enum machine_mode)); |
| static void etoibm PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *, |
| enum machine_mode)); |
| static void toibm PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *, |
| enum machine_mode)); |
| static void make_nan PROTO((unsigned EMUSHORT *, int, enum machine_mode)); |
| static void uditoe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void ditoe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void etoudi PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void etodi PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| static void esqrt PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); |
| |
| /* Copy 32-bit numbers obtained from array containing 16-bit numbers, |
| swapping ends if required, into output array of longs. The |
| result is normally passed to fprintf by the ASM_OUTPUT_ macros. */ |
| |
| static void |
| endian (e, x, mode) |
| unsigned EMUSHORT e[]; |
| long x[]; |
| enum machine_mode mode; |
| { |
| unsigned long th, t; |
| |
| if (REAL_WORDS_BIG_ENDIAN) |
| { |
| switch (mode) |
| { |
| |
| case TFmode: |
| /* Swap halfwords in the fourth long. */ |
| th = (unsigned long) e[6] & 0xffff; |
| t = (unsigned long) e[7] & 0xffff; |
| t |= th << 16; |
| x[3] = (long) t; |
| |
| case XFmode: |
| |
| /* Swap halfwords in the third long. */ |
| th = (unsigned long) e[4] & 0xffff; |
| t = (unsigned long) e[5] & 0xffff; |
| t |= th << 16; |
| x[2] = (long) t; |
| /* fall into the double case */ |
| |
| case DFmode: |
| |
| /* swap halfwords in the second word */ |
| th = (unsigned long) e[2] & 0xffff; |
| t = (unsigned long) e[3] & 0xffff; |
| t |= th << 16; |
| x[1] = (long) t; |
| /* fall into the float case */ |
| |
| case HFmode: |
| case SFmode: |
| |
| /* swap halfwords in the first word */ |
| th = (unsigned long) e[0] & 0xffff; |
| t = (unsigned long) e[1] & 0xffff; |
| t |= th << 16; |
| x[0] = (long) t; |
| break; |
| |
| default: |
| abort (); |
| } |
| } |
| else |
| { |
| /* Pack the output array without swapping. */ |
| |
| switch (mode) |
| { |
| |
| case TFmode: |
| |
| /* Pack the fourth long. */ |
| th = (unsigned long) e[7] & 0xffff; |
| t = (unsigned long) e[6] & 0xffff; |
| t |= th << 16; |
| x[3] = (long) t; |
| |
| case XFmode: |
| |
| /* Pack the third long. |
| Each element of the input REAL_VALUE_TYPE array has 16 useful bits |
| in it. */ |
| th = (unsigned long) e[5] & 0xffff; |
| t = (unsigned long) e[4] & 0xffff; |
| t |= th << 16; |
| x[2] = (long) t; |
| /* fall into the double case */ |
| |
| case DFmode: |
| |
| /* pack the second long */ |
| th = (unsigned long) e[3] & 0xffff; |
| t = (unsigned long) e[2] & 0xffff; |
| t |= th << 16; |
| x[1] = (long) t; |
| /* fall into the float case */ |
| |
| case HFmode: |
| case SFmode: |
| |
| /* pack the first long */ |
| th = (unsigned long) e[1] & 0xffff; |
| t = (unsigned long) e[0] & 0xffff; |
| t |= th << 16; |
| x[0] = (long) t; |
| break; |
| |
| default: |
| abort (); |
| } |
| } |
| } |
| |
| |
| /* 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, eisneg (d1) ^ eisneg (d2)); |
| 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; |
| HOST_WIDE_INT 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 HOST_WIDE_INT 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 HFmode: |
| 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; |
| case TFmode: |
| asctoe113 (s, tem); |
| e113toe (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); |
| eneg (e); |
| PUT_REAL (e, &r); |
| return (r); |
| } |
| |
| |
| /* Round real toward zero to HOST_WIDE_INT; |
| implements REAL_VALUE_FIX (x). */ |
| |
| HOST_WIDE_INT |
| efixi (x) |
| REAL_VALUE_TYPE x; |
| { |
| unsigned EMUSHORT f[NE], g[NE]; |
| HOST_WIDE_INT l; |
| |
| GET_REAL (&x, f); |
| #ifdef NANS |
| if (eisnan (f)) |
| { |
| warning ("conversion from NaN to int"); |
| return (-1); |
| } |
| #endif |
| eifrac (f, &l, g); |
| return l; |
| } |
| |
| /* Round real toward zero to unsigned HOST_WIDE_INT |
| implements REAL_VALUE_UNSIGNED_FIX (x). |
| Negative input returns zero. */ |
| |
| unsigned HOST_WIDE_INT |
| efixui (x) |
| REAL_VALUE_TYPE x; |
| { |
| unsigned EMUSHORT f[NE], g[NE]; |
| unsigned HOST_WIDE_INT l; |
| |
| GET_REAL (&x, f); |
| #ifdef NANS |
| if (eisnan (f)) |
| { |
| warning ("conversion from NaN to unsigned int"); |
| return (-1); |
| } |
| #endif |
| euifrac (f, &l, g); |
| return l; |
| } |
| |
| |
| /* REAL_VALUE_FROM_INT macro. */ |
| |
| void |
| ereal_from_int (d, i, j, mode) |
| REAL_VALUE_TYPE *d; |
| HOST_WIDE_INT i, j; |
| enum machine_mode mode; |
| { |
| unsigned EMUSHORT df[NE], dg[NE]; |
| HOST_WIDE_INT low, high; |
| int sign; |
| |
| if (GET_MODE_CLASS (mode) != MODE_FLOAT) |
| abort (); |
| sign = 0; |
| low = i; |
| if ((high = j) < 0) |
| { |
| sign = 1; |
| /* complement and add 1 */ |
| high = ~high; |
| if (low) |
| low = -low; |
| else |
| high += 1; |
| } |
| eldexp (eone, HOST_BITS_PER_WIDE_INT, df); |
| ultoe ((unsigned HOST_WIDE_INT *) &high, dg); |
| emul (dg, df, dg); |
| ultoe ((unsigned HOST_WIDE_INT *) &low, df); |
| eadd (df, dg, dg); |
| if (sign) |
| eneg (dg); |
| |
| /* A REAL_VALUE_TYPE may not be wide enough to hold the two HOST_WIDE_INTS. |
| Avoid double-rounding errors later by rounding off now from the |
| extra-wide internal format to the requested precision. */ |
| switch (GET_MODE_BITSIZE (mode)) |
| { |
| case 32: |
| etoe24 (dg, df); |
| e24toe (df, dg); |
| break; |
| |
| case 64: |
| etoe53 (dg, df); |
| e53toe (df, dg); |
| break; |
| |
| case 96: |
| etoe64 (dg, df); |
| e64toe (df, dg); |
| break; |
| |
| case 128: |
| etoe113 (dg, df); |
| e113toe (df, dg); |
| break; |
| |
| default: |
| abort (); |
| } |
| |
| PUT_REAL (dg, d); |
| } |
| |
| |
| /* REAL_VALUE_FROM_UNSIGNED_INT macro. */ |
| |
| void |
| ereal_from_uint (d, i, j, mode) |
| REAL_VALUE_TYPE *d; |
| unsigned HOST_WIDE_INT i, j; |
| enum machine_mode mode; |
| { |
| unsigned EMUSHORT df[NE], dg[NE]; |
| unsigned HOST_WIDE_INT low, high; |
| |
| if (GET_MODE_CLASS (mode) != MODE_FLOAT) |
| abort (); |
| low = i; |
| high = j; |
| eldexp (eone, HOST_BITS_PER_WIDE_INT, df); |
| ultoe (&high, dg); |
| emul (dg, df, dg); |
| ultoe (&low, df); |
| eadd (df, dg, dg); |
| |
| /* A REAL_VALUE_TYPE may not be wide enough to hold the two HOST_WIDE_INTS. |
| Avoid double-rounding errors later by rounding off now from the |
| extra-wide internal format to the requested precision. */ |
| switch (GET_MODE_BITSIZE (mode)) |
| { |
| case 32: |
| etoe24 (dg, df); |
| e24toe (df, dg); |
| break; |
| |
| case 64: |
| etoe53 (dg, df); |
| e53toe (df, dg); |
| break; |
| |
| case 96: |
| etoe64 (dg, df); |
| e64toe (df, dg); |
| break; |
| |
| case 128: |
| etoe113 (dg, df); |
| e113toe (df, dg); |
| break; |
| |
| default: |
| abort (); |
| } |
| |
| PUT_REAL (dg, d); |
| } |
| |
| |
| /* REAL_VALUE_TO_INT macro. */ |
| |
| void |
| ereal_to_int (low, high, rr) |
| HOST_WIDE_INT *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_WIDE_INT, df); |
| ediv (df, d, dg); /* dg = d / 2^32 is the high word */ |
| euifrac (dg, (unsigned HOST_WIDE_INT *) high, dh); |
| emul (df, dh, dg); /* fractional part is the low word */ |
| euifrac (dg, (unsigned HOST_WIDE_INT *)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 just checks the sign bit, so that -0 counts as negative. */ |
| |
| int |
| target_negative (x) |
| REAL_VALUE_TYPE x; |
| { |
| return ereal_isneg (x); |
| } |
| |
| /* 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 TFmode: |
| etoe113 (e, t); |
| e113toe (t, t); |
| break; |
| |
| case XFmode: |
| etoe64 (e, t); |
| e64toe (t, t); |
| break; |
| |
| case DFmode: |
| etoe53 (e, t); |
| e53toe (t, t); |
| break; |
| |
| case HFmode: |
| case SFmode: |
| etoe24 (e, t); |
| e24toe (t, t); |
| break; |
| |
| case SImode: |
| r = etrunci (arg); |
| return (r); |
| |
| /* If an unsupported type was requested, presume that |
| the machine files know something useful to do with |
| the unmodified value. */ |
| |
| default: |
| return (arg); |
| } |
| PUT_REAL (t, &r); |
| return (r); |
| } |
| |
| /* Try to change R into its exact multiplicative inverse in machine mode |
| MODE. Return nonzero function value if successful. */ |
| |
| int |
| exact_real_inverse (mode, r) |
| enum machine_mode mode; |
| REAL_VALUE_TYPE *r; |
| { |
| unsigned EMUSHORT e[NE], einv[NE]; |
| REAL_VALUE_TYPE rinv; |
| int i; |
| |
| GET_REAL (r, e); |
| |
| /* Test for input in range. Don't transform IEEE special values. */ |
| if (eisinf (e) || eisnan (e) || (ecmp (e, ezero) == 0)) |
| return 0; |
| |
| /* Test for a power of 2: all significand bits zero except the MSB. |
| We are assuming the target has binary (or hex) arithmetic. */ |
| if (e[NE - 2] != 0x8000) |
| return 0; |
| |
| for (i = 0; i < NE - 2; i++) |
| { |
| if (e[i] != 0) |
| return 0; |
| } |
| |
| /* Compute the inverse and truncate it to the required mode. */ |
| ediv (e, eone, einv); |
| PUT_REAL (einv, &rinv); |
| rinv = real_value_truncate (mode, rinv); |
| |
| #ifdef CHECK_FLOAT_VALUE |
| /* This check is not redundant. It may, for example, flush |
| a supposedly IEEE denormal value to zero. */ |
| i = 0; |
| if (CHECK_FLOAT_VALUE (mode, rinv, i)) |
| return 0; |
| #endif |
| GET_REAL (&rinv, einv); |
| |
| /* Check the bits again, because the truncation might have |
| generated an arbitrary saturation value on overflow. */ |
| if (einv[NE - 2] != 0x8000) |
| return 0; |
| |
| for (i = 0; i < NE - 2; i++) |
| { |
| if (einv[i] != 0) |
| return 0; |
| } |
| |
| /* Fail if the computed inverse is out of range. */ |
| if (eisinf (einv) || eisnan (einv) || (ecmp (einv, ezero) == 0)) |
| return 0; |
| |
| /* Output the reciprocal and return success flag. */ |
| PUT_REAL (einv, r); |
| return 1; |
| } |
| #endif /* REAL_ARITHMETIC defined */ |
| |
| /* Used for debugging--print the value of R in human-readable format |
| on stderr. */ |
| |
| void |
| debug_real (r) |
| REAL_VALUE_TYPE r; |
| { |
| char dstr[30]; |
| |
| REAL_VALUE_TO_DECIMAL (r, "%.20g", dstr); |
| fprintf (stderr, "%s", dstr); |
| } |
| |
| |
| /* The following routines convert REAL_VALUE_TYPE to the various floating |
| point formats that are meaningful to supported computers. |
| |
| The results are returned in 32-bit pieces, each piece stored in a `long'. |
| This is so they can be printed by statements like |
| |
| fprintf (file, "%lx, %lx", L[0], L[1]); |
| |
| that will work on both narrow- and wide-word host computers. */ |
| |
| /* Convert R to a 128-bit long double precision value. The output array L |
| contains four 32-bit pieces of the result, in the order they would appear |
| in memory. */ |
| |
| void |
| etartdouble (r, l) |
| REAL_VALUE_TYPE r; |
| long l[]; |
| { |
| unsigned EMUSHORT e[NE]; |
| |
| GET_REAL (&r, e); |
| etoe113 (e, e); |
| endian (e, l, TFmode); |
| } |
| |
| /* Convert R to a double extended precision value. The output array L |
| contains three 32-bit pieces of the result, in the order they would |
| appear in memory. */ |
| |
| void |
| etarldouble (r, l) |
| REAL_VALUE_TYPE r; |
| long l[]; |
| { |
| unsigned EMUSHORT e[NE]; |
| |
| GET_REAL (&r, e); |
| etoe64 (e, e); |
| endian (e, l, XFmode); |
| } |
| |
| /* Convert R to a double precision value. The output array L contains two |
| 32-bit pieces of the result, in the order they would appear in memory. */ |
| |
| void |
| etardouble (r, l) |
| REAL_VALUE_TYPE r; |
| long l[]; |
| { |
| unsigned EMUSHORT e[NE]; |
| |
| GET_REAL (&r, e); |
| etoe53 (e, e); |
| endian (e, l, DFmode); |
| } |
| |
| /* Convert R to a single precision float value stored in the least-significant |
| bits of a `long'. */ |
| |
| long |
| etarsingle (r) |
| REAL_VALUE_TYPE r; |
| { |
| unsigned EMUSHORT e[NE]; |
| long l; |
| |
| GET_REAL (&r, e); |
| etoe24 (e, e); |
| endian (e, &l, SFmode); |
| return ((long) l); |
| } |
| |
| /* Convert X to a decimal ASCII string S for output to an assembly |
| language file. Note, there is no standard way to spell infinity or |
| a NaN, so these values may require special treatment in the tm.h |
| macros. */ |
| |
| void |
| ereal_to_decimal (x, s) |
| REAL_VALUE_TYPE x; |
| char *s; |
| { |
| unsigned EMUSHORT e[NE]; |
| |
| GET_REAL (&x, e); |
| etoasc (e, s, 20); |
| } |
| |
| /* Compare X and Y. Return 1 if X > Y, 0 if X == Y, -1 if X < Y, |
| or -2 if either is a NaN. */ |
| |
| int |
| ereal_cmp (x, y) |
| REAL_VALUE_TYPE x, y; |
| { |
| unsigned EMUSHORT ex[NE], ey[NE]; |
| |
| GET_REAL (&x, ex); |
| GET_REAL (&y, ey); |
| return (ecmp (ex, ey)); |
| } |
| |
| /* Return 1 if the sign bit of X is set, else return 0. */ |
| |
| int |
| ereal_isneg (x) |
| REAL_VALUE_TYPE x; |
| { |
| unsigned EMUSHORT ex[NE]; |
| |
| GET_REAL (&x, ex); |
| return (eisneg (ex)); |
| } |
| |
| /* End of REAL_ARITHMETIC interface */ |
| |
| /* |
| 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, similar to 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 exploded e-type data structure of a number (a "word" is 16 bits): |
| |
| ei[0] sign word (0 for positive, 0xffff for negative) |
| ei[1] biased exponent (value = EXONE for the number 1.0) |
| 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 e-type 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 |
| e113toe (&d, e) 128-bit 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 HOST_WIDE_INT and e type fraction |
| euifrac (e, &l, frac) e to unsigned HOST_WIDE_INT 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) 80-bit long double to ASCII string |
| e113toasc (&d, str, n) 128-bit 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) HOST_WIDE_INT to e type |
| ultoe (&l, e) unsigned HOST_WIDE_INT 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 exploded e-type 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 |
| eiisneg (ai) 1 if sign bit of ai != 0, else 0 |
| 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). */ |
| |
| /* Definitions for error codes that are passed to the common error handling |
| routine mtherr. |
| |
| 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 IEEE should be defined. |
| These numbers have 53-bit significands. In this mode, constants |
| are provided as arrays of hexadecimal 16 bit integers. |
| The endian-ness of generated values is controlled by |
| REAL_WORDS_BIG_ENDIAN. |
| |
| To accommodate other types of computer arithmetic, all |
| constants are also provided in a normal decimal radix |
| 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 */ |
| |
| #if LONG_DOUBLE_TYPE_SIZE == 128 |
| /* 0.0 */ |
| unsigned EMUSHORT ezero[NE] = |
| {0x0000, 0x0000, 0x0000, 0x0000, |
| 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000,}; |
| extern unsigned EMUSHORT ezero[]; |
| |
| /* 5.0E-1 */ |
| unsigned EMUSHORT ehalf[NE] = |
| {0x0000, 0x0000, 0x0000, 0x0000, |
| 0x0000, 0x0000, 0x0000, 0x0000, 0x8000, 0x3ffe,}; |
| extern unsigned EMUSHORT ehalf[]; |
| |
| /* 1.0E0 */ |
| unsigned EMUSHORT eone[NE] = |
| {0x0000, 0x0000, 0x0000, 0x0000, |
| 0x0000, 0x0000, 0x0000, 0x0000, 0x8000, 0x3fff,}; |
| extern unsigned EMUSHORT eone[]; |
| |
| /* 2.0E0 */ |
| unsigned EMUSHORT etwo[NE] = |
| {0x0000, 0x0000, 0x0000, 0x0000, |
| 0x0000, 0x0000, 0x0000, 0x0000, 0x8000, 0x4000,}; |
| extern unsigned EMUSHORT etwo[]; |
| |
| /* 3.2E1 */ |
| unsigned EMUSHORT e32[NE] = |
| {0x0000, 0x0000, 0x0000, 0x0000, |
| 0x0000, 0x0000, 0x0000, 0x0000, 0x8000, 0x4004,}; |
| extern unsigned EMUSHORT e32[]; |
| |
| /* 6.93147180559945309417232121458176568075500134360255E-1 */ |
| unsigned EMUSHORT elog2[NE] = |
| {0x40f3, 0xf6af, 0x03f2, 0xb398, |
| 0xc9e3, 0x79ab, 0150717, 0013767, 0130562, 0x3ffe,}; |
| extern unsigned EMUSHORT elog2[]; |
| |
| /* 1.41421356237309504880168872420969807856967187537695E0 */ |
| unsigned EMUSHORT esqrt2[NE] = |
| {0x1d6f, 0xbe9f, 0x754a, 0x89b3, |
| 0x597d, 0x6484, 0174736, 0171463, 0132404, 0x3fff,}; |
| extern unsigned EMUSHORT esqrt2[]; |
| |
| /* 3.14159265358979323846264338327950288419716939937511E0 */ |
| unsigned EMUSHORT epi[NE] = |
| {0x2902, 0x1cd1, 0x80dc, 0x628b, |
| 0xc4c6, 0xc234, 0020550, 0155242, 0144417, 0040000,}; |
| extern unsigned EMUSHORT epi[]; |
| |
| #else |
| /* LONG_DOUBLE_TYPE_SIZE is other than 128 */ |
| unsigned EMUSHORT ezero[NE] = |
| {0, 0000000, 0000000, 0000000, 0000000, 0000000,}; |
| unsigned EMUSHORT ehalf[NE] = |
| {0, 0000000, 0000000, 0000000, 0100000, 0x3ffe,}; |
| unsigned EMUSHORT eone[NE] = |
| {0, 0000000, 0000000, 0000000, 0100000, 0x3fff,}; |
| unsigned EMUSHORT etwo[NE] = |
| {0, 0000000, 0000000, 0000000, 0100000, 0040000,}; |
| unsigned EMUSHORT e32[NE] = |
| {0, 0000000, 0000000, 0000000, 0100000, 0040004,}; |
| unsigned EMUSHORT elog2[NE] = |
| {0xc9e4, 0x79ab, 0150717, 0013767, 0130562, 0x3ffe,}; |
| unsigned EMUSHORT esqrt2[NE] = |
| {0x597e, 0x6484, 0174736, 0171463, 0132404, 0x3fff,}; |
| unsigned EMUSHORT epi[NE] = |
| {0xc4c6, 0xc234, 0020550, 0155242, 0144417, 0040000,}; |
| #endif |
| |
| /* Control register for rounding precision. |
| This can be set to 113 (if NE=10), 80 (if NE=6), 64, 56, 53, or 24 bits. */ |
| |
| int rndprc = NBITS; |
| extern int rndprc; |
| |
| /* Clear out entire e-type number X. */ |
| |
| static void |
| eclear (x) |
| register unsigned EMUSHORT *x; |
| { |
| register int i; |
| |
| for (i = 0; i < NE; i++) |
| *x++ = 0; |
| } |
| |
| /* Move e-type number from A to B. */ |
| |
| static void |
| emov (a, b) |
| register unsigned EMUSHORT *a, *b; |
| { |
| register int i; |
| |
| for (i = 0; i < NE; i++) |
| *b++ = *a++; |
| } |
| |
| |
| /* Absolute value of e-type X. */ |
| |
| static void |
| eabs (x) |
| unsigned EMUSHORT x[]; |
| { |
| /* sign is top bit of last word of external format */ |
| x[NE - 1] &= 0x7fff; |
| } |
| |
| /* Negate the e-type number X. */ |
| |
| static void |
| eneg (x) |
| unsigned EMUSHORT x[]; |
| { |
| |
| x[NE - 1] ^= 0x8000; /* Toggle the sign bit */ |
| } |
| |
| /* Return 1 if sign bit of e-type number X is nonzero, else zero. */ |
| |
| static int |
| eisneg (x) |
| unsigned EMUSHORT x[]; |
| { |
| |
| if (x[NE - 1] & 0x8000) |
| return (1); |
| else |
| return (0); |
| } |
| |
| /* Return 1 if e-type number X is infinity, else return zero. */ |
| |
| static 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. */ |
| |
| static 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 e-type number X with infinity pattern (IEEE) |
| or largest possible number (non-IEEE). */ |
| |
| static 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 == 113) |
| { |
| *(x - 9) = 0; |
| *(x - 8) = 0; |
| } |
| 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. */ |
| |
| static void |
| enan (x, sign) |
| register unsigned EMUSHORT *x; |
| int sign; |
| { |
| register int i; |
| |
| for (i = 0; i < NE - 2; i++) |
| *x++ = 0; |
| *x++ = 0xc000; |
| *x = (sign << 15) | 0x7fff; |
| } |
| |
| /* Move in an e-type number A, converting it to exploded e-type B. */ |
| |
| static 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 out exploded e-type number A, converting it to e type B. */ |
| |
| static void |
| emovo (a, b) |
| unsigned EMUSHORT *a, *b; |
| { |
| register unsigned EMUSHORT *p, *q; |
| unsigned EMUSHORT i; |
| int j; |
| |
| 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, eiisneg (a)); |
| return; |
| } |
| #endif |
| einfin (b); |
| return; |
| } |
| #endif |
| /* skip over guard word */ |
| ++p; |
| /* move the significand */ |
| for (j = 0; j < NE - 1; j++) |
| *q-- = *p++; |
| } |
| |
| /* Clear out exploded e-type number XI. */ |
| |
| static void |
| ecleaz (xi) |
| register unsigned EMUSHORT *xi; |
| { |
| register int i; |
| |
| for (i = 0; i < NI; i++) |
| *xi++ = 0; |
| } |
| |
| /* Clear out exploded e-type XI, but don't touch the sign. */ |
| |
| static void |
| ecleazs (xi) |
| register unsigned EMUSHORT *xi; |
| { |
| register int i; |
| |
| ++xi; |
| for (i = 0; i < NI - 1; i++) |
| *xi++ = 0; |
| } |
| |
| /* Move exploded e-type number from A to B. */ |
| |
| static 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 exploded e-type NaN. |
| The explicit pattern for this is maximum exponent and |
| top two significant bits set. */ |
| |
| static void |
| einan (x) |
| unsigned EMUSHORT x[]; |
| { |
| |
| ecleaz (x); |
| x[E] = 0x7fff; |
| x[M + 1] = 0xc000; |
| } |
| |
| /* Return nonzero if exploded e-type X is a NaN. */ |
| |
| static 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); |
| } |
| |
| /* Return nonzero if sign of exploded e-type X is nonzero. */ |
| |
| static int |
| eiisneg (x) |
| unsigned EMUSHORT x[]; |
| { |
| |
| return x[0] != 0; |
| } |
| |
| /* Fill exploded e-type X with infinity pattern. |
| This has maximum exponent and significand all zeros. */ |
| |
| static void |
| eiinfin (x) |
| unsigned EMUSHORT x[]; |
| { |
| |
| ecleaz (x); |
| x[E] = 0x7fff; |
| } |
| |
| /* Return nonzero if exploded e-type X is infinite. */ |
| |
| static 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 exploded e-type format. |
| Guard words are included in the comparison. |
| |
| Returns +1 if a > b |
| 0 if a == b |
| -1 if a < b */ |
| |
| static 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 of exploded e-type X down by 1 bit. */ |
| |
| static 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 of exploded e-type X up by 1 bit. */ |
| |
| static 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 of exploded e-type X down by 8 bits. */ |
| |
| static 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 of exploded e-type X up by 8 bits. */ |
| |
| static 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 of exploded e-type X up by 16 bits. */ |
| |
| static 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 of exploded e-type X down by 16 bits. */ |
| |
| static 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 of exploded e-type X and Y. X + Y replaces Y. */ |
| |
| static 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 of exploded e-type X and Y. Y - X replaces Y. */ |
| |
| static 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; |
| } |
| } |
| |
| |
| static unsigned EMUSHORT equot[NI]; |
| |
| |
| #if 0 |
| /* Radix 2 shift-and-add versions of multiply and divide */ |
| |
| |
| /* Divide significands */ |
| |
| 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 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); |
| } |
| |
| #else |
| |
| /* Radix 65536 versions of multiply and divide. */ |
| |
| /* Multiply significand of e-type number B |
| by 16-bit quantity A, return e-type result to C. */ |
| |
| static void |
| m16m (a, b, c) |
| unsigned int a; |
| unsigned EMUSHORT b[], c[]; |
| { |
| register unsigned EMUSHORT *pp; |
| register unsigned EMULONG carry; |
| unsigned EMUSHORT *ps; |
| unsigned EMUSHORT p[NI]; |
| unsigned EMULONG aa, m; |
| int i; |
| |
| aa = a; |
| pp = &p[NI-2]; |
| *pp++ = 0; |
| *pp = 0; |
| ps = &b[NI-1]; |
| |
| for (i=M+1; i<NI; i++) |
| { |
| if (*ps == 0) |
| { |
| --ps; |
| --pp; |
| *(pp-1) = 0; |
| } |
| else |
| { |
| m = (unsigned EMULONG) aa * *ps--; |
| carry = (m & 0xffff) + *pp; |
| *pp-- = (unsigned EMUSHORT)carry; |
| carry = (carry >> 16) + (m >> 16) + *pp; |
| *pp = (unsigned EMUSHORT)carry; |
| *(pp-1) = carry >> 16; |
| } |
| } |
| for (i=M; i<NI; i++) |
| c[i] = p[i]; |
| } |
| |
| /* Divide significands of exploded e-types NUM / DEN. Neither the |
| numerator NUM nor the denominator DEN is permitted to have its high guard |
| word nonzero. */ |
| |
| static int |
| edivm (den, num) |
| unsigned EMUSHORT den[], num[]; |
| { |
| int i; |
| register unsigned EMUSHORT *p; |
| unsigned EMULONG tnum; |
| unsigned EMUSHORT j, tdenm, tquot; |
| unsigned EMUSHORT tprod[NI+1]; |
| |
| p = &equot[0]; |
| *p++ = num[0]; |
| *p++ = num[1]; |
| |
| for (i=M; i<NI; i++) |
| { |
| *p++ = 0; |
| } |
| eshdn1 (num); |
| tdenm = den[M+1]; |
| for (i=M; i<NI; i++) |
| { |
| /* Find trial quotient digit (the radix is 65536). */ |
| tnum = (((unsigned EMULONG) num[M]) << 16) + num[M+1]; |
| |
| /* Do not execute the divide instruction if it will overflow. */ |
| if ((tdenm * 0xffffL) < tnum) |
| tquot = 0xffff; |
| else |
| tquot = tnum / tdenm; |
| /* Multiply denominator by trial quotient digit. */ |
| m16m ((unsigned int)tquot, den, tprod); |
| /* The quotient digit may have been overestimated. */ |
| if (ecmpm (tprod, num) > 0) |
| { |
| tquot -= 1; |
| esubm (den, tprod); |
| if (ecmpm (tprod, num) > 0) |
| { |
| tquot -= 1; |
| esubm (den, tprod); |
| } |
| } |
| esubm (tprod, num); |
| equot[i] = tquot; |
| eshup6(num); |
| } |
| /* 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 of exploded e-type A and B, result in B. */ |
| |
| static int |
| emulm (a, b) |
| unsigned EMUSHORT a[], b[]; |
| { |
| unsigned EMUSHORT *p, *q; |
| unsigned EMUSHORT pprod[NI]; |
| unsigned EMUSHORT j; |
| int i; |
| |
| equot[0] = b[0]; |
| equot[1] = b[1]; |
| for (i=M; i<NI; i++) |
| equot[i] = 0; |
| |
| j = 0; |
| p = &a[NI-1]; |
| q = &equot[NI-1]; |
| for (i=M+1; i<NI; i++) |
| { |
| if (*p == 0) |
| { |
| --p; |
| } |
| else |
| { |
| m16m ((unsigned int) *p--, b, pprod); |
| eaddm(pprod, equot); |
| } |
| j |= *q; |
| eshdn6(equot); |
| } |
| |
| for (i=0; i<NI; i++) |
| b[i] = equot[i]; |
| |
| /* return flag for lost nonzero bits */ |
| return ((int)j); |
| } |
| #endif |
| |
| |
| /* Normalize and round off. |
| |
| The internal format number to be rounded is S. |
| Input LOST is 0 if the value is exact. This is the so-called sticky bit. |
| |
| Input SUBFLG indicates whether the number was obtained |
| by a subtraction operation. In that case if LOST is nonzero |
| 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. If it is nonzero, the |
| returned value will be rounded to RNDPRC bits. |
| |
| For future reference: In order for emdnorm to round off denormal |
| significands at the right point, the input exponent must be |
| adjusted to be the actual value it would have after conversion to |
| the final floating point type. This adjustment has been |
| implemented for all type conversions (etoe53, etc.) and decimal |
| conversions, but not for the arithmetic functions (eadd, etc.). |
| Data types having standard 15-bit exponents are not affected by |
| this, but SFmode and DFmode are affected. For example, ediv with |
| rndprc = 24 will not round correctly to 24-bit precision if the |
| result is denormal. */ |
| |
| 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]; |
| |
| static 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; |
| re = rw - 1; |
| rebit = 1; |
| break; |
| case 113: |
| rw = 10; |
| rmsk = 0x7fff; |
| rmbit = 0x4000; |
| rebit = 0x8000; |
| re = rw; |
| break; |
| case 64: |
| rw = 7; |
| rmsk = 0xffff; |
| rmbit = 0x8000; |
| re = rw - 1; |
| rebit = 1; |
| break; |
| /* For DEC or IBM arithmetic */ |
| case 56: |
| rw = 6; |
| rmsk = 0xff; |
| rmbit = 0x80; |
| rebit = 0x100; |
| re = rw; |
| break; |
| case 53: |
| rw = 6; |
| rmsk = 0x7ff; |
| rmbit = 0x0400; |
| rebit = 0x800; |
| re = rw; |
| break; |
| case 24: |
| rw = 4; |
| rmsk = 0xff; |
| rmbit = 0x80; |
| rebit = 0x100; |
| re = rw; |
| break; |
| } |
| rbit[re] = rebit; |
| rlast = rndprc; |
| } |
| |
| /* Shift down 1 temporarily if the data structure has an implied |
| most significant bit and the number is denormal. |
| Intel long double denormals also lose one bit of precision. */ |
| if ((exp <= 0) && (rndprc != NBITS) |
| && ((rndprc != 64) || ((rndprc == 64) && ! REAL_WORDS_BIG_ENDIAN))) |
| { |
| lost |= s[NI - 1] & 1; |
| eshdn1 (s); |
| } |
| /* Clear out all bits below the rounding bit, |
| remembering in r if any were nonzero. */ |
| r = s[rw] & rmsk; |
| if (rndprc < NBITS) |
| { |
| i = rw + 1; |
| while (i < NI) |
| { |
| if (s[i]) |
| r |= 1; |
| s[i] = 0; |
| ++i; |
| } |
| } |
| 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: |
| /* Undo the temporary shift for denormal values. */ |
| if ((exp <= 0) && (rndprc != NBITS) |
| && ((rndprc != 64) || ((rndprc == 64) && ! REAL_WORDS_BIG_ENDIAN))) |
| { |
| 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) || (rndprc == 113)) |
| { |
| 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. C = B - A, all e type numbers. */ |
| |
| static int subflg = 0; |
| |
| static 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, 0); |
| return; |
| } |
| #endif |
| subflg = 1; |
| eadd1 (a, b, c); |
| } |
| |
| /* Add. C = A + B, all e type. */ |
| |
| static 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, 0); |
| return; |
| } |
| #endif |
| subflg = 0; |
| eadd1 (a, b, c); |
| } |
| |
| /* Arithmetic common to both addition and subtraction. */ |
| |
| static 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 denormalized 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) |
| { |
| ltb += 1; |
| if (ltb >= 0x7fff) |
| { |
| eclear (c); |
| if (ai[0] != 0) |
| eneg (c); |
| einfin (c); |
| return; |
| } |
| 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: C = B/A, all e type. */ |
| |
| static void |
| ediv (a, b, c) |
| unsigned EMUSHORT *a, *b, *c; |
| { |
| unsigned EMUSHORT ai[NI], bi[NI]; |
| int i, sign; |
| EMULONG lt, lta, ltb; |
| |
| /* IEEE says if result is not a NaN, the sign is "-" if and only if |
| operands have opposite signs -- but flush -0 to 0 later if not IEEE. */ |
| sign = eisneg(a) ^ eisneg(b); |
| |
| #ifdef NANS |
| /* Return any NaN input. */ |
| 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, sign); |
| return; |
| } |
| #endif |
| /* Infinity over anything else is infinity. */ |
| #ifdef INFINITY |
| if (eisinf (b)) |
| { |
| einfin (c); |
| goto divsign; |
| } |
| /* Anything else over infinity is zero. */ |
| if (eisinf (a)) |
| { |
| eclear (c); |
| goto divsign; |
| } |
| #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); |
| goto divsign; |
| } |
| 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; |
| } |
| } |
| /* Divide by zero is not an invalid operation. |
| It is a divide-by-zero operation! */ |
| einfin (c); |
| mtherr ("ediv", SING); |
| goto divsign; |
| } |
| dnzro2: |
| |
| i = edivm (ai, bi); |
| /* calculate exponent */ |
| lt = ltb - lta + EXONE; |
| emdnorm (bi, i, 0, lt, 64); |
| emovo (bi, c); |
| |
| divsign: |
| |
| if (sign |
| #ifndef IEEE |
| && (ecmp (c, ezero) != 0) |
| #endif |
| ) |
| *(c+(NE-1)) |= 0x8000; |
| else |
| *(c+(NE-1)) &= ~0x8000; |
| } |
| |
| /* Multiply e-types A and B, return e-type product C. */ |
| |
| static void |
| emul (a, b, c) |
| unsigned EMUSHORT *a, *b, *c; |
| { |
| unsigned EMUSHORT ai[NI], bi[NI]; |
| int i, j, sign; |
| EMULONG lt, lta, ltb; |
| |
| /* IEEE says if result is not a NaN, the sign is "-" if and only if |
| operands have opposite signs -- but flush -0 to 0 later if not IEEE. */ |
| sign = eisneg(a) ^ eisneg(b); |
| |
| #ifdef NANS |
| /* NaN times anything is the same NaN. */ |
| 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, sign); |
| return; |
| } |
| #endif |
| /* Infinity times anything else is infinity. */ |
| #ifdef INFINITY |
| if (eisinf (a) || eisinf (b)) |
| { |
| einfin (c); |
| goto mulsign; |
| } |
| #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); |
| goto mulsign; |
| } |
| mnzer1: |
| |
| if (bi[E] == 0) |
| { |
| for (i = 1; i < NI - 1; i++) |
| { |
| if (bi[i] != 0) |
| { |
| ltb -= enormlz (bi); |
| goto mnzer2; |
| } |
| } |
| eclear (c); |
| goto mulsign; |
| } |
| mnzer2: |
| |
| /* Multiply significands */ |
| j = emulm (ai, bi); |
| /* calculate exponent */ |
| lt = lta + ltb - (EXONE - 1); |
| emdnorm (bi, j, 0, lt, 64); |
| emovo (bi, c); |
| |
| mulsign: |
| |
| if (sign |
| #ifndef IEEE |
| && (ecmp (c, ezero) != 0) |
| #endif |
| ) |
| *(c+(NE-1)) |= 0x8000; |
| else |
| *(c+(NE-1)) &= ~0x8000; |
| } |
| |
| /* Convert double precision PE to e-type Y. */ |
| |
| static void |
| e53toe (pe, y) |
| unsigned EMUSHORT *pe, *y; |
| { |
| #ifdef DEC |
| |
| dectoe (pe, y); |
| |
| #else |
| #ifdef IBM |
| |
| ibmtoe (pe, y, DFmode); |
| |
| #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); |
| if (! REAL_WORDS_BIG_ENDIAN) |
| e += 3; |
| 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 |
| if (! REAL_WORDS_BIG_ENDIAN) |
| { |
| if (((pe[3] & 0xf) != 0) || (pe[2] != 0) |
| || (pe[1] != 0) || (pe[0] != 0)) |
| { |
| enan (y, yy[0] != 0); |
| return; |
| } |
| } |
| else |
| { |
| if (((pe[0] & 0xf) != 0) || (pe[1] != 0) |
| || (pe[2] != 0) || (pe[3] != 0)) |
| { |
| enan (y, yy[0] != 0); |
| return; |
| } |
| } |
| #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 IEEE |
| if (! REAL_WORDS_BIG_ENDIAN) |
| { |
| *p++ = *(--e); |
| *p++ = *(--e); |
| *p++ = *(--e); |
| } |
| else |
| { |
| ++e; |
| *p++ = *e++; |
| *p++ = *e++; |
| *p++ = *e++; |
| } |
| #endif |
| eshift (yy, -5); |
| if (denorm) |
| { /* 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 IBM */ |
| #endif /* not DEC */ |
| } |
| |
| /* Convert double extended precision float PE to e type Y. */ |
| |
| static 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; |
| /* This precision is not ordinarily supported on DEC or IBM. */ |
| #ifdef DEC |
| for (i = 0; i < 5; i++) |
| *p++ = *e++; |
| #endif |
| #ifdef IBM |
| p = &yy[0] + (NE - 1); |
| *p-- = *e++; |
| ++e; |
| for (i = 0; i < 5; i++) |
| *p-- = *e++; |
| #endif |
| #ifdef IEEE |
| if (! REAL_WORDS_BIG_ENDIAN) |
| { |
| for (i = 0; i < 5; i++) |
| *p++ = *e++; |
| |
| /* For denormal long double Intel format, shift significand up one |
| -- but only if the top significand bit is zero. A top bit of 1 |
| is "pseudodenormal" when the exponent is zero. */ |
| if((yy[NE-1] & 0x7fff) == 0 && (yy[NE-2] & 0x8000) == 0) |
| { |
| unsigned EMUSHORT temp[NI]; |
| |
| emovi(yy, temp); |
| eshup1(temp); |
| emovo(temp,y); |
| return; |
| } |
| } |
| else |
| { |
| p = &yy[0] + (NE - 1); |
| #ifdef ARM_EXTENDED_IEEE_FORMAT |
| /* For ARMs, the exponent is in the lowest 15 bits of the word. */ |
| *p-- = (e[0] & 0x8000) | (e[1] & 0x7ffff); |
| e += 2; |
| #else |
| *p-- = *e++; |
| ++e; |
| #endif |
| for (i = 0; i < 4; i++) |
| *p-- = *e++; |
| } |
| #endif |
| #ifdef INFINITY |
| /* Point to the exponent field and check max exponent cases. */ |
| p = &yy[NE - 1]; |
| if ((*p & 0x7fff) == 0x7fff) |
| { |
| #ifdef NANS |
| if (! REAL_WORDS_BIG_ENDIAN) |
| { |
| for (i = 0; i < 4; i++) |
| { |
| if ((i != 3 && pe[i] != 0) |
| /* Anything but 0x8000 here, including 0, is a NaN. */ |
| || (i == 3 && pe[i] != 0x8000)) |
| { |
| enan (y, (*p & 0x8000) != 0); |
| return; |
| } |
| } |
| } |
| else |
| { |
| #ifdef ARM_EXTENDED_IEEE_FORMAT |
| for (i = 2; i <= 5; i++) |
| { |
| if (pe[i] != 0) |
| { |
| enan (y, (*p & 0x8000) != 0); |
| return; |
| } |
| } |
| #else /* not ARM */ |
| /* In Motorola extended precision format, the most significant |
| bit of an infinity mantissa could be either 1 or 0. It is |
| the lower order bits that tell whether the value is a NaN. */ |
| if ((pe[2] & 0x7fff) != 0) |
| goto bigend_nan; |
| |
| for (i = 3; i <= 5; i++) |
| { |
| if (pe[i] != 0) |
| { |
| bigend_nan: |
| enan (y, (*p & 0x8000) != 0); |
| return; |
| } |
| } |
| #endif /* not ARM */ |
| } |
| #endif /* NANS */ |
| eclear (y); |
| einfin (y); |
| if (*p & 0x8000) |
| eneg (y); |
| return; |
| } |
| #endif /* INFINITY */ |
| p = yy; |
| q = y; |
| for (i = 0; i < NE; i++) |
| *q++ = *p++; |
| } |
| |
| /* Convert 128-bit long double precision float PE to e type Y. */ |
| |
| static void |
| e113toe (pe, y) |
| unsigned EMUSHORT *pe, *y; |
| { |
| register unsigned EMUSHORT r; |
| unsigned EMUSHORT *e, *p; |
| unsigned EMUSHORT yy[NI]; |
| int denorm, i; |
| |
| e = pe; |
| denorm = 0; |
| ecleaz (yy); |
| #ifdef IEEE |
| if (! REAL_WORDS_BIG_ENDIAN) |
| e += 7; |
| #endif |
| r = *e; |
| yy[0] = 0; |
| if (r & 0x8000) |
| yy[0] = 0xffff; |
| r &= 0x7fff; |
| #ifdef INFINITY |
| if (r == 0x7fff) |
| { |
| #ifdef NANS |
| if (! REAL_WORDS_BIG_ENDIAN) |
| { |
| for (i = 0; i < 7; i++) |
| { |
| if (pe[i] != 0) |
| { |
| enan (y, yy[0] != 0); |
| return; |
| } |
| } |
| } |
| else |
| { |
| for (i = 1; i < 8; i++) |
| { |
| if (pe[i] != 0) |
| { |
| enan (y, yy[0] != 0); |
| return; |
| } |
| } |
| } |
| #endif /* NANS */ |
| eclear (y); |
| einfin (y); |
| if (yy[0]) |
| eneg (y); |
| return; |
| } |
| #endif /* INFINITY */ |
| yy[E] = r; |
| p = &yy[M + 1]; |
| #ifdef IEEE |
| if (! REAL_WORDS_BIG_ENDIAN) |
| { |
| for (i = 0; i < 7; i++) |
| *p++ = *(--e |