real.c revision 52284
1184610Salfred/* real.c - implementation of REAL_ARITHMETIC, REAL_VALUE_ATOF, 2184610Salfred and support for XFmode IEEE extended real floating point arithmetic. 3184610Salfred Copyright (C) 1993, 94-98, 1999 Free Software Foundation, Inc. 4184610Salfred Contributed by Stephen L. Moshier (moshier@world.std.com). 5184610Salfred 6184610SalfredThis file is part of GNU CC. 7184610Salfred 8184610SalfredGNU CC is free software; you can redistribute it and/or modify 9184610Salfredit under the terms of the GNU General Public License as published by 10184610Salfredthe Free Software Foundation; either version 2, or (at your option) 11184610Salfredany later version. 12184610Salfred 13184610SalfredGNU CC is distributed in the hope that it will be useful, 14184610Salfredbut WITHOUT ANY WARRANTY; without even the implied warranty of 15184610SalfredMERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 16184610SalfredGNU General Public License for more details. 17184610Salfred 18184610SalfredYou should have received a copy of the GNU General Public License 19184610Salfredalong with GNU CC; see the file COPYING. If not, write to 20184610Salfredthe Free Software Foundation, 59 Temple Place - Suite 330, 21184610SalfredBoston, MA 02111-1307, USA. */ 22184610Salfred 23184610Salfred#include "config.h" 24184610Salfred#include "system.h" 25184610Salfred#include "tree.h" 26184610Salfred#include "toplev.h" 27184610Salfred 28184610Salfred/* To enable support of XFmode extended real floating point, define 29184610SalfredLONG_DOUBLE_TYPE_SIZE 96 in the tm.h file (m68k.h or i386.h). 30184610Salfred 31184610SalfredTo support cross compilation between IEEE, VAX and IBM floating 32184610Salfredpoint formats, define REAL_ARITHMETIC in the tm.h file. 33184610Salfred 34184610SalfredIn either case the machine files (tm.h) must not contain any code 35184610Salfredthat tries to use host floating point arithmetic to convert 36184610SalfredREAL_VALUE_TYPEs from `double' to `float', pass them to fprintf, 37184610Salfredetc. In cross-compile situations a REAL_VALUE_TYPE may not 38184610Salfredbe intelligible to the host computer's native arithmetic. 39184610Salfred 40184610SalfredThe emulator defaults to the host's floating point format so that 41264149Sianits decimal conversion functions can be used if desired (see 42264149Sianreal.h). 43264149Sian 44264149SianThe first part of this file interfaces gcc to a floating point 45264149Sianarithmetic suite that was not written with gcc in mind. Avoid 46264149Sianchanging the low-level arithmetic routines unless you have suitable 47264149Siantest programs available. A special version of the PARANOIA floating 48264149Sianpoint arithmetic tester, modified for this purpose, can be found on 49184610Salfredusc.edu: /pub/C-numanal/ieeetest.zoo. Other tests, and libraries of 50184610SalfredXFmode and TFmode transcendental functions, can be obtained by ftp from 51194677Sthompsanetlib.att.com: netlib/cephes. */ 52194677Sthompsa 53194677Sthompsa/* Type of computer arithmetic. 54194677Sthompsa Only one of DEC, IBM, IEEE, C4X, or UNK should get defined. 55194677Sthompsa 56194677Sthompsa `IEEE', when REAL_WORDS_BIG_ENDIAN is non-zero, refers generically 57194677Sthompsa to big-endian IEEE floating-point data structure. This definition 58194677Sthompsa should work in SFmode `float' type and DFmode `double' type on 59194677Sthompsa virtually all big-endian IEEE machines. If LONG_DOUBLE_TYPE_SIZE 60194677Sthompsa has been defined to be 96, then IEEE also invokes the particular 61194677Sthompsa XFmode (`long double' type) data structure used by the Motorola 62194677Sthompsa 680x0 series processors. 63194677Sthompsa 64194677Sthompsa `IEEE', when REAL_WORDS_BIG_ENDIAN is zero, refers generally to 65194677Sthompsa little-endian IEEE machines. In this case, if LONG_DOUBLE_TYPE_SIZE 66194677Sthompsa has been defined to be 96, then IEEE also invokes the particular 67194677Sthompsa XFmode `long double' data structure used by the Intel 80x86 series 68194677Sthompsa processors. 69194677Sthompsa 70194677Sthompsa `DEC' refers specifically to the Digital Equipment Corp PDP-11 71194677Sthompsa and VAX floating point data structure. This model currently 72194677Sthompsa supports no type wider than DFmode. 73264031Sian 74264149Sian `IBM' refers specifically to the IBM System/370 and compatible 75188746Sthompsa floating point data structure. This model currently supports 76184610Salfred no type wider than DFmode. The IBM conversions were contributed by 77184610Salfred frank@atom.ansto.gov.au (Frank Crawford). 78188942Sthompsa 79188942Sthompsa `C4X' refers specifically to the floating point format used on 80184610Salfred Texas Instruments TMS320C3x and TMS320C4x digital signal 81188942Sthompsa processors. This supports QFmode (32-bit float, double) and HFmode 82188942Sthompsa (40-bit long double) where BITS_PER_BYTE is 32. Unlike IEEE 83264149Sian floats, C4x floats are not rounded to be even. The C4x conversions 84184610Salfred were contributed by m.hayes@elec.canterbury.ac.nz (Michael Hayes) and 85264923Sian Haj.Ten.Brugge@net.HCC.nl (Herman ten Brugge). 86264923Sian 87207077Sthompsa If LONG_DOUBLE_TYPE_SIZE = 64 (the default, unless tm.h defines it) 88184610Salfred then `long double' and `double' are both implemented, but they 89276701Shselasky both mean DFmode. In this case, the software floating-point 90184610Salfred support available here is activated by writing 91184610Salfred #define REAL_ARITHMETIC 92184610Salfred in tm.h. 93184610Salfred 94184610Salfred The case LONG_DOUBLE_TYPE_SIZE = 128 activates TFmode support 95264031Sian and may deactivate XFmode since `long double' is used to refer 96264031Sian to both modes. 97264031Sian 98264031Sian The macros FLOAT_WORDS_BIG_ENDIAN, HOST_FLOAT_WORDS_BIG_ENDIAN, 99264031Sian contributed by Richard Earnshaw <Richard.Earnshaw@cl.cam.ac.uk>, 100264031Sian separate the floating point unit's endian-ness from that of 101264031Sian the integer addressing. This permits one to define a big-endian 102264031Sian FPU on a little-endian machine (e.g., ARM). An extension to 103264031Sian BYTES_BIG_ENDIAN may be required for some machines in the future. 104264031Sian These optional macros may be defined in tm.h. In real.h, they 105264031Sian default to WORDS_BIG_ENDIAN, etc., so there is no need to define 106264031Sian them for any normal host or target machine on which the floats 107264031Sian and the integers have the same endian-ness. */ 108264031Sian 109264031Sian 110264031Sian/* The following converts gcc macros into the ones used by this file. */ 111264031Sian 112264031Sian/* REAL_ARITHMETIC defined means that macros in real.h are 113264031Sian defined to call emulator functions. */ 114264031Sian#ifdef REAL_ARITHMETIC 115264031Sian 116264031Sian#if TARGET_FLOAT_FORMAT == VAX_FLOAT_FORMAT 117264031Sian/* PDP-11, Pro350, VAX: */ 118264031Sian#define DEC 1 119264031Sian#else /* it's not VAX */ 120264031Sian#if TARGET_FLOAT_FORMAT == IBM_FLOAT_FORMAT 121264031Sian/* IBM System/370 style */ 122264031Sian#define IBM 1 123184610Salfred#else /* it's also not an IBM */ 124187259Sthompsa#if TARGET_FLOAT_FORMAT == C4X_FLOAT_FORMAT 125187259Sthompsa/* TMS320C3x/C4x style */ 126187259Sthompsa#define C4X 1 127188413Sthompsa#else /* it's also not a C4X */ 128187259Sthompsa#if TARGET_FLOAT_FORMAT == IEEE_FLOAT_FORMAT 129187259Sthompsa#define IEEE 130264010Sian#else /* it's not IEEE either */ 131264010Sian/* UNKnown arithmetic. We don't support this and can't go on. */ 132264010Sianunknown arithmetic type 133264010Sian#define UNK 1 134264010Sian#endif /* not IEEE */ 135264010Sian#endif /* not C4X */ 136264010Sian#endif /* not IBM */ 137264010Sian#endif /* not VAX */ 138264010Sian 139264010Sian#define REAL_WORDS_BIG_ENDIAN FLOAT_WORDS_BIG_ENDIAN 140264010Sian 141264010Sian#else 142264010Sian/* REAL_ARITHMETIC not defined means that the *host's* data 143264010Sian structure will be used. It may differ by endian-ness from the 144264010Sian target machine's structure and will get its ends swapped 145184610Salfred accordingly (but not here). Probably only the decimal <-> binary 146192984Sthompsa functions in this file will actually be used in this case. */ 147192984Sthompsa 148184610Salfred#if HOST_FLOAT_FORMAT == VAX_FLOAT_FORMAT 149192984Sthompsa#define DEC 1 150192984Sthompsa#else /* it's not VAX */ 151184610Salfred#if HOST_FLOAT_FORMAT == IBM_FLOAT_FORMAT 152189265Sthompsa/* IBM System/370 style */ 153184610Salfred#define IBM 1 154184610Salfred#else /* it's also not an IBM */ 155184610Salfred#if HOST_FLOAT_FORMAT == IEEE_FLOAT_FORMAT 156184610Salfred#define IEEE 157264149Sian#else /* it's not IEEE either */ 158184610Salfredunknown arithmetic type 159264010Sian#define UNK 1 160264010Sian#endif /* not IEEE */ 161184610Salfred#endif /* not IBM */ 162184610Salfred#endif /* not VAX */ 163184610Salfred 164286385Sian#define REAL_WORDS_BIG_ENDIAN HOST_FLOAT_WORDS_BIG_ENDIAN 165184610Salfred 166184610Salfred#endif /* REAL_ARITHMETIC not defined */ 167184610Salfred 168264010Sian/* Define INFINITY for support of infinity. 169264010Sian Define NANS for support of Not-a-Number's (NaN's). */ 170184610Salfred#if !defined(DEC) && !defined(IBM) && !defined(C4X) 171184610Salfred#define INFINITY 172184610Salfred#define NANS 173184610Salfred#endif 174184610Salfred 175184610Salfred/* Support of NaNs requires support of infinity. */ 176184610Salfred#ifdef NANS 177184610Salfred#ifndef INFINITY 178184610Salfred#define INFINITY 179184610Salfred#endif 180184610Salfred#endif 181239299Shselasky 182184610Salfred/* Find a host integer type that is at least 16 bits wide, 183193045Sthompsa and another type at least twice whatever that size is. */ 184193045Sthompsa 185184610Salfred#if HOST_BITS_PER_CHAR >= 16 186239180Shselasky#define EMUSHORT char 187192984Sthompsa#define EMUSHORT_SIZE HOST_BITS_PER_CHAR 188264149Sian#define EMULONG_SIZE (2 * HOST_BITS_PER_CHAR) 189192984Sthompsa#else 190192984Sthompsa#if HOST_BITS_PER_SHORT >= 16 191192984Sthompsa#define EMUSHORT short 192264010Sian#define EMUSHORT_SIZE HOST_BITS_PER_SHORT 193264010Sian#define EMULONG_SIZE (2 * HOST_BITS_PER_SHORT) 194192984Sthompsa#else 195192984Sthompsa#if HOST_BITS_PER_INT >= 16 196192984Sthompsa#define EMUSHORT int 197185948Sthompsa#define EMUSHORT_SIZE HOST_BITS_PER_INT 198264149Sian#define EMULONG_SIZE (2 * HOST_BITS_PER_INT) 199264149Sian#else 200286385Sian#if HOST_BITS_PER_LONG >= 16 201264149Sian#define EMUSHORT long 202264149Sian#define EMUSHORT_SIZE HOST_BITS_PER_LONG 203264149Sian#define EMULONG_SIZE (2 * HOST_BITS_PER_LONG) 204264149Sian#else 205264149Sian/* You will have to modify this program to have a smaller unit size. */ 206264149Sian#define EMU_NON_COMPILE 207192984Sthompsa#endif 208192984Sthompsa#endif 209192984Sthompsa#endif 210192984Sthompsa#endif 211197570Sthompsa 212184610Salfred#if HOST_BITS_PER_SHORT >= EMULONG_SIZE 213192984Sthompsa#define EMULONG short 214184610Salfred#else 215187259Sthompsa#if HOST_BITS_PER_INT >= EMULONG_SIZE 216184610Salfred#define EMULONG int 217184610Salfred#else 218184610Salfred#if HOST_BITS_PER_LONG >= EMULONG_SIZE 219190734Sthompsa#define EMULONG long 220200308Sthompsa#else 221190734Sthompsa#if HOST_BITS_PER_LONGLONG >= EMULONG_SIZE 222184610Salfred#define EMULONG long long int 223184610Salfred#else 224187259Sthompsa/* You will have to modify this program to have a smaller unit size. */ 225184610Salfred#define EMU_NON_COMPILE 226184610Salfred#endif 227184610Salfred#endif 228264031Sian#endif 229190734Sthompsa#endif 230190734Sthompsa 231184610Salfred 232184610Salfred/* The host interface doesn't work if no 16-bit size exists. */ 233184610Salfred#if EMUSHORT_SIZE != 16 234192984Sthompsa#define EMU_NON_COMPILE 235194228Sthompsa#endif 236194228Sthompsa 237194228Sthompsa/* OK to continue compilation. */ 238194228Sthompsa#ifndef EMU_NON_COMPILE 239194228Sthompsa 240194228Sthompsa/* Construct macros to translate between REAL_VALUE_TYPE and e type. 241264149Sian In GET_REAL and PUT_REAL, r and e are pointers. 242194228Sthompsa A REAL_VALUE_TYPE is guaranteed to occupy contiguous locations 243264149Sian in memory, with no holes. */ 244194228Sthompsa 245194228Sthompsa#if LONG_DOUBLE_TYPE_SIZE == 96 246194228Sthompsa/* Number of 16 bit words in external e type format */ 247194228Sthompsa#define NE 6 248197570Sthompsa#define MAXDECEXP 4932 249239180Shselasky#define MINDECEXP -4956 250184610Salfred#define GET_REAL(r,e) bcopy ((char *) r, (char *) e, 2*NE) 251184610Salfred#define PUT_REAL(e,r) \ 252184610Salfreddo { \ 253184610Salfred if (2*NE < sizeof(*r)) \ 254184610Salfred bzero((char *)r, sizeof(*r)); \ 255184610Salfred bcopy ((char *) e, (char *) r, 2*NE); \ 256184610Salfred} while (0) 257237102Smarius#else /* no XFmode */ 258184610Salfred#if LONG_DOUBLE_TYPE_SIZE == 128 259184610Salfred#define NE 10 260184610Salfred#define MAXDECEXP 4932 261184610Salfred#define MINDECEXP -4977 262184610Salfred#define GET_REAL(r,e) bcopy ((char *) r, (char *) e, 2*NE) 263184610Salfred#define PUT_REAL(e,r) bcopy ((char *) e, (char *) r, 2*NE) 264184610Salfred#else 265184610Salfred#define NE 6 266184610Salfred#define MAXDECEXP 4932 267184610Salfred#define MINDECEXP -4956 268237102Smarius#ifdef REAL_ARITHMETIC 269237102Smarius/* Emulator uses target format internally 270237102Smarius but host stores it in host endian-ness. */ 271264923Sian 272264923Sian#define GET_REAL(r,e) \ 273264923Siando { \ 274264923Sian if (HOST_FLOAT_WORDS_BIG_ENDIAN == REAL_WORDS_BIG_ENDIAN) \ 275264923Sian e53toe ((unsigned EMUSHORT *) (r), (e)); \ 276264923Sian else \ 277264923Sian { \ 278264923Sian unsigned EMUSHORT w[4]; \ 279264923Sian bcopy (((EMUSHORT *) r), &w[3], sizeof (EMUSHORT)); \ 280264923Sian bcopy (((EMUSHORT *) r) + 1, &w[2], sizeof (EMUSHORT)); \ 281264923Sian bcopy (((EMUSHORT *) r) + 2, &w[1], sizeof (EMUSHORT)); \ 282264923Sian bcopy (((EMUSHORT *) r) + 3, &w[0], sizeof (EMUSHORT)); \ 283264923Sian e53toe (w, (e)); \ 284264923Sian } \ 285264923Sian } while (0) 286264923Sian 287264923Sian#define PUT_REAL(e,r) \ 288264923Siando { \ 289264923Sian if (HOST_FLOAT_WORDS_BIG_ENDIAN == REAL_WORDS_BIG_ENDIAN) \ 290264923Sian etoe53 ((e), (unsigned EMUSHORT *) (r)); \ 291264923Sian else \ 292264923Sian { \ 293264923Sian unsigned EMUSHORT w[4]; \ 294264923Sian etoe53 ((e), w); \ 295264923Sian bcopy (&w[3], ((EMUSHORT *) r), sizeof (EMUSHORT)); \ 296264923Sian bcopy (&w[2], ((EMUSHORT *) r) + 1, sizeof (EMUSHORT)); \ 297273181Sjoerg bcopy (&w[1], ((EMUSHORT *) r) + 2, sizeof (EMUSHORT)); \ 298264923Sian bcopy (&w[0], ((EMUSHORT *) r) + 3, sizeof (EMUSHORT)); \ 299264923Sian } \ 300264923Sian } while (0) 301264923Sian 302264923Sian#else /* not REAL_ARITHMETIC */ 303264923Sian 304264923Sian/* emulator uses host format */ 305264923Sian#define GET_REAL(r,e) e53toe ((unsigned EMUSHORT *) (r), (e)) 306264923Sian#define PUT_REAL(e,r) etoe53 ((e), (unsigned EMUSHORT *) (r)) 307264923Sian 308264923Sian#endif /* not REAL_ARITHMETIC */ 309264923Sian#endif /* not TFmode */ 310264923Sian#endif /* not XFmode */ 311264923Sian 312264923Sian 313264923Sian/* Number of 16 bit words in internal format */ 314264923Sian#define NI (NE+3) 315264923Sian 316264923Sian/* Array offset to exponent */ 317264923Sian#define E 1 318264923Sian 319264923Sian/* Array offset to high guard word */ 320264923Sian#define M 2 321264923Sian 322264923Sian/* Number of bits of precision */ 323264923Sian#define NBITS ((NI-4)*16) 324264923Sian 325264923Sian/* Maximum number of decimal digits in ASCII conversion 326264923Sian * = NBITS*log10(2) 327264923Sian */ 328264923Sian#define NDEC (NBITS*8/27) 329264923Sian 330264923Sian/* The exponent of 1.0 */ 331264923Sian#define EXONE (0x3fff) 332264923Sian 333264923Sianextern int extra_warnings; 334264923Sianextern unsigned EMUSHORT ezero[], ehalf[], eone[], etwo[]; 335264923Sianextern unsigned EMUSHORT elog2[], esqrt2[]; 336264923Sian 337264923Sianstatic void endian PROTO((unsigned EMUSHORT *, long *, 338264923Sian enum machine_mode)); 339264923Sianstatic void eclear PROTO((unsigned EMUSHORT *)); 340264923Sianstatic void emov PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 341264923Sian#if 0 342264923Sianstatic void eabs PROTO((unsigned EMUSHORT *)); 343264923Sian#endif 344264923Sianstatic void eneg PROTO((unsigned EMUSHORT *)); 345264923Sianstatic int eisneg PROTO((unsigned EMUSHORT *)); 346264923Sianstatic int eisinf PROTO((unsigned EMUSHORT *)); 347264923Sianstatic int eisnan PROTO((unsigned EMUSHORT *)); 348264923Sianstatic void einfin PROTO((unsigned EMUSHORT *)); 349264923Sianstatic void enan PROTO((unsigned EMUSHORT *, int)); 350264923Sianstatic void emovi PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 351264923Sianstatic void emovo PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 352264923Sianstatic void ecleaz PROTO((unsigned EMUSHORT *)); 353264923Sianstatic void ecleazs PROTO((unsigned EMUSHORT *)); 354264923Sianstatic void emovz PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 355264923Sianstatic void einan PROTO((unsigned EMUSHORT *)); 356264923Sianstatic int eiisnan PROTO((unsigned EMUSHORT *)); 357264923Sianstatic int eiisneg PROTO((unsigned EMUSHORT *)); 358264923Sian#if 0 359264923Sianstatic void eiinfin PROTO((unsigned EMUSHORT *)); 360264923Sian#endif 361264923Sianstatic int eiisinf PROTO((unsigned EMUSHORT *)); 362264923Sianstatic int ecmpm PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 363264923Sianstatic void eshdn1 PROTO((unsigned EMUSHORT *)); 364264923Sianstatic void eshup1 PROTO((unsigned EMUSHORT *)); 365264923Sianstatic void eshdn8 PROTO((unsigned EMUSHORT *)); 366264923Sianstatic void eshup8 PROTO((unsigned EMUSHORT *)); 367264923Sianstatic void eshup6 PROTO((unsigned EMUSHORT *)); 368264923Sianstatic void eshdn6 PROTO((unsigned EMUSHORT *)); 369264923Sianstatic void eaddm PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 370264923Sianstatic void esubm PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 371264923Sianstatic void m16m PROTO((unsigned int, unsigned short *, 372264923Sian unsigned short *)); 373264923Sianstatic int edivm PROTO((unsigned short *, unsigned short *)); 374264923Sianstatic int emulm PROTO((unsigned short *, unsigned short *)); 375264923Sianstatic void emdnorm PROTO((unsigned EMUSHORT *, int, int, EMULONG, int)); 376264923Sianstatic void esub PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *, 377264923Sian unsigned EMUSHORT *)); 378264923Sianstatic void eadd PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *, 379264923Sian unsigned EMUSHORT *)); 380264923Sianstatic void eadd1 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *, 381264923Sian unsigned EMUSHORT *)); 382264923Sianstatic void ediv PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *, 383264923Sian unsigned EMUSHORT *)); 384264923Sianstatic void emul PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *, 385264923Sian unsigned EMUSHORT *)); 386264923Sianstatic void e53toe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 387264923Sianstatic void e64toe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 388264923Sianstatic void e113toe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 389264923Sianstatic void e24toe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 390264923Sianstatic void etoe113 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 391264923Sianstatic void toe113 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 392264923Sianstatic void etoe64 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 393264923Sianstatic void toe64 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 394264923Sianstatic void etoe53 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 395264923Sianstatic void toe53 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 396264923Sianstatic void etoe24 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 397264923Sianstatic void toe24 PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 398264923Sianstatic int ecmp PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 399264923Sian#if 0 400264923Sianstatic void eround PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 401264923Sian#endif 402264923Sianstatic void ltoe PROTO((HOST_WIDE_INT *, unsigned EMUSHORT *)); 403264923Sianstatic void ultoe PROTO((unsigned HOST_WIDE_INT *, unsigned EMUSHORT *)); 404264923Sianstatic void eifrac PROTO((unsigned EMUSHORT *, HOST_WIDE_INT *, 405264923Sian unsigned EMUSHORT *)); 406264923Sianstatic void euifrac PROTO((unsigned EMUSHORT *, unsigned HOST_WIDE_INT *, 407264923Sian unsigned EMUSHORT *)); 408264923Sianstatic int eshift PROTO((unsigned EMUSHORT *, int)); 409264923Sianstatic int enormlz PROTO((unsigned EMUSHORT *)); 410264923Sian#if 0 411264923Sianstatic void e24toasc PROTO((unsigned EMUSHORT *, char *, int)); 412264923Sianstatic void e53toasc PROTO((unsigned EMUSHORT *, char *, int)); 413264923Sianstatic void e64toasc PROTO((unsigned EMUSHORT *, char *, int)); 414264923Sianstatic void e113toasc PROTO((unsigned EMUSHORT *, char *, int)); 415264923Sian#endif /* 0 */ 416264923Sianstatic void etoasc PROTO((unsigned EMUSHORT *, char *, int)); 417264923Sianstatic void asctoe24 PROTO((const char *, unsigned EMUSHORT *)); 418264923Sianstatic void asctoe53 PROTO((const char *, unsigned EMUSHORT *)); 419264923Sianstatic void asctoe64 PROTO((const char *, unsigned EMUSHORT *)); 420264923Sianstatic void asctoe113 PROTO((const char *, unsigned EMUSHORT *)); 421264923Sianstatic void asctoe PROTO((const char *, unsigned EMUSHORT *)); 422264923Sianstatic void asctoeg PROTO((const char *, unsigned EMUSHORT *, int)); 423264923Sianstatic void efloor PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 424264923Sian#if 0 425264923Sianstatic void efrexp PROTO((unsigned EMUSHORT *, int *, 426264923Sian unsigned EMUSHORT *)); 427264923Sian#endif 428264923Sianstatic void eldexp PROTO((unsigned EMUSHORT *, int, unsigned EMUSHORT *)); 429264923Sian#if 0 430264923Sianstatic void eremain PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *, 431264923Sian unsigned EMUSHORT *)); 432264923Sian#endif 433264923Sianstatic void eiremain PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 434264923Sianstatic void mtherr PROTO((const char *, int)); 435264923Sian#ifdef DEC 436264923Sianstatic void dectoe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 437264923Sianstatic void etodec PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 438264923Sianstatic void todec PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 439264923Sian#endif 440264923Sian#ifdef IBM 441264923Sianstatic void ibmtoe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *, 442264923Sian enum machine_mode)); 443264923Sianstatic void etoibm PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *, 444264923Sian enum machine_mode)); 445264923Sianstatic void toibm PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *, 446264923Sian enum machine_mode)); 447264923Sian#endif 448264923Sian#ifdef C4X 449264923Sianstatic void c4xtoe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *, 450264923Sian enum machine_mode)); 451264923Sianstatic void etoc4x PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *, 452264923Sian enum machine_mode)); 453264923Sianstatic void toc4x PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *, 454264923Sian enum machine_mode)); 455264923Sian#endif 456264923Sianstatic void make_nan PROTO((unsigned EMUSHORT *, int, enum machine_mode)); 457264923Sian#if 0 458264923Sianstatic void uditoe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 459264923Sianstatic void ditoe PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 460264923Sianstatic void etoudi PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 461264923Sianstatic void etodi PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 462264923Sianstatic void esqrt PROTO((unsigned EMUSHORT *, unsigned EMUSHORT *)); 463264923Sian#endif 464264923Sian 465264923Sian/* Copy 32-bit numbers obtained from array containing 16-bit numbers, 466264923Sian swapping ends if required, into output array of longs. The 467264923Sian result is normally passed to fprintf by the ASM_OUTPUT_ macros. */ 468264923Sian 469264923Sianstatic void 470264923Sianendian (e, x, mode) 471264923Sian unsigned EMUSHORT e[]; 472264923Sian long x[]; 473264923Sian enum machine_mode mode; 474264923Sian{ 475264923Sian unsigned long th, t; 476264923Sian 477264923Sian if (REAL_WORDS_BIG_ENDIAN) 478264923Sian { 479264923Sian switch (mode) 480264923Sian { 481264923Sian case TFmode: 482264923Sian /* Swap halfwords in the fourth long. */ 483264923Sian th = (unsigned long) e[6] & 0xffff; 484264923Sian t = (unsigned long) e[7] & 0xffff; 485264923Sian t |= th << 16; 486264923Sian x[3] = (long) t; 487264923Sian 488264923Sian case XFmode: 489264923Sian /* Swap halfwords in the third long. */ 490264923Sian th = (unsigned long) e[4] & 0xffff; 491264923Sian t = (unsigned long) e[5] & 0xffff; 492264923Sian t |= th << 16; 493264923Sian x[2] = (long) t; 494264923Sian /* fall into the double case */ 495264923Sian 496282505Shselasky case DFmode: 497264923Sian /* Swap halfwords in the second word. */ 498264923Sian th = (unsigned long) e[2] & 0xffff; 499264923Sian t = (unsigned long) e[3] & 0xffff; 500264923Sian t |= th << 16; 501264923Sian x[1] = (long) t; 502264923Sian /* fall into the float case */ 503264923Sian 504264923Sian case SFmode: 505264923Sian case HFmode: 506264923Sian /* Swap halfwords in the first word. */ 507264923Sian th = (unsigned long) e[0] & 0xffff; 508264923Sian t = (unsigned long) e[1] & 0xffff; 509264923Sian t |= th << 16; 510264923Sian x[0] = (long) t; 511264923Sian break; 512264923Sian 513264923Sian default: 514264923Sian abort (); 515264923Sian } 516264923Sian } 517264923Sian else 518264923Sian { 519264923Sian /* Pack the output array without swapping. */ 520264923Sian 521264923Sian switch (mode) 522264923Sian { 523264923Sian case TFmode: 524264923Sian /* Pack the fourth long. */ 525264923Sian th = (unsigned long) e[7] & 0xffff; 526264923Sian t = (unsigned long) e[6] & 0xffff; 527264923Sian t |= th << 16; 528264923Sian x[3] = (long) t; 529264923Sian 530264923Sian case XFmode: 531264923Sian /* Pack the third long. 532264923Sian Each element of the input REAL_VALUE_TYPE array has 16 useful bits 533264923Sian in it. */ 534264923Sian th = (unsigned long) e[5] & 0xffff; 535264923Sian t = (unsigned long) e[4] & 0xffff; 536264923Sian t |= th << 16; 537264934Sian x[2] = (long) t; 538264934Sian /* fall into the double case */ 539264923Sian 540264923Sian case DFmode: 541264923Sian /* Pack the second long */ 542264923Sian th = (unsigned long) e[3] & 0xffff; 543264923Sian t = (unsigned long) e[2] & 0xffff; 544264923Sian t |= th << 16; 545264923Sian x[1] = (long) t; 546264923Sian /* fall into the float case */ 547264923Sian 548264923Sian case SFmode: 549264923Sian case HFmode: 550264923Sian /* Pack the first long */ 551264923Sian th = (unsigned long) e[1] & 0xffff; 552264923Sian t = (unsigned long) e[0] & 0xffff; 553264923Sian t |= th << 16; 554264923Sian x[0] = (long) t; 555264923Sian break; 556264923Sian 557264923Sian default: 558264923Sian abort (); 559264923Sian } 560264923Sian } 561264923Sian} 562264923Sian 563264923Sian 564264923Sian/* This is the implementation of the REAL_ARITHMETIC macro. */ 565264923Sian 566292366Sianvoid 567264923Sianearith (value, icode, r1, r2) 568264923Sian REAL_VALUE_TYPE *value; 569264923Sian int icode; 570264923Sian REAL_VALUE_TYPE *r1; 571264923Sian REAL_VALUE_TYPE *r2; 572264923Sian{ 573264923Sian unsigned EMUSHORT d1[NE], d2[NE], v[NE]; 574264923Sian enum tree_code code; 575264923Sian 576264923Sian GET_REAL (r1, d1); 577264923Sian GET_REAL (r2, d2); 578264923Sian#ifdef NANS 579264923Sian/* Return NaN input back to the caller. */ 580264923Sian if (eisnan (d1)) 581264923Sian { 582264923Sian PUT_REAL (d1, value); 583264923Sian return; 584264923Sian } 585264923Sian if (eisnan (d2)) 586264923Sian { 587264923Sian PUT_REAL (d2, value); 588264923Sian return; 589264923Sian } 590264923Sian#endif 591264923Sian code = (enum tree_code) icode; 592264923Sian switch (code) 593264923Sian { 594264923Sian case PLUS_EXPR: 595264923Sian eadd (d2, d1, v); 596264923Sian break; 597264923Sian 598264923Sian case MINUS_EXPR: 599264923Sian esub (d2, d1, v); /* d1 - d2 */ 600264923Sian break; 601264923Sian 602264923Sian case MULT_EXPR: 603264923Sian emul (d2, d1, v); 604264923Sian break; 605264923Sian 606264923Sian case RDIV_EXPR: 607264923Sian#ifndef REAL_INFINITY 608264923Sian if (ecmp (d2, ezero) == 0) 609264923Sian { 610264923Sian#ifdef NANS 611264923Sian enan (v, eisneg (d1) ^ eisneg (d2)); 612264923Sian break; 613264923Sian#else 614264923Sian abort (); 615264923Sian#endif 616264923Sian } 617264923Sian#endif 618264923Sian ediv (d2, d1, v); /* d1/d2 */ 619264923Sian break; 620264923Sian 621264923Sian case MIN_EXPR: /* min (d1,d2) */ 622264923Sian if (ecmp (d1, d2) < 0) 623264923Sian emov (d1, v); 624264923Sian else 625264923Sian emov (d2, v); 626264923Sian break; 627264923Sian 628264923Sian case MAX_EXPR: /* max (d1,d2) */ 629264923Sian if (ecmp (d1, d2) > 0) 630264923Sian emov (d1, v); 631264923Sian else 632264923Sian emov (d2, v); 633264923Sian break; 634264923Sian default: 635264923Sian emov (ezero, v); 636264923Sian break; 637264923Sian } 638264923SianPUT_REAL (v, value); 639264923Sian} 640264923Sian 641264923Sian 642264923Sian/* Truncate REAL_VALUE_TYPE toward zero to signed HOST_WIDE_INT. 643264923Sian implements REAL_VALUE_RNDZINT (x) (etrunci (x)). */ 644264923Sian 645264923SianREAL_VALUE_TYPE 646264923Sianetrunci (x) 647264923Sian REAL_VALUE_TYPE x; 648264923Sian{ 649264923Sian unsigned EMUSHORT f[NE], g[NE]; 650264923Sian REAL_VALUE_TYPE r; 651264923Sian HOST_WIDE_INT l; 652264923Sian 653264923Sian GET_REAL (&x, g); 654264923Sian#ifdef NANS 655264923Sian if (eisnan (g)) 656264923Sian return (x); 657264923Sian#endif 658264923Sian eifrac (g, &l, f); 659264923Sian ltoe (&l, g); 660264923Sian PUT_REAL (g, &r); 661264923Sian return (r); 662264923Sian} 663264923Sian 664264923Sian 665264923Sian/* Truncate REAL_VALUE_TYPE toward zero to unsigned HOST_WIDE_INT; 666264923Sian implements REAL_VALUE_UNSIGNED_RNDZINT (x) (etruncui (x)). */ 667264923Sian 668264923SianREAL_VALUE_TYPE 669264923Sianetruncui (x) 670264923Sian REAL_VALUE_TYPE x; 671264923Sian{ 672264923Sian unsigned EMUSHORT f[NE], g[NE]; 673264923Sian REAL_VALUE_TYPE r; 674264923Sian unsigned HOST_WIDE_INT l; 675264923Sian 676264923Sian GET_REAL (&x, g); 677264923Sian#ifdef NANS 678264923Sian if (eisnan (g)) 679264923Sian return (x); 680264923Sian#endif 681264923Sian euifrac (g, &l, f); 682264923Sian ultoe (&l, g); 683264923Sian PUT_REAL (g, &r); 684264923Sian return (r); 685264923Sian} 686264923Sian 687264923Sian 688264923Sian/* This is the REAL_VALUE_ATOF function. It converts a decimal or hexadecimal 689264923Sian string to binary, rounding off as indicated by the machine_mode argument. 690264923Sian Then it promotes the rounded value to REAL_VALUE_TYPE. */ 691264923Sian 692264923SianREAL_VALUE_TYPE 693264923Sianereal_atof (s, t) 694264923Sian const char *s; 695264923Sian enum machine_mode t; 696264923Sian{ 697264923Sian unsigned EMUSHORT tem[NE], e[NE]; 698264923Sian REAL_VALUE_TYPE r; 699264923Sian 700264923Sian switch (t) 701264923Sian { 702264923Sian#ifdef C4X 703264923Sian case QFmode: 704264923Sian case HFmode: 705264923Sian asctoe53 (s, tem); 706264923Sian e53toe (tem, e); 707264923Sian break; 708264923Sian#else 709264923Sian case HFmode: 710264923Sian#endif 711264923Sian 712264923Sian case SFmode: 713264923Sian asctoe24 (s, tem); 714264923Sian e24toe (tem, e); 715264923Sian break; 716264923Sian 717264923Sian case DFmode: 718264923Sian asctoe53 (s, tem); 719264923Sian e53toe (tem, e); 720264923Sian break; 721264923Sian 722264923Sian case XFmode: 723264923Sian asctoe64 (s, tem); 724264923Sian e64toe (tem, e); 725264923Sian break; 726264923Sian 727264923Sian case TFmode: 728264923Sian asctoe113 (s, tem); 729264923Sian e113toe (tem, e); 730264923Sian break; 731264923Sian 732264923Sian default: 733264923Sian asctoe (s, e); 734264923Sian } 735264923Sian PUT_REAL (e, &r); 736264923Sian return (r); 737264923Sian} 738264923Sian 739264923Sian 740264923Sian/* Expansion of REAL_NEGATE. */ 741264923Sian 742264923SianREAL_VALUE_TYPE 743264923Sianereal_negate (x) 744264923Sian REAL_VALUE_TYPE x; 745264923Sian{ 746264923Sian unsigned EMUSHORT e[NE]; 747264923Sian REAL_VALUE_TYPE r; 748264923Sian 749264923Sian GET_REAL (&x, e); 750264923Sian eneg (e); 751264923Sian PUT_REAL (e, &r); 752264923Sian return (r); 753264923Sian} 754264923Sian 755264923Sian 756264923Sian/* Round real toward zero to HOST_WIDE_INT; 757264923Sian implements REAL_VALUE_FIX (x). */ 758264923Sian 759264923SianHOST_WIDE_INT 760264923Sianefixi (x) 761264923Sian REAL_VALUE_TYPE x; 762264923Sian{ 763264923Sian unsigned EMUSHORT f[NE], g[NE]; 764264923Sian HOST_WIDE_INT l; 765264923Sian 766264923Sian GET_REAL (&x, f); 767264923Sian#ifdef NANS 768264923Sian if (eisnan (f)) 769264923Sian { 770264923Sian warning ("conversion from NaN to int"); 771264923Sian return (-1); 772264923Sian } 773264923Sian#endif 774264923Sian eifrac (f, &l, g); 775264923Sian return l; 776264923Sian} 777264923Sian 778264923Sian/* Round real toward zero to unsigned HOST_WIDE_INT 779264923Sian implements REAL_VALUE_UNSIGNED_FIX (x). 780264923Sian Negative input returns zero. */ 781264923Sian 782264923Sianunsigned HOST_WIDE_INT 783264923Sianefixui (x) 784264923Sian REAL_VALUE_TYPE x; 785264923Sian{ 786264923Sian unsigned EMUSHORT f[NE], g[NE]; 787264923Sian unsigned HOST_WIDE_INT l; 788264923Sian 789264923Sian GET_REAL (&x, f); 790264923Sian#ifdef NANS 791264923Sian if (eisnan (f)) 792264923Sian { 793264923Sian warning ("conversion from NaN to unsigned int"); 794264923Sian return (-1); 795264923Sian } 796264923Sian#endif 797264923Sian euifrac (f, &l, g); 798264923Sian return l; 799264923Sian} 800264923Sian 801264923Sian 802264923Sian/* REAL_VALUE_FROM_INT macro. */ 803264923Sian 804264923Sianvoid 805264923Sianereal_from_int (d, i, j, mode) 806264923Sian REAL_VALUE_TYPE *d; 807264923Sian HOST_WIDE_INT i, j; 808264923Sian enum machine_mode mode; 809264923Sian{ 810264923Sian unsigned EMUSHORT df[NE], dg[NE]; 811264923Sian HOST_WIDE_INT low, high; 812264923Sian int sign; 813264923Sian 814264923Sian if (GET_MODE_CLASS (mode) != MODE_FLOAT) 815264923Sian abort (); 816264923Sian sign = 0; 817264923Sian low = i; 818264923Sian if ((high = j) < 0) 819264923Sian { 820264923Sian sign = 1; 821264923Sian /* complement and add 1 */ 822264923Sian high = ~high; 823264923Sian if (low) 824264923Sian low = -low; 825264923Sian else 826264923Sian high += 1; 827264923Sian } 828264923Sian eldexp (eone, HOST_BITS_PER_WIDE_INT, df); 829264923Sian ultoe ((unsigned HOST_WIDE_INT *) &high, dg); 830264923Sian emul (dg, df, dg); 831264923Sian ultoe ((unsigned HOST_WIDE_INT *) &low, df); 832264923Sian eadd (df, dg, dg); 833264923Sian if (sign) 834264923Sian eneg (dg); 835264923Sian 836264923Sian /* A REAL_VALUE_TYPE may not be wide enough to hold the two HOST_WIDE_INTS. 837264923Sian Avoid double-rounding errors later by rounding off now from the 838264923Sian extra-wide internal format to the requested precision. */ 839264923Sian switch (GET_MODE_BITSIZE (mode)) 840264923Sian { 841264923Sian case 32: 842264923Sian etoe24 (dg, df); 843264923Sian e24toe (df, dg); 844264923Sian break; 845264923Sian 846264923Sian case 64: 847264923Sian etoe53 (dg, df); 848264923Sian e53toe (df, dg); 849264923Sian break; 850264923Sian 851264923Sian case 96: 852264923Sian etoe64 (dg, df); 853264923Sian e64toe (df, dg); 854264923Sian break; 855264923Sian 856264923Sian case 128: 857264923Sian etoe113 (dg, df); 858264923Sian e113toe (df, dg); 859264923Sian break; 860264923Sian 861264923Sian default: 862264923Sian abort (); 863264923Sian } 864264923Sian 865264923Sian PUT_REAL (dg, d); 866264923Sian} 867264923Sian 868264923Sian 869264923Sian/* REAL_VALUE_FROM_UNSIGNED_INT macro. */ 870264923Sian 871264923Sianvoid 872264923Sianereal_from_uint (d, i, j, mode) 873264923Sian REAL_VALUE_TYPE *d; 874264923Sian unsigned HOST_WIDE_INT i, j; 875264923Sian enum machine_mode mode; 876264923Sian{ 877264923Sian unsigned EMUSHORT df[NE], dg[NE]; 878264923Sian unsigned HOST_WIDE_INT low, high; 879264923Sian 880264923Sian if (GET_MODE_CLASS (mode) != MODE_FLOAT) 881264923Sian abort (); 882264923Sian low = i; 883264923Sian high = j; 884264923Sian eldexp (eone, HOST_BITS_PER_WIDE_INT, df); 885264923Sian ultoe (&high, dg); 886264923Sian emul (dg, df, dg); 887264923Sian ultoe (&low, df); 888264923Sian eadd (df, dg, dg); 889264923Sian 890264923Sian /* A REAL_VALUE_TYPE may not be wide enough to hold the two HOST_WIDE_INTS. 891264923Sian Avoid double-rounding errors later by rounding off now from the 892264923Sian extra-wide internal format to the requested precision. */ 893264923Sian switch (GET_MODE_BITSIZE (mode)) 894264923Sian { 895264923Sian case 32: 896264923Sian etoe24 (dg, df); 897264923Sian e24toe (df, dg); 898264923Sian break; 899264923Sian 900264923Sian case 64: 901264923Sian etoe53 (dg, df); 902264923Sian e53toe (df, dg); 903264923Sian break; 904264923Sian 905264923Sian case 96: 906264923Sian etoe64 (dg, df); 907264923Sian e64toe (df, dg); 908264923Sian break; 909201028Sthompsa 910184610Salfred case 128: 911184610Salfred etoe113 (dg, df); 912292080Simp e113toe (df, dg); 913292080Simp break; 914292080Simp 915292080Simp default: 916292080Simp abort (); 917292080Simp } 918264010Sian 919264923Sian PUT_REAL (dg, d); 920264923Sian} 921264923Sian 922264923Sian 923264923Sian/* REAL_VALUE_TO_INT macro. */ 924264923Sian 925264923Sianvoid 926264923Sianereal_to_int (low, high, rr) 927264923Sian HOST_WIDE_INT *low, *high; 928264923Sian REAL_VALUE_TYPE rr; 929264923Sian{ 930264923Sian unsigned EMUSHORT d[NE], df[NE], dg[NE], dh[NE]; 931264923Sian int s; 932264923Sian 933264923Sian GET_REAL (&rr, d); 934264923Sian#ifdef NANS 935264923Sian if (eisnan (d)) 936264923Sian { 937264923Sian warning ("conversion from NaN to int"); 938264923Sian *low = -1; 939264923Sian *high = -1; 940264923Sian return; 941264923Sian } 942264923Sian#endif 943264923Sian /* convert positive value */ 944264923Sian s = 0; 945264923Sian if (eisneg (d)) 946264923Sian { 947264923Sian eneg (d); 948264923Sian s = 1; 949264923Sian } 950264923Sian eldexp (eone, HOST_BITS_PER_WIDE_INT, df); 951264923Sian ediv (df, d, dg); /* dg = d / 2^32 is the high word */ 952264923Sian euifrac (dg, (unsigned HOST_WIDE_INT *) high, dh); 953264923Sian emul (df, dh, dg); /* fractional part is the low word */ 954264923Sian euifrac (dg, (unsigned HOST_WIDE_INT *)low, dh); 955264923Sian if (s) 956264923Sian { 957264923Sian /* complement and add 1 */ 958264923Sian *high = ~(*high); 959264923Sian if (*low) 960264923Sian *low = -(*low); 961264923Sian else 962264923Sian *high += 1; 963264923Sian } 964264923Sian} 965264923Sian 966264923Sian 967264923Sian/* REAL_VALUE_LDEXP macro. */ 968264923Sian 969264923SianREAL_VALUE_TYPE 970264923Sianereal_ldexp (x, n) 971264923Sian REAL_VALUE_TYPE x; 972264923Sian int n; 973264923Sian{ 974264010Sian unsigned EMUSHORT e[NE], y[NE]; 975264010Sian REAL_VALUE_TYPE r; 976264010Sian 977264010Sian GET_REAL (&x, e); 978264010Sian#ifdef NANS 979264010Sian if (eisnan (e)) 980264010Sian return (x); 981264010Sian#endif 982264010Sian eldexp (e, n, y); 983264010Sian PUT_REAL (y, &r); 984264010Sian return (r); 985264010Sian} 986264010Sian 987264010Sian/* These routines are conditionally compiled because functions 988264010Sian of the same names may be defined in fold-const.c. */ 989264149Sian 990264149Sian#ifdef REAL_ARITHMETIC 991264010Sian 992264010Sian/* Check for infinity in a REAL_VALUE_TYPE. */ 993264010Sian 994264010Sianint 995264010Siantarget_isinf (x) 996264010Sian REAL_VALUE_TYPE x; 997264010Sian{ 998264010Sian unsigned EMUSHORT e[NE]; 999264010Sian 1000264010Sian#ifdef INFINITY 1001264010Sian GET_REAL (&x, e); 1002264010Sian return (eisinf (e)); 1003264010Sian#else 1004264010Sian return 0; 1005264010Sian#endif 1006264010Sian} 1007264010Sian 1008264010Sian/* Check whether a REAL_VALUE_TYPE item is a NaN. */ 1009264010Sian 1010264010Sianint 1011264010Siantarget_isnan (x) 1012264010Sian REAL_VALUE_TYPE x; 1013264010Sian{ 1014264010Sian unsigned EMUSHORT e[NE]; 1015264010Sian 1016264010Sian#ifdef NANS 1017264010Sian GET_REAL (&x, e); 1018264010Sian return (eisnan (e)); 1019264010Sian#else 1020264010Sian return (0); 1021264010Sian#endif 1022264010Sian} 1023264010Sian 1024264010Sian 1025264010Sian/* Check for a negative REAL_VALUE_TYPE number. 1026264010Sian This just checks the sign bit, so that -0 counts as negative. */ 1027264010Sian 1028264010Sianint 1029264010Siantarget_negative (x) 1030264010Sian REAL_VALUE_TYPE x; 1031264010Sian{ 1032264010Sian return ereal_isneg (x); 1033264010Sian} 1034264010Sian 1035264010Sian/* Expansion of REAL_VALUE_TRUNCATE. 1036264010Sian The result is in floating point, rounded to nearest or even. */ 1037264010Sian 1038264010SianREAL_VALUE_TYPE 1039264010Sianreal_value_truncate (mode, arg) 1040264010Sian enum machine_mode mode; 1041269569Sn_hibma REAL_VALUE_TYPE arg; 1042264010Sian{ 1043264010Sian unsigned EMUSHORT e[NE], t[NE]; 1044264010Sian REAL_VALUE_TYPE r; 1045264010Sian 1046264010Sian GET_REAL (&arg, e); 1047264010Sian#ifdef NANS 1048264010Sian if (eisnan (e)) 1049184610Salfred return (arg); 1050184610Salfred#endif 1051184610Salfred eclear (t); 1052192984Sthompsa switch (mode) 1053237102Smarius { 1054184610Salfred case TFmode: 1055192499Sthompsa etoe113 (e, t); 1056184610Salfred e113toe (t, t); 1057184610Salfred break; 1058184610Salfred 1059184610Salfred case XFmode: 1060184610Salfred etoe64 (e, t); 1061184610Salfred e64toe (t, t); 1062237102Smarius break; 1063237102Smarius 1064237102Smarius case DFmode: 1065237102Smarius etoe53 (e, t); 1066237102Smarius e53toe (t, t); 1067237102Smarius break; 1068237236Smarius 1069237102Smarius case SFmode: 1070264923Sian#ifndef C4X 1071264923Sian case HFmode: 1072264923Sian#endif 1073264923Sian etoe24 (e, t); 1074264923Sian e24toe (t, t); 1075264923Sian break; 1076237236Smarius 1077237236Smarius#ifdef C4X 1078237236Smarius case HFmode: 1079237102Smarius case QFmode: 1080237102Smarius etoe53 (e, t); 1081184610Salfred e53toe (t, t); 1082184610Salfred break; 1083184610Salfred#endif 1084184610Salfred 1085184610Salfred case SImode: 1086192984Sthompsa r = etrunci (arg); 1087184610Salfred return (r); 1088184610Salfred 1089184610Salfred /* If an unsupported type was requested, presume that 1090264010Sian the machine files know something useful to do with 1091264010Sian the unmodified value. */ 1092184610Salfred 1093184610Salfred default: 1094184610Salfred return (arg); 1095286385Sian } 1096184610Salfred PUT_REAL (t, &r); 1097194228Sthompsa return (r); 1098189265Sthompsa} 1099239180Shselasky 1100184610Salfred/* Try to change R into its exact multiplicative inverse in machine mode 1101184610Salfred MODE. Return nonzero function value if successful. */ 1102264010Sian 1103184610Salfredint 1104194228Sthompsaexact_real_inverse (mode, r) 1105264010Sian enum machine_mode mode; 1106189265Sthompsa REAL_VALUE_TYPE *r; 1107184610Salfred{ 1108184610Salfred unsigned EMUSHORT e[NE], einv[NE]; 1109184610Salfred REAL_VALUE_TYPE rinv; 1110199816Sthompsa int i; 1111184610Salfred 1112184610Salfred GET_REAL (r, e); 1113184610Salfred 1114189265Sthompsa /* Test for input in range. Don't transform IEEE special values. */ 1115194677Sthompsa if (eisinf (e) || eisnan (e) || (ecmp (e, ezero) == 0)) 1116194677Sthompsa return 0; 1117189265Sthompsa 1118184610Salfred /* Test for a power of 2: all significand bits zero except the MSB. 1119184610Salfred We are assuming the target has binary (or hex) arithmetic. */ 1120184610Salfred if (e[NE - 2] != 0x8000) 1121184610Salfred return 0; 1122184610Salfred 1123184610Salfred for (i = 0; i < NE - 2; i++) 1124184610Salfred { 1125184610Salfred if (e[i] != 0) 1126194228Sthompsa return 0; 1127189265Sthompsa } 1128184610Salfred 1129184610Salfred /* Compute the inverse and truncate it to the required mode. */ 1130184610Salfred ediv (e, eone, einv); 1131214843Sn_hibma PUT_REAL (einv, &rinv); 1132214843Sn_hibma rinv = real_value_truncate (mode, rinv); 1133184610Salfred 1134184610Salfred#ifdef CHECK_FLOAT_VALUE 1135184610Salfred /* This check is not redundant. It may, for example, flush 1136184610Salfred a supposedly IEEE denormal value to zero. */ 1137184610Salfred i = 0; 1138184610Salfred if (CHECK_FLOAT_VALUE (mode, rinv, i)) 1139184610Salfred return 0; 1140184610Salfred#endif 1141184610Salfred GET_REAL (&rinv, einv); 1142184610Salfred 1143184610Salfred /* Check the bits again, because the truncation might have 1144184610Salfred generated an arbitrary saturation value on overflow. */ 1145214761Sn_hibma if (einv[NE - 2] != 0x8000) 1146194228Sthompsa return 0; 1147184610Salfred 1148239299Shselasky for (i = 0; i < NE - 2; i++) 1149239299Shselasky { 1150239299Shselasky if (einv[i] != 0) 1151239299Shselasky return 0; 1152184610Salfred } 1153184610Salfred 1154184610Salfred /* Fail if the computed inverse is out of range. */ 1155239180Shselasky if (eisinf (einv) || eisnan (einv) || (ecmp (einv, ezero) == 0)) 1156239180Shselasky return 0; 1157184610Salfred 1158239299Shselasky /* Output the reciprocal and return success flag. */ 1159239180Shselasky PUT_REAL (einv, r); 1160239180Shselasky return 1; 1161239299Shselasky} 1162239299Shselasky#endif /* REAL_ARITHMETIC defined */ 1163239180Shselasky 1164239180Shselasky/* Used for debugging--print the value of R in human-readable format 1165239180Shselasky on stderr. */ 1166239180Shselasky 1167239180Shselaskyvoid 1168239180Shselaskydebug_real (r) 1169239299Shselasky REAL_VALUE_TYPE r; 1170239180Shselasky{ 1171239180Shselasky char dstr[30]; 1172239180Shselasky 1173192984Sthompsa REAL_VALUE_TO_DECIMAL (r, "%.20g", dstr); 1174184610Salfred fprintf (stderr, "%s", dstr); 1175184610Salfred} 1176264149Sian 1177264149Sian 1178264149Sian/* The following routines convert REAL_VALUE_TYPE to the various floating 1179264149Sian point formats that are meaningful to supported computers. 1180264149Sian 1181297583Sian The results are returned in 32-bit pieces, each piece stored in a `long'. 1182264149Sian This is so they can be printed by statements like 1183184610Salfred 1184264149Sian fprintf (file, "%lx, %lx", L[0], L[1]); 1185264149Sian 1186264149Sian that will work on both narrow- and wide-word host computers. */ 1187184610Salfred 1188184610Salfred/* Convert R to a 128-bit long double precision value. The output array L 1189264149Sian contains four 32-bit pieces of the result, in the order they would appear 1190264149Sian in memory. */ 1191264149Sian 1192184610Salfredvoid 1193297583Sianetartdouble (r, l) 1194184610Salfred REAL_VALUE_TYPE r; 1195184610Salfred long l[]; 1196184610Salfred{ 1197194677Sthompsa unsigned EMUSHORT e[NE]; 1198184610Salfred 1199194677Sthompsa GET_REAL (&r, e); 1200194677Sthompsa etoe113 (e, e); 1201264031Sian endian (e, l, TFmode); 1202264031Sian} 1203184610Salfred 1204184610Salfred/* Convert R to a double extended precision value. The output array L 1205297583Sian contains three 32-bit pieces of the result, in the order they would 1206297583Sian appear in memory. */ 1207184610Salfred 1208264031Sianvoid 1209264031Sianetarldouble (r, l) 1210264031Sian REAL_VALUE_TYPE r; 1211264031Sian long l[]; 1212264031Sian{ 1213264031Sian unsigned EMUSHORT e[NE]; 1214184610Salfred 1215184610Salfred GET_REAL (&r, e); 1216264031Sian etoe64 (e, e); 1217264031Sian endian (e, l, XFmode); 1218264031Sian} 1219264031Sian 1220264031Sian/* Convert R to a double precision value. The output array L contains two 1221264031Sian 32-bit pieces of the result, in the order they would appear in memory. */ 1222264800Shselasky 1223264800Shselaskyvoid 1224264800Shselaskyetardouble (r, l) 1225264800Shselasky REAL_VALUE_TYPE r; 1226264800Shselasky long l[]; 1227264031Sian{ 1228194677Sthompsa unsigned EMUSHORT e[NE]; 1229264031Sian 1230264800Shselasky GET_REAL (&r, e); 1231264800Shselasky etoe53 (e, e); 1232264800Shselasky endian (e, l, DFmode); 1233264031Sian} 1234264031Sian 1235264031Sian/* Convert R to a single precision float value stored in the least-significant 1236264031Sian bits of a `long'. */ 1237264031Sian 1238264031Sianlong 1239264031Sianetarsingle (r) 1240264031Sian REAL_VALUE_TYPE r; 1241264031Sian{ 1242264031Sian unsigned EMUSHORT e[NE]; 1243184610Salfred long l; 1244264031Sian 1245264031Sian GET_REAL (&r, e); 1246264031Sian etoe24 (e, e); 1247194228Sthompsa endian (e, &l, SFmode); 1248184610Salfred return ((long) l); 1249264031Sian} 1250184610Salfred 1251184610Salfred/* Convert X to a decimal ASCII string S for output to an assembly 1252184610Salfred language file. Note, there is no standard way to spell infinity or 1253184610Salfred a NaN, so these values may require special treatment in the tm.h 1254194677Sthompsa macros. */ 1255184610Salfred 1256194677Sthompsavoid 1257194677Sthompsaereal_to_decimal (x, s) 1258184610Salfred REAL_VALUE_TYPE x; 1259186730Salfred char *s; 1260184610Salfred{ 1261184610Salfred unsigned EMUSHORT e[NE]; 1262264031Sian 1263264031Sian GET_REAL (&x, e); 1264264031Sian etoasc (e, s, 20); 1265264031Sian} 1266184610Salfred 1267297583Sian/* Compare X and Y. Return 1 if X > Y, 0 if X == Y, -1 if X < Y, 1268297583Sian or -2 if either is a NaN. */ 1269264031Sian 1270194677Sthompsaint 1271184610Salfredereal_cmp (x, y) 1272184610Salfred REAL_VALUE_TYPE x, y; 1273264031Sian{ 1274184610Salfred unsigned EMUSHORT ex[NE], ey[NE]; 1275264031Sian 1276264031Sian GET_REAL (&x, ex); 1277264031Sian GET_REAL (&y, ey); 1278264031Sian return (ecmp (ex, ey)); 1279264031Sian} 1280264031Sian 1281264031Sian/* Return 1 if the sign bit of X is set, else return 0. */ 1282264031Sian 1283264031Sianint 1284264031Sianereal_isneg (x) 1285264031Sian REAL_VALUE_TYPE x; 1286264031Sian{ 1287264031Sian unsigned EMUSHORT ex[NE]; 1288264031Sian 1289264031Sian GET_REAL (&x, ex); 1290264031Sian return (eisneg (ex)); 1291264031Sian} 1292264031Sian 1293264031Sian/* End of REAL_ARITHMETIC interface */ 1294264031Sian 1295264031Sian/* 1296264031Sian Extended precision IEEE binary floating point arithmetic routines 1297264031Sian 1298264031Sian Numbers are stored in C language as arrays of 16-bit unsigned 1299264031Sian short integers. The arguments of the routines are pointers to 1300264031Sian the arrays. 1301264031Sian 1302184610Salfred External e type data structure, similar to Intel 8087 chip 1303186730Salfred temporary real format but possibly with a larger significand: 1304184610Salfred 1305186730Salfred NE-1 significand words (least significant word first, 1306186730Salfred most significant bit is normally set) 1307186730Salfred exponent (value = EXONE for 1.0, 1308186730Salfred top bit is the sign) 1309186730Salfred 1310186730Salfred 1311186730Salfred Internal exploded e-type data structure of a number (a "word" is 16 bits): 1312186730Salfred 1313186730Salfred ei[0] sign word (0 for positive, 0xffff for negative) 1314184610Salfred ei[1] biased exponent (value = EXONE for the number 1.0) 1315184610Salfred ei[2] high guard word (always zero after normalization) 1316184610Salfred ei[3] 1317184610Salfred to ei[NI-2] significand (NI-4 significand words, 1318184610Salfred most significant word first, 1319184610Salfred most significant bit is set) 1320184610Salfred ei[NI-1] low guard word (0x8000 bit is rounding place) 1321184610Salfred 1322184610Salfred 1323194228Sthompsa 1324184610Salfred Routines for external format e-type numbers 1325264031Sian 1326184610Salfred asctoe (string, e) ASCII string to extended double e type 1327184610Salfred asctoe64 (string, &d) ASCII string to long double 1328194677Sthompsa asctoe53 (string, &d) ASCII string to double 1329194228Sthompsa asctoe24 (string, &f) ASCII string to single 1330184610Salfred asctoeg (string, e, prec) ASCII string to specified precision 1331184610Salfred e24toe (&f, e) IEEE single precision to e type 1332184610Salfred e53toe (&d, e) IEEE double precision to e type 1333194677Sthompsa e64toe (&d, e) IEEE long double precision to e type 1334188413Sthompsa e113toe (&d, e) 128-bit long double precision to e type 1335194677Sthompsa#if 0 1336188413Sthompsa eabs (e) absolute value 1337184610Salfred#endif 1338184610Salfred eadd (a, b, c) c = b + a 1339184610Salfred eclear (e) e = 0 1340184610Salfred ecmp (a, b) Returns 1 if a > b, 0 if a == b, 1341184610Salfred -1 if a < b, -2 if either a or b is a NaN. 1342184610Salfred ediv (a, b, c) c = b / a 1343192984Sthompsa efloor (a, b) truncate to integer, toward -infinity 1344184610Salfred efrexp (a, exp, s) extract exponent and significand 1345184610Salfred eifrac (e, &l, frac) e to HOST_WIDE_INT and e type fraction 1346184610Salfred euifrac (e, &l, frac) e to unsigned HOST_WIDE_INT and e type fraction 1347184610Salfred einfin (e) set e to infinity, leaving its sign alone 1348192984Sthompsa eldexp (a, n, b) multiply by 2**n 1349184610Salfred emov (a, b) b = a 1350297583Sian emul (a, b, c) c = b * a 1351297583Sian eneg (e) e = -e 1352184610Salfred#if 0 1353184610Salfred eround (a, b) b = nearest integer value to a 1354184610Salfred#endif 1355184610Salfred esub (a, b, c) c = b - a 1356184610Salfred#if 0 1357184610Salfred e24toasc (&f, str, n) single to ASCII string, n digits after decimal 1358184610Salfred e53toasc (&d, str, n) double to ASCII string, n digits after decimal 1359194228Sthompsa e64toasc (&d, str, n) 80-bit long double to ASCII string 1360188413Sthompsa e113toasc (&d, str, n) 128-bit long double to ASCII string 1361184610Salfred#endif 1362184610Salfred etoasc (e, str, n) e to ASCII string, n digits after decimal 1363184610Salfred etoe24 (e, &f) convert e type to IEEE single precision 1364192984Sthompsa etoe53 (e, &d) convert e type to IEEE double precision 1365184610Salfred etoe64 (e, &d) convert e type to IEEE long double precision 1366184610Salfred ltoe (&l, e) HOST_WIDE_INT to e type 1367184610Salfred ultoe (&l, e) unsigned HOST_WIDE_INT to e type 1368184610Salfred eisneg (e) 1 if sign bit of e != 0, else 0 1369192984Sthompsa eisinf (e) 1 if e has maximum exponent (non-IEEE) 1370184610Salfred or is infinite (IEEE) 1371297583Sian eisnan (e) 1 if e is a NaN 1372297583Sian 1373184610Salfred 1374184610Salfred Routines for internal format exploded e-type numbers 1375184610Salfred 1376184610Salfred eaddm (ai, bi) add significands, bi = bi + ai 1377184610Salfred ecleaz (ei) ei = 0 1378184610Salfred ecleazs (ei) set ei = 0 but leave its sign alone 1379184610Salfred ecmpm (ai, bi) compare significands, return 1, 0, or -1 1380194228Sthompsa edivm (ai, bi) divide significands, bi = bi / ai 1381188413Sthompsa emdnorm (ai,l,s,exp) normalize and round off 1382184610Salfred emovi (a, ai) convert external a to internal ai 1383184610Salfred emovo (ai, a) convert internal ai to external a 1384184610Salfred emovz (ai, bi) bi = ai, low guard word of bi = 0 1385192984Sthompsa emulm (ai, bi) multiply significands, bi = bi * ai 1386184610Salfred enormlz (ei) left-justify the significand 1387184610Salfred eshdn1 (ai) shift significand and guards down 1 bit 1388184610Salfred eshdn8 (ai) shift down 8 bits 1389184610Salfred eshdn6 (ai) shift down 16 bits 1390192984Sthompsa eshift (ai, n) shift ai n bits up (or down if n < 0) 1391184610Salfred eshup1 (ai) shift significand and guards up 1 bit 1392297583Sian eshup8 (ai) shift up 8 bits 1393297583Sian eshup6 (ai) shift up 16 bits 1394184610Salfred esubm (ai, bi) subtract significands, bi = bi - ai 1395184610Salfred eiisinf (ai) 1 if infinite 1396184610Salfred eiisnan (ai) 1 if a NaN 1397184610Salfred eiisneg (ai) 1 if sign bit of ai != 0, else 0 1398184610Salfred einan (ai) set ai = NaN 1399184610Salfred#if 0 1400184610Salfred eiinfin (ai) set ai = infinity 1401184610Salfred#endif 1402184610Salfred 1403184610Salfred The result is always normalized and rounded to NI-4 word precision 1404184610Salfred after each arithmetic operation. 1405184610Salfred 1406184610Salfred Exception flags are NOT fully supported. 1407194228Sthompsa 1408188413Sthompsa Signaling NaN's are NOT supported; they are treated the same 1409184610Salfred as quiet NaN's. 1410184610Salfred 1411264010Sian Define INFINITY for support of infinity; otherwise a 1412264010Sian saturation arithmetic is implemented. 1413264010Sian 1414264010Sian Define NANS for support of Not-a-Number items; otherwise the 1415264010Sian arithmetic will never produce a NaN output, and might be confused 1416264010Sian by a NaN input. 1417264010Sian If NaN's are supported, the output of `ecmp (a,b)' is -2 if 1418264010Sian either a or b is a NaN. This means asking `if (ecmp (a,b) < 0)' 1419264010Sian may not be legitimate. Use `if (ecmp (a,b) == -1)' for `less than' 1420264010Sian if in doubt. 1421264010Sian 1422184610Salfred Denormals are always supported here where appropriate (e.g., not 1423264010Sian for conversion to DEC numbers). */ 1424264010Sian 1425184610Salfred/* Definitions for error codes that are passed to the common error handling 1426264010Sian routine mtherr. 1427264010Sian 1428264010Sian For Digital Equipment PDP-11 and VAX computers, certain 1429264010Sian IBM systems, and others that use numbers with a 56-bit 1430237102Smarius significand, the symbol DEC should be defined. In this 1431264010Sian mode, most floating point constants are given as arrays 1432264010Sian of octal integers to eliminate decimal to binary conversion 1433264010Sian errors that might be introduced by the compiler. 1434264010Sian 1435264010Sian For computers, such as IBM PC, that follow the IEEE 1436264010Sian Standard for Binary Floating Point Arithmetic (ANSI/IEEE 1437264010Sian Std 754-1985), the symbol IEEE should be defined. 1438264010Sian These numbers have 53-bit significands. In this mode, constants 1439264010Sian are provided as arrays of hexadecimal 16 bit integers. 1440184610Salfred The endian-ness of generated values is controlled by 1441264010Sian REAL_WORDS_BIG_ENDIAN. 1442264010Sian 1443264010Sian To accommodate other types of computer arithmetic, all 1444184610Salfred constants are also provided in a normal decimal radix 1445264010Sian which one can hope are correctly converted to a suitable 1446264010Sian format by the available C language compiler. To invoke 1447264010Sian this mode, the symbol UNK is defined. 1448264010Sian 1449264010Sian An important difference among these modes is a predefined 1450264010Sian set of machine arithmetic constants for each. The numbers 1451264010Sian MACHEP (the machine roundoff error), MAXNUM (largest number 1452264010Sian represented), and several other parameters are preset by 1453264010Sian the configuration symbol. Check the file const.c to 1454264010Sian ensure that these values are correct for your computer. 1455264010Sian 1456264010Sian For ANSI C compatibility, define ANSIC equal to 1. Currently 1457264010Sian this affects only the atan2 function and others that use it. */ 1458264010Sian 1459264010Sian/* Constant definitions for math error conditions. */ 1460264010Sian 1461264010Sian#define DOMAIN 1 /* argument domain error */ 1462264010Sian#define SING 2 /* argument singularity */ 1463264010Sian#define OVERFLOW 3 /* overflow range error */ 1464264010Sian#define UNDERFLOW 4 /* underflow range error */ 1465264010Sian#define TLOSS 5 /* total loss of precision */ 1466184610Salfred#define PLOSS 6 /* partial loss of precision */ 1467184610Salfred#define INVALID 7 /* NaN-producing operation */ 1468264010Sian 1469264010Sian/* e type constants used by high precision check routines */ 1470264010Sian 1471264010Sian#if LONG_DOUBLE_TYPE_SIZE == 128 1472264010Sian/* 0.0 */ 1473264010Sianunsigned EMUSHORT ezero[NE] = 1474264010Sian {0x0000, 0x0000, 0x0000, 0x0000, 1475264010Sian 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000,}; 1476264010Sianextern unsigned EMUSHORT ezero[]; 1477264010Sian 1478264010Sian/* 5.0E-1 */ 1479264010Sianunsigned EMUSHORT ehalf[NE] = 1480264010Sian {0x0000, 0x0000, 0x0000, 0x0000, 1481264010Sian 0x0000, 0x0000, 0x0000, 0x0000, 0x8000, 0x3ffe,}; 1482264010Sianextern unsigned EMUSHORT ehalf[]; 1483264010Sian 1484264010Sian/* 1.0E0 */ 1485264010Sianunsigned EMUSHORT eone[NE] = 1486264010Sian {0x0000, 0x0000, 0x0000, 0x0000, 1487264010Sian 0x0000, 0x0000, 0x0000, 0x0000, 0x8000, 0x3fff,}; 1488264010Sianextern unsigned EMUSHORT eone[]; 1489264010Sian 1490264010Sian/* 2.0E0 */ 1491264010Sianunsigned EMUSHORT etwo[NE] = 1492264010Sian {0x0000, 0x0000, 0x0000, 0x0000, 1493264010Sian 0x0000, 0x0000, 0x0000, 0x0000, 0x8000, 0x4000,}; 1494264010Sianextern unsigned EMUSHORT etwo[]; 1495264010Sian 1496264010Sian/* 3.2E1 */ 1497264010Sianunsigned EMUSHORT e32[NE] = 1498264800Shselasky {0x0000, 0x0000, 0x0000, 0x0000, 1499264010Sian 0x0000, 0x0000, 0x0000, 0x0000, 0x8000, 0x4004,}; 1500264010Sianextern unsigned EMUSHORT e32[]; 1501264010Sian 1502264010Sian/* 6.93147180559945309417232121458176568075500134360255E-1 */ 1503264010Sianunsigned EMUSHORT elog2[NE] = 1504264010Sian {0x40f3, 0xf6af, 0x03f2, 0xb398, 1505264010Sian 0xc9e3, 0x79ab, 0150717, 0013767, 0130562, 0x3ffe,}; 1506264010Sianextern unsigned EMUSHORT elog2[]; 1507264010Sian 1508264010Sian/* 1.41421356237309504880168872420969807856967187537695E0 */ 1509264010Sianunsigned EMUSHORT esqrt2[NE] = 1510264010Sian {0x1d6f, 0xbe9f, 0x754a, 0x89b3, 1511264010Sian 0x597d, 0x6484, 0174736, 0171463, 0132404, 0x3fff,}; 1512264010Sianextern unsigned EMUSHORT esqrt2[]; 1513264010Sian 1514264010Sian/* 3.14159265358979323846264338327950288419716939937511E0 */ 1515264010Sianunsigned EMUSHORT epi[NE] = 1516264010Sian {0x2902, 0x1cd1, 0x80dc, 0x628b, 1517264010Sian 0xc4c6, 0xc234, 0020550, 0155242, 0144417, 0040000,}; 1518264010Sianextern unsigned EMUSHORT epi[]; 1519264010Sian 1520264010Sian#else 1521264010Sian/* LONG_DOUBLE_TYPE_SIZE is other than 128 */ 1522264010Sianunsigned EMUSHORT ezero[NE] = 1523264010Sian {0, 0000000, 0000000, 0000000, 0000000, 0000000,}; 1524264010Sianunsigned EMUSHORT ehalf[NE] = 1525264010Sian {0, 0000000, 0000000, 0000000, 0100000, 0x3ffe,}; 1526264010Sianunsigned EMUSHORT eone[NE] = 1527264010Sian {0, 0000000, 0000000, 0000000, 0100000, 0x3fff,}; 1528264010Sianunsigned EMUSHORT etwo[NE] = 1529264010Sian {0, 0000000, 0000000, 0000000, 0100000, 0040000,}; 1530264010Sianunsigned EMUSHORT e32[NE] = 1531264010Sian {0, 0000000, 0000000, 0000000, 0100000, 0040004,}; 1532264010Sianunsigned EMUSHORT elog2[NE] = 1533264010Sian {0xc9e4, 0x79ab, 0150717, 0013767, 0130562, 0x3ffe,}; 1534264010Sianunsigned EMUSHORT esqrt2[NE] = 1535264010Sian {0x597e, 0x6484, 0174736, 0171463, 0132404, 0x3fff,}; 1536264010Sianunsigned EMUSHORT epi[NE] = 1537264010Sian {0xc4c6, 0xc234, 0020550, 0155242, 0144417, 0040000,}; 1538264010Sian#endif 1539264010Sian 1540264010Sian/* Control register for rounding precision. 1541264010Sian This can be set to 113 (if NE=10), 80 (if NE=6), 64, 56, 53, or 24 bits. */ 1542264010Sian 1543264010Sianint rndprc = NBITS; 1544264010Sianextern int rndprc; 1545264010Sian 1546264010Sian/* Clear out entire e-type number X. */ 1547264010Sian 1548264010Sianstatic void 1549264010Sianeclear (x) 1550264010Sian register unsigned EMUSHORT *x; 1551264010Sian{ 1552264010Sian register int i; 1553264010Sian 1554264010Sian for (i = 0; i < NE; i++) 1555264010Sian *x++ = 0; 1556264010Sian} 1557264010Sian 1558264010Sian/* Move e-type number from A to B. */ 1559264010Sian 1560264010Sianstatic void 1561264010Sianemov (a, b) 1562264010Sian register unsigned EMUSHORT *a, *b; 1563264010Sian{ 1564264010Sian register int i; 1565264010Sian 1566264010Sian for (i = 0; i < NE; i++) 1567184610Salfred *b++ = *a++; 1568184610Salfred} 1569184610Salfred 1570184610Salfred 1571184610Salfred#if 0 1572184610Salfred/* Absolute value of e-type X. */ 1573184610Salfred 1574184610Salfredstatic void 1575184610Salfredeabs (x) 1576184610Salfred unsigned EMUSHORT x[]; 1577184610Salfred{ 1578184610Salfred /* sign is top bit of last word of external format */ 1579184610Salfred x[NE - 1] &= 0x7fff; 1580184610Salfred} 1581184610Salfred#endif /* 0 */ 1582184610Salfred 1583184610Salfred/* Negate the e-type number X. */ 1584184610Salfred 1585184610Salfredstatic void 1586184610Salfredeneg (x) 1587184610Salfred unsigned EMUSHORT x[]; 1588184610Salfred{ 1589184610Salfred 1590184610Salfred x[NE - 1] ^= 0x8000; /* Toggle the sign bit */ 1591184610Salfred} 1592184610Salfred 1593184610Salfred/* Return 1 if sign bit of e-type number X is nonzero, else zero. */ 1594184610Salfred 1595184610Salfredstatic int 1596184610Salfredeisneg (x) 1597184610Salfred unsigned EMUSHORT x[]; 1598184610Salfred{ 1599184610Salfred 1600184610Salfred if (x[NE - 1] & 0x8000) 1601184610Salfred return (1); 1602184610Salfred else 1603184610Salfred return (0); 1604184610Salfred} 1605184610Salfred 1606184610Salfred/* Return 1 if e-type number X is infinity, else return zero. */ 1607184610Salfred 1608184610Salfredstatic int 1609184610Salfredeisinf (x) 1610184610Salfred unsigned EMUSHORT x[]; 1611184610Salfred{ 1612184610Salfred 1613184610Salfred#ifdef NANS 1614192984Sthompsa if (eisnan (x)) 1615184610Salfred return (0); 1616184610Salfred#endif 1617184610Salfred if ((x[NE - 1] & 0x7fff) == 0x7fff) 1618184610Salfred return (1); 1619184610Salfred else 1620264010Sian return (0); 1621184610Salfred} 1622184610Salfred 1623184610Salfred/* Check if e-type number is not a number. The bit pattern is one that we 1624192984Sthompsa defined, so we know for sure how to detect it. */ 1625184610Salfred 1626184610Salfredstatic int 1627184610Salfredeisnan (x) 1628184610Salfred unsigned EMUSHORT x[]; 1629192984Sthompsa{ 1630184610Salfred#ifdef NANS 1631297583Sian int i; 1632297583Sian 1633264010Sian /* NaN has maximum exponent */ 1634184610Salfred if ((x[NE - 1] & 0x7fff) != 0x7fff) 1635184610Salfred return (0); 1636184610Salfred /* ... and non-zero significand field. */ 1637184610Salfred for (i = 0; i < NE - 1; i++) 1638184610Salfred { 1639184610Salfred if (*x++ != 0) 1640184610Salfred return (1); 1641264010Sian } 1642264010Sian#endif 1643184610Salfred 1644194228Sthompsa return (0); 1645188413Sthompsa} 1646184610Salfred 1647184610Salfred/* Fill e-type number X with infinity pattern (IEEE) 1648184610Salfred or largest possible number (non-IEEE). */ 1649184610Salfred 1650184610Salfredstatic void 1651184610Salfredeinfin (x) 1652194228Sthompsa register unsigned EMUSHORT *x; 1653188413Sthompsa{ 1654184610Salfred register int i; 1655184610Salfred 1656184610Salfred#ifdef INFINITY 1657184610Salfred for (i = 0; i < NE - 1; i++) 1658184610Salfred *x++ = 0; 1659184610Salfred *x |= 32767; 1660194228Sthompsa#else 1661188413Sthompsa for (i = 0; i < NE - 1; i++) 1662184610Salfred *x++ = 0xffff; 1663184610Salfred *x |= 32766; 1664184610Salfred if (rndprc < NBITS) 1665192984Sthompsa { 1666184610Salfred if (rndprc == 113) 1667184610Salfred { 1668184610Salfred *(x - 9) = 0; 1669297583Sian *(x - 8) = 0; 1670184610Salfred } 1671184610Salfred if (rndprc == 64) 1672184610Salfred { 1673184610Salfred *(x - 5) = 0; 1674184610Salfred } 1675264149Sian if (rndprc == 53) 1676264149Sian { 1677264149Sian *(x - 4) = 0xf800; 1678264149Sian } 1679264149Sian else 1680264149Sian { 1681297583Sian *(x - 4) = 0; 1682297583Sian *(x - 3) = 0; 1683264149Sian *(x - 2) = 0xff00; 1684264149Sian } 1685264149Sian } 1686264149Sian#endif 1687264149Sian} 1688264149Sian 1689264149Sian/* Output an e-type NaN. 1690264149Sian This generates Intel's quiet NaN pattern for extended real. 1691264149Sian The exponent is 7fff, the leading mantissa word is c000. */ 1692264149Sian 1693264149Sianstatic void 1694264149Sianenan (x, sign) 1695264149Sian register unsigned EMUSHORT *x; 1696264149Sian int sign; 1697264149Sian{ 1698286385Sian register int i; 1699264149Sian 1700297583Sian for (i = 0; i < NE - 2; i++) 1701297583Sian *x++ = 0; 1702264149Sian *x++ = 0xc000; 1703264149Sian *x = (sign << 15) | 0x7fff; 1704264149Sian} 1705264149Sian 1706264149Sian/* Move in an e-type number A, converting it to exploded e-type B. */ 1707264149Sian 1708264149Sianstatic void 1709264149Sianemovi (a, b) 1710264149Sian unsigned EMUSHORT *a, *b; 1711264149Sian{ 1712264149Sian register unsigned EMUSHORT *p, *q; 1713286385Sian int i; 1714286385Sian 1715286385Sian q = b; 1716286385Sian p = a + (NE - 1); /* point to last word of external number */ 1717286385Sian /* get the sign bit */ 1718264149Sian if (*p & 0x8000) 1719264149Sian *q++ = 0xffff; 1720264149Sian else 1721286385Sian *q++ = 0; 1722264149Sian /* get the exponent */ 1723264149Sian *q = *p--; 1724264149Sian *q++ &= 0x7fff; /* delete the sign bit */ 1725264149Sian#ifdef INFINITY 1726297583Sian if ((*(q - 1) & 0x7fff) == 0x7fff) 1727297583Sian { 1728281367Sian#ifdef NANS 1729264149Sian if (eisnan (a)) 1730264149Sian { 1731264149Sian *q++ = 0; 1732264149Sian for (i = 3; i < NI; i++) 1733264149Sian *q++ = *p--; 1734264149Sian return; 1735286385Sian } 1736264149Sian#endif 1737264149Sian 1738264149Sian for (i = 2; i < NI; i++) 1739264149Sian *q++ = 0; 1740264149Sian return; 1741264149Sian } 1742264149Sian#endif 1743264149Sian 1744264149Sian /* clear high guard word */ 1745297583Sian *q++ = 0; 1746297583Sian /* move in the significand */ 1747264149Sian for (i = 0; i < NE - 1; i++) 1748264149Sian *q++ = *p--; 1749264149Sian /* clear low guard word */ 1750264149Sian *q = 0; 1751264149Sian} 1752264149Sian 1753264149Sian/* Move out exploded e-type number A, converting it to e type B. */ 1754264149Sian 1755264149Sianstatic void 1756264149Sianemovo (a, b) 1757264149Sian unsigned EMUSHORT *a, *b; 1758264149Sian{ 1759264149Sian register unsigned EMUSHORT *p, *q; 1760264149Sian unsigned EMUSHORT i; 1761264149Sian int j; 1762264149Sian 1763264149Sian p = a; 1764264149Sian q = b + (NE - 1); /* point to output exponent */ 1765264149Sian /* combine sign and exponent */ 1766264149Sian i = *p++; 1767264149Sian if (i) 1768297583Sian *q-- = *p++ | 0x8000; 1769297583Sian else 1770281367Sian *q-- = *p++; 1771264149Sian#ifdef INFINITY 1772264149Sian if (*(p - 1) == 0x7fff) 1773264149Sian { 1774264149Sian#ifdef NANS 1775264149Sian if (eiisnan (a)) 1776264149Sian { 1777264149Sian enan (b, eiisneg (a)); 1778264149Sian return; 1779264149Sian } 1780264149Sian#endif 1781264149Sian einfin (b); 1782264149Sian return; 1783264149Sian } 1784264149Sian#endif 1785264149Sian /* skip over guard word */ 1786264149Sian ++p; 1787264149Sian /* move the significand */ 1788264149Sian for (j = 0; j < NE - 1; j++) 1789264149Sian *q-- = *p++; 1790297583Sian} 1791297583Sian 1792264149Sian/* Clear out exploded e-type number XI. */ 1793264149Sian 1794264149Sianstatic void 1795264149Sianecleaz (xi) 1796264149Sian register unsigned EMUSHORT *xi; 1797264149Sian{ 1798264149Sian register int i; 1799264149Sian 1800264149Sian for (i = 0; i < NI; i++) 1801264149Sian *xi++ = 0; 1802264149Sian} 1803264149Sian 1804264149Sian/* Clear out exploded e-type XI, but don't touch the sign. */ 1805264149Sian 1806264149Sianstatic void 1807264149Sianecleazs (xi) 1808264149Sian register unsigned EMUSHORT *xi; 1809264149Sian{ 1810264149Sian register int i; 1811297583Sian 1812297583Sian ++xi; 1813264149Sian for (i = 0; i < NI - 1; i++) 1814264149Sian *xi++ = 0; 1815264149Sian} 1816264149Sian 1817264149Sian/* Move exploded e-type number from A to B. */ 1818264149Sian 1819264149Sianstatic void 1820264149Sianemovz (a, b) 1821264149Sian register unsigned EMUSHORT *a, *b; 1822264149Sian{ 1823264149Sian register int i; 1824264149Sian 1825264149Sian for (i = 0; i < NI - 1; i++) 1826286382Sian *b++ = *a++; 1827286382Sian /* clear low guard word */ 1828286382Sian *b = 0; 1829286382Sian} 1830286382Sian 1831286382Sian/* Generate exploded e-type NaN. 1832286382Sian The explicit pattern for this is maximum exponent and 1833297583Sian top two significant bits set. */ 1834297583Sian 1835286382Sianstatic void 1836286382Sianeinan (x) 1837286382Sian unsigned EMUSHORT x[]; 1838286382Sian{ 1839286382Sian 1840286382Sian ecleaz (x); 1841286382Sian x[E] = 0x7fff; 1842286382Sian x[M + 1] = 0xc000; 1843286382Sian} 1844286382Sian 1845286382Sian/* Return nonzero if exploded e-type X is a NaN. */ 1846286382Sian 1847286382Sianstatic int 1848286382Sianeiisnan (x) 1849286382Sian unsigned EMUSHORT x[]; 1850286382Sian{ 1851286382Sian int i; 1852286382Sian 1853286382Sian if ((x[E] & 0x7fff) == 0x7fff) 1854286382Sian { 1855286382Sian for (i = M + 1; i < NI; i++) 1856286382Sian { 1857286382Sian if (x[i] != 0) 1858286382Sian return (1); 1859286382Sian } 1860286382Sian } 1861286382Sian return (0); 1862286382Sian} 1863297583Sian 1864297583Sian/* Return nonzero if sign of exploded e-type X is nonzero. */ 1865286382Sian 1866286382Sianstatic int 1867286382Sianeiisneg (x) 1868286382Sian unsigned EMUSHORT x[]; 1869286382Sian{ 1870286382Sian 1871286382Sian return x[0] != 0; 1872286382Sian} 1873286382Sian 1874286382Sian#if 0 1875286382Sian/* Fill exploded e-type X with infinity pattern. 1876286382Sian This has maximum exponent and significand all zeros. */ 1877286382Sian 1878286382Sianstatic void 1879286382Sianeiinfin (x) 1880286382Sian unsigned EMUSHORT x[]; 1881286382Sian{ 1882286382Sian 1883286382Sian ecleaz (x); 1884286382Sian x[E] = 0x7fff; 1885286382Sian} 1886286382Sian#endif /* 0 */ 1887286382Sian 1888286382Sian/* Return nonzero if exploded e-type X is infinite. */ 1889286382Sian 1890286382Sianstatic int 1891297583Sianeiisinf (x) 1892297583Sian unsigned EMUSHORT x[]; 1893286382Sian{ 1894286382Sian 1895286382Sian#ifdef NANS 1896286382Sian if (eiisnan (x)) 1897286382Sian return (0); 1898286382Sian#endif 1899286382Sian if ((x[E] & 0x7fff) == 0x7fff) 1900286382Sian return (1); 1901286382Sian return (0); 1902286382Sian} 1903286382Sian 1904286382Sian 1905286382Sian/* Compare significands of numbers in internal exploded e-type format. 1906286382Sian Guard words are included in the comparison. 1907286382Sian 1908264149Sian Returns +1 if a > b 1909264149Sian 0 if a == b 1910264149Sian -1 if a < b */ 1911264149Sian 1912264149Sianstatic int 1913264149Sianecmpm (a, b) 1914264149Sian register unsigned EMUSHORT *a, *b; 1915264149Sian{ 1916264149Sian int i; 1917264149Sian 1918264149Sian a += M; /* skip up to significand area */ 1919264149Sian b += M; 1920264149Sian for (i = M; i < NI; i++) 1921264149Sian { 1922264149Sian if (*a++ != *b++) 1923264149Sian goto difrnt; 1924264149Sian } 1925264149Sian return (0); 1926264149Sian 1927264149Sian difrnt: 1928264149Sian if (*(--a) > *(--b)) 1929264149Sian return (1); 1930286385Sian else 1931264149Sian return (-1); 1932264149Sian} 1933264149Sian 1934264149Sian/* Shift significand of exploded e-type X down by 1 bit. */ 1935264149Sian 1936264149Sianstatic void 1937264149Sianeshdn1 (x) 1938264149Sian register unsigned EMUSHORT *x; 1939264800Shselasky{ 1940264800Shselasky register unsigned EMUSHORT bits; 1941264800Shselasky int i; 1942264149Sian 1943264149Sian x += M; /* point to significand area */ 1944264149Sian 1945264149Sian bits = 0; 1946264149Sian for (i = M; i < NI; i++) 1947264149Sian { 1948286382Sian if (*x & 1) 1949286382Sian bits |= 1; 1950286382Sian *x >>= 1; 1951286382Sian if (bits & 2) 1952286382Sian *x |= 0x8000; 1953286382Sian bits <<= 1; 1954286382Sian ++x; 1955286382Sian } 1956286382Sian} 1957264149Sian 1958264149Sian/* Shift significand of exploded e-type X up by 1 bit. */ 1959264149Sian 1960264149Sianstatic void 1961264149Sianeshup1 (x) 1962264149Sian register unsigned EMUSHORT *x; 1963264149Sian{ 1964264149Sian register unsigned EMUSHORT bits; 1965184610Salfred int i; 1966192984Sthompsa 1967184610Salfred x += NI - 1; 1968184610Salfred bits = 0; 1969184610Salfred 1970194228Sthompsa for (i = M; i < NI; i++) 1971184610Salfred { 1972184610Salfred if (*x & 0x8000) 1973184610Salfred bits |= 1; 1974192984Sthompsa *x <<= 1; 1975184610Salfred if (bits & 2) 1976184610Salfred *x |= 1; 1977184610Salfred bits <<= 1; 1978194228Sthompsa --x; 1979184610Salfred } 1980184610Salfred} 1981184610Salfred 1982192984Sthompsa 1983184610Salfred/* Shift significand of exploded e-type X down by 8 bits. */ 1984184610Salfred 1985184610Salfredstatic void 1986194228Sthompsaeshdn8 (x) 1987184610Salfred register unsigned EMUSHORT *x; 1988184610Salfred{ 1989184610Salfred register unsigned EMUSHORT newbyt, oldbyt; 1990192984Sthompsa int i; 1991184610Salfred 1992184610Salfred x += M; 1993184610Salfred oldbyt = 0; 1994194228Sthompsa for (i = M; i < NI; i++) 1995184610Salfred { 1996184610Salfred newbyt = *x << 8; 1997197570Sthompsa *x >>= 8; 1998197570Sthompsa *x |= oldbyt; 1999197570Sthompsa oldbyt = newbyt; 2000197570Sthompsa ++x; 2001237102Smarius } 2002197570Sthompsa} 2003197570Sthompsa 2004/* Shift significand of exploded e-type X up by 8 bits. */ 2005 2006static void 2007eshup8 (x) 2008 register unsigned EMUSHORT *x; 2009{ 2010 int i; 2011 register unsigned EMUSHORT newbyt, oldbyt; 2012 2013 x += NI - 1; 2014 oldbyt = 0; 2015 2016 for (i = M; i < NI; i++) 2017 { 2018 newbyt = *x >> 8; 2019 *x <<= 8; 2020 *x |= oldbyt; 2021 oldbyt = newbyt; 2022 --x; 2023 } 2024} 2025 2026/* Shift significand of exploded e-type X up by 16 bits. */ 2027 2028static void 2029eshup6 (x) 2030 register unsigned EMUSHORT *x; 2031{ 2032 int i; 2033 register unsigned EMUSHORT *p; 2034 2035 p = x + M; 2036 x += M + 1; 2037 2038 for (i = M; i < NI - 1; i++) 2039 *p++ = *x++; 2040 2041 *p = 0; 2042} 2043 2044/* Shift significand of exploded e-type X down by 16 bits. */ 2045 2046static void 2047eshdn6 (x) 2048 register unsigned EMUSHORT *x; 2049{ 2050 int i; 2051 register unsigned EMUSHORT *p; 2052 2053 x += NI - 1; 2054 p = x + 1; 2055 2056 for (i = M; i < NI - 1; i++) 2057 *(--p) = *(--x); 2058 2059 *(--p) = 0; 2060} 2061 2062/* Add significands of exploded e-type X and Y. X + Y replaces Y. */ 2063 2064static void 2065eaddm (x, y) 2066 unsigned EMUSHORT *x, *y; 2067{ 2068 register unsigned EMULONG a; 2069 int i; 2070 unsigned int carry; 2071 2072 x += NI - 1; 2073 y += NI - 1; 2074 carry = 0; 2075 for (i = M; i < NI; i++) 2076 { 2077 a = (unsigned EMULONG) (*x) + (unsigned EMULONG) (*y) + carry; 2078 if (a & 0x10000) 2079 carry = 1; 2080 else 2081 carry = 0; 2082 *y = (unsigned EMUSHORT) a; 2083 --x; 2084 --y; 2085 } 2086} 2087 2088/* Subtract significands of exploded e-type X and Y. Y - X replaces Y. */ 2089 2090static void 2091esubm (x, y) 2092 unsigned EMUSHORT *x, *y; 2093{ 2094 unsigned EMULONG a; 2095 int i; 2096 unsigned int carry; 2097 2098 x += NI - 1; 2099 y += NI - 1; 2100 carry = 0; 2101 for (i = M; i < NI; i++) 2102 { 2103 a = (unsigned EMULONG) (*y) - (unsigned EMULONG) (*x) - carry; 2104 if (a & 0x10000) 2105 carry = 1; 2106 else 2107 carry = 0; 2108 *y = (unsigned EMUSHORT) a; 2109 --x; 2110 --y; 2111 } 2112} 2113 2114 2115static unsigned EMUSHORT equot[NI]; 2116 2117 2118#if 0 2119/* Radix 2 shift-and-add versions of multiply and divide */ 2120 2121 2122/* Divide significands */ 2123 2124int 2125edivm (den, num) 2126 unsigned EMUSHORT den[], num[]; 2127{ 2128 int i; 2129 register unsigned EMUSHORT *p, *q; 2130 unsigned EMUSHORT j; 2131 2132 p = &equot[0]; 2133 *p++ = num[0]; 2134 *p++ = num[1]; 2135 2136 for (i = M; i < NI; i++) 2137 { 2138 *p++ = 0; 2139 } 2140 2141 /* Use faster compare and subtraction if denominator has only 15 bits of 2142 significance. */ 2143 2144 p = &den[M + 2]; 2145 if (*p++ == 0) 2146 { 2147 for (i = M + 3; i < NI; i++) 2148 { 2149 if (*p++ != 0) 2150 goto fulldiv; 2151 } 2152 if ((den[M + 1] & 1) != 0) 2153 goto fulldiv; 2154 eshdn1 (num); 2155 eshdn1 (den); 2156 2157 p = &den[M + 1]; 2158 q = &num[M + 1]; 2159 2160 for (i = 0; i < NBITS + 2; i++) 2161 { 2162 if (*p <= *q) 2163 { 2164 *q -= *p; 2165 j = 1; 2166 } 2167 else 2168 { 2169 j = 0; 2170 } 2171 eshup1 (equot); 2172 equot[NI - 2] |= j; 2173 eshup1 (num); 2174 } 2175 goto divdon; 2176 } 2177 2178 /* The number of quotient bits to calculate is NBITS + 1 scaling guard 2179 bit + 1 roundoff bit. */ 2180 2181 fulldiv: 2182 2183 p = &equot[NI - 2]; 2184 for (i = 0; i < NBITS + 2; i++) 2185 { 2186 if (ecmpm (den, num) <= 0) 2187 { 2188 esubm (den, num); 2189 j = 1; /* quotient bit = 1 */ 2190 } 2191 else 2192 j = 0; 2193 eshup1 (equot); 2194 *p |= j; 2195 eshup1 (num); 2196 } 2197 2198 divdon: 2199 2200 eshdn1 (equot); 2201 eshdn1 (equot); 2202 2203 /* test for nonzero remainder after roundoff bit */ 2204 p = &num[M]; 2205 j = 0; 2206 for (i = M; i < NI; i++) 2207 { 2208 j |= *p++; 2209 } 2210 if (j) 2211 j = 1; 2212 2213 2214 for (i = 0; i < NI; i++) 2215 num[i] = equot[i]; 2216 return ((int) j); 2217} 2218 2219 2220/* Multiply significands */ 2221 2222int 2223emulm (a, b) 2224 unsigned EMUSHORT a[], b[]; 2225{ 2226 unsigned EMUSHORT *p, *q; 2227 int i, j, k; 2228 2229 equot[0] = b[0]; 2230 equot[1] = b[1]; 2231 for (i = M; i < NI; i++) 2232 equot[i] = 0; 2233 2234 p = &a[NI - 2]; 2235 k = NBITS; 2236 while (*p == 0) /* significand is not supposed to be zero */ 2237 { 2238 eshdn6 (a); 2239 k -= 16; 2240 } 2241 if ((*p & 0xff) == 0) 2242 { 2243 eshdn8 (a); 2244 k -= 8; 2245 } 2246 2247 q = &equot[NI - 1]; 2248 j = 0; 2249 for (i = 0; i < k; i++) 2250 { 2251 if (*p & 1) 2252 eaddm (b, equot); 2253 /* remember if there were any nonzero bits shifted out */ 2254 if (*q & 1) 2255 j |= 1; 2256 eshdn1 (a); 2257 eshdn1 (equot); 2258 } 2259 2260 for (i = 0; i < NI; i++) 2261 b[i] = equot[i]; 2262 2263 /* return flag for lost nonzero bits */ 2264 return (j); 2265} 2266 2267#else 2268 2269/* Radix 65536 versions of multiply and divide. */ 2270 2271/* Multiply significand of e-type number B 2272 by 16-bit quantity A, return e-type result to C. */ 2273 2274static void 2275m16m (a, b, c) 2276 unsigned int a; 2277 unsigned EMUSHORT b[], c[]; 2278{ 2279 register unsigned EMUSHORT *pp; 2280 register unsigned EMULONG carry; 2281 unsigned EMUSHORT *ps; 2282 unsigned EMUSHORT p[NI]; 2283 unsigned EMULONG aa, m; 2284 int i; 2285 2286 aa = a; 2287 pp = &p[NI-2]; 2288 *pp++ = 0; 2289 *pp = 0; 2290 ps = &b[NI-1]; 2291 2292 for (i=M+1; i<NI; i++) 2293 { 2294 if (*ps == 0) 2295 { 2296 --ps; 2297 --pp; 2298 *(pp-1) = 0; 2299 } 2300 else 2301 { 2302 m = (unsigned EMULONG) aa * *ps--; 2303 carry = (m & 0xffff) + *pp; 2304 *pp-- = (unsigned EMUSHORT)carry; 2305 carry = (carry >> 16) + (m >> 16) + *pp; 2306 *pp = (unsigned EMUSHORT)carry; 2307 *(pp-1) = carry >> 16; 2308 } 2309 } 2310 for (i=M; i<NI; i++) 2311 c[i] = p[i]; 2312} 2313 2314/* Divide significands of exploded e-types NUM / DEN. Neither the 2315 numerator NUM nor the denominator DEN is permitted to have its high guard 2316 word nonzero. */ 2317 2318static int 2319edivm (den, num) 2320 unsigned EMUSHORT den[], num[]; 2321{ 2322 int i; 2323 register unsigned EMUSHORT *p; 2324 unsigned EMULONG tnum; 2325 unsigned EMUSHORT j, tdenm, tquot; 2326 unsigned EMUSHORT tprod[NI+1]; 2327 2328 p = &equot[0]; 2329 *p++ = num[0]; 2330 *p++ = num[1]; 2331 2332 for (i=M; i<NI; i++) 2333 { 2334 *p++ = 0; 2335 } 2336 eshdn1 (num); 2337 tdenm = den[M+1]; 2338 for (i=M; i<NI; i++) 2339 { 2340 /* Find trial quotient digit (the radix is 65536). */ 2341 tnum = (((unsigned EMULONG) num[M]) << 16) + num[M+1]; 2342 2343 /* Do not execute the divide instruction if it will overflow. */ 2344 if ((tdenm * (unsigned long)0xffff) < tnum) 2345 tquot = 0xffff; 2346 else 2347 tquot = tnum / tdenm; 2348 /* Multiply denominator by trial quotient digit. */ 2349 m16m ((unsigned int)tquot, den, tprod); 2350 /* The quotient digit may have been overestimated. */ 2351 if (ecmpm (tprod, num) > 0) 2352 { 2353 tquot -= 1; 2354 esubm (den, tprod); 2355 if (ecmpm (tprod, num) > 0) 2356 { 2357 tquot -= 1; 2358 esubm (den, tprod); 2359 } 2360 } 2361 esubm (tprod, num); 2362 equot[i] = tquot; 2363 eshup6(num); 2364 } 2365 /* test for nonzero remainder after roundoff bit */ 2366 p = &num[M]; 2367 j = 0; 2368 for (i=M; i<NI; i++) 2369 { 2370 j |= *p++; 2371 } 2372 if (j) 2373 j = 1; 2374 2375 for (i=0; i<NI; i++) 2376 num[i] = equot[i]; 2377 2378 return ((int)j); 2379} 2380 2381/* Multiply significands of exploded e-type A and B, result in B. */ 2382 2383static int 2384emulm (a, b) 2385 unsigned EMUSHORT a[], b[]; 2386{ 2387 unsigned EMUSHORT *p, *q; 2388 unsigned EMUSHORT pprod[NI]; 2389 unsigned EMUSHORT j; 2390 int i; 2391 2392 equot[0] = b[0]; 2393 equot[1] = b[1]; 2394 for (i=M; i<NI; i++) 2395 equot[i] = 0; 2396 2397 j = 0; 2398 p = &a[NI-1]; 2399 q = &equot[NI-1]; 2400 for (i=M+1; i<NI; i++) 2401 { 2402 if (*p == 0) 2403 { 2404 --p; 2405 } 2406 else 2407 { 2408 m16m ((unsigned int) *p--, b, pprod); 2409 eaddm(pprod, equot); 2410 } 2411 j |= *q; 2412 eshdn6(equot); 2413 } 2414 2415 for (i=0; i<NI; i++) 2416 b[i] = equot[i]; 2417 2418 /* return flag for lost nonzero bits */ 2419 return ((int)j); 2420} 2421#endif 2422 2423 2424/* Normalize and round off. 2425 2426 The internal format number to be rounded is S. 2427 Input LOST is 0 if the value is exact. This is the so-called sticky bit. 2428 2429 Input SUBFLG indicates whether the number was obtained 2430 by a subtraction operation. In that case if LOST is nonzero 2431 then the number is slightly smaller than indicated. 2432 2433 Input EXP is the biased exponent, which may be negative. 2434 the exponent field of S is ignored but is replaced by 2435 EXP as adjusted by normalization and rounding. 2436 2437 Input RCNTRL is the rounding control. If it is nonzero, the 2438 returned value will be rounded to RNDPRC bits. 2439 2440 For future reference: In order for emdnorm to round off denormal 2441 significands at the right point, the input exponent must be 2442 adjusted to be the actual value it would have after conversion to 2443 the final floating point type. This adjustment has been 2444 implemented for all type conversions (etoe53, etc.) and decimal 2445 conversions, but not for the arithmetic functions (eadd, etc.). 2446 Data types having standard 15-bit exponents are not affected by 2447 this, but SFmode and DFmode are affected. For example, ediv with 2448 rndprc = 24 will not round correctly to 24-bit precision if the 2449 result is denormal. */ 2450 2451static int rlast = -1; 2452static int rw = 0; 2453static unsigned EMUSHORT rmsk = 0; 2454static unsigned EMUSHORT rmbit = 0; 2455static unsigned EMUSHORT rebit = 0; 2456static int re = 0; 2457static unsigned EMUSHORT rbit[NI]; 2458 2459static void 2460emdnorm (s, lost, subflg, exp, rcntrl) 2461 unsigned EMUSHORT s[]; 2462 int lost; 2463 int subflg; 2464 EMULONG exp; 2465 int rcntrl; 2466{ 2467 int i, j; 2468 unsigned EMUSHORT r; 2469 2470 /* Normalize */ 2471 j = enormlz (s); 2472 2473 /* a blank significand could mean either zero or infinity. */ 2474#ifndef INFINITY 2475 if (j > NBITS) 2476 { 2477 ecleazs (s); 2478 return; 2479 } 2480#endif 2481 exp -= j; 2482#ifndef INFINITY 2483 if (exp >= 32767L) 2484 goto overf; 2485#else 2486 if ((j > NBITS) && (exp < 32767)) 2487 { 2488 ecleazs (s); 2489 return; 2490 } 2491#endif 2492 if (exp < 0L) 2493 { 2494 if (exp > (EMULONG) (-NBITS - 1)) 2495 { 2496 j = (int) exp; 2497 i = eshift (s, j); 2498 if (i) 2499 lost = 1; 2500 } 2501 else 2502 { 2503 ecleazs (s); 2504 return; 2505 } 2506 } 2507 /* Round off, unless told not to by rcntrl. */ 2508 if (rcntrl == 0) 2509 goto mdfin; 2510 /* Set up rounding parameters if the control register changed. */ 2511 if (rndprc != rlast) 2512 { 2513 ecleaz (rbit); 2514 switch (rndprc) 2515 { 2516 default: 2517 case NBITS: 2518 rw = NI - 1; /* low guard word */ 2519 rmsk = 0xffff; 2520 rmbit = 0x8000; 2521 re = rw - 1; 2522 rebit = 1; 2523 break; 2524 2525 case 113: 2526 rw = 10; 2527 rmsk = 0x7fff; 2528 rmbit = 0x4000; 2529 rebit = 0x8000; 2530 re = rw; 2531 break; 2532 2533 case 64: 2534 rw = 7; 2535 rmsk = 0xffff; 2536 rmbit = 0x8000; 2537 re = rw - 1; 2538 rebit = 1; 2539 break; 2540 2541 /* For DEC or IBM arithmetic */ 2542 case 56: 2543 rw = 6; 2544 rmsk = 0xff; 2545 rmbit = 0x80; 2546 rebit = 0x100; 2547 re = rw; 2548 break; 2549 2550 case 53: 2551 rw = 6; 2552 rmsk = 0x7ff; 2553 rmbit = 0x0400; 2554 rebit = 0x800; 2555 re = rw; 2556 break; 2557 2558 /* For C4x arithmetic */ 2559 case 32: 2560 rw = 5; 2561 rmsk = 0xffff; 2562 rmbit = 0x8000; 2563 rebit = 1; 2564 re = rw - 1; 2565 break; 2566 2567 case 24: 2568 rw = 4; 2569 rmsk = 0xff; 2570 rmbit = 0x80; 2571 rebit = 0x100; 2572 re = rw; 2573 break; 2574 } 2575 rbit[re] = rebit; 2576 rlast = rndprc; 2577 } 2578 2579 /* Shift down 1 temporarily if the data structure has an implied 2580 most significant bit and the number is denormal. 2581 Intel long double denormals also lose one bit of precision. */ 2582 if ((exp <= 0) && (rndprc != NBITS) 2583 && ((rndprc != 64) || ((rndprc == 64) && ! REAL_WORDS_BIG_ENDIAN))) 2584 { 2585 lost |= s[NI - 1] & 1; 2586 eshdn1 (s); 2587 } 2588 /* Clear out all bits below the rounding bit, 2589 remembering in r if any were nonzero. */ 2590 r = s[rw] & rmsk; 2591 if (rndprc < NBITS) 2592 { 2593 i = rw + 1; 2594 while (i < NI) 2595 { 2596 if (s[i]) 2597 r |= 1; 2598 s[i] = 0; 2599 ++i; 2600 } 2601 } 2602 s[rw] &= ~rmsk; 2603 if ((r & rmbit) != 0) 2604 { 2605#ifndef C4X 2606 if (r == rmbit) 2607 { 2608 if (lost == 0) 2609 { /* round to even */ 2610 if ((s[re] & rebit) == 0) 2611 goto mddone; 2612 } 2613 else 2614 { 2615 if (subflg != 0) 2616 goto mddone; 2617 } 2618 } 2619#endif 2620 eaddm (rbit, s); 2621 } 2622 mddone: 2623/* Undo the temporary shift for denormal values. */ 2624 if ((exp <= 0) && (rndprc != NBITS) 2625 && ((rndprc != 64) || ((rndprc == 64) && ! REAL_WORDS_BIG_ENDIAN))) 2626 { 2627 eshup1 (s); 2628 } 2629 if (s[2] != 0) 2630 { /* overflow on roundoff */ 2631 eshdn1 (s); 2632 exp += 1; 2633 } 2634 mdfin: 2635 s[NI - 1] = 0; 2636 if (exp >= 32767L) 2637 { 2638#ifndef INFINITY 2639 overf: 2640#endif 2641#ifdef INFINITY 2642 s[1] = 32767; 2643 for (i = 2; i < NI - 1; i++) 2644 s[i] = 0; 2645 if (extra_warnings) 2646 warning ("floating point overflow"); 2647#else 2648 s[1] = 32766; 2649 s[2] = 0; 2650 for (i = M + 1; i < NI - 1; i++) 2651 s[i] = 0xffff; 2652 s[NI - 1] = 0; 2653 if ((rndprc < 64) || (rndprc == 113)) 2654 { 2655 s[rw] &= ~rmsk; 2656 if (rndprc == 24) 2657 { 2658 s[5] = 0; 2659 s[6] = 0; 2660 } 2661 } 2662#endif 2663 return; 2664 } 2665 if (exp < 0) 2666 s[1] = 0; 2667 else 2668 s[1] = (unsigned EMUSHORT) exp; 2669} 2670 2671/* Subtract. C = B - A, all e type numbers. */ 2672 2673static int subflg = 0; 2674 2675static void 2676esub (a, b, c) 2677 unsigned EMUSHORT *a, *b, *c; 2678{ 2679 2680#ifdef NANS 2681 if (eisnan (a)) 2682 { 2683 emov (a, c); 2684 return; 2685 } 2686 if (eisnan (b)) 2687 { 2688 emov (b, c); 2689 return; 2690 } 2691/* Infinity minus infinity is a NaN. 2692 Test for subtracting infinities of the same sign. */ 2693 if (eisinf (a) && eisinf (b) 2694 && ((eisneg (a) ^ eisneg (b)) == 0)) 2695 { 2696 mtherr ("esub", INVALID); 2697 enan (c, 0); 2698 return; 2699 } 2700#endif 2701 subflg = 1; 2702 eadd1 (a, b, c); 2703} 2704 2705/* Add. C = A + B, all e type. */ 2706 2707static void 2708eadd (a, b, c) 2709 unsigned EMUSHORT *a, *b, *c; 2710{ 2711 2712#ifdef NANS 2713/* NaN plus anything is a NaN. */ 2714 if (eisnan (a)) 2715 { 2716 emov (a, c); 2717 return; 2718 } 2719 if (eisnan (b)) 2720 { 2721 emov (b, c); 2722 return; 2723 } 2724/* Infinity minus infinity is a NaN. 2725 Test for adding infinities of opposite signs. */ 2726 if (eisinf (a) && eisinf (b) 2727 && ((eisneg (a) ^ eisneg (b)) != 0)) 2728 { 2729 mtherr ("esub", INVALID); 2730 enan (c, 0); 2731 return; 2732 } 2733#endif 2734 subflg = 0; 2735 eadd1 (a, b, c); 2736} 2737 2738/* Arithmetic common to both addition and subtraction. */ 2739 2740static void 2741eadd1 (a, b, c) 2742 unsigned EMUSHORT *a, *b, *c; 2743{ 2744 unsigned EMUSHORT ai[NI], bi[NI], ci[NI]; 2745 int i, lost, j, k; 2746 EMULONG lt, lta, ltb; 2747 2748#ifdef INFINITY 2749 if (eisinf (a)) 2750 { 2751 emov (a, c); 2752 if (subflg) 2753 eneg (c); 2754 return; 2755 } 2756 if (eisinf (b)) 2757 { 2758 emov (b, c); 2759 return; 2760 } 2761#endif 2762 emovi (a, ai); 2763 emovi (b, bi); 2764 if (subflg) 2765 ai[0] = ~ai[0]; 2766 2767 /* compare exponents */ 2768 lta = ai[E]; 2769 ltb = bi[E]; 2770 lt = lta - ltb; 2771 if (lt > 0L) 2772 { /* put the larger number in bi */ 2773 emovz (bi, ci); 2774 emovz (ai, bi); 2775 emovz (ci, ai); 2776 ltb = bi[E]; 2777 lt = -lt; 2778 } 2779 lost = 0; 2780 if (lt != 0L) 2781 { 2782 if (lt < (EMULONG) (-NBITS - 1)) 2783 goto done; /* answer same as larger addend */ 2784 k = (int) lt; 2785 lost = eshift (ai, k); /* shift the smaller number down */ 2786 } 2787 else 2788 { 2789 /* exponents were the same, so must compare significands */ 2790 i = ecmpm (ai, bi); 2791 if (i == 0) 2792 { /* the numbers are identical in magnitude */ 2793 /* if different signs, result is zero */ 2794 if (ai[0] != bi[0]) 2795 { 2796 eclear (c); 2797 return; 2798 } 2799 /* if same sign, result is double */ 2800 /* double denormalized tiny number */ 2801 if ((bi[E] == 0) && ((bi[3] & 0x8000) == 0)) 2802 { 2803 eshup1 (bi); 2804 goto done; 2805 } 2806 /* add 1 to exponent unless both are zero! */ 2807 for (j = 1; j < NI - 1; j++) 2808 { 2809 if (bi[j] != 0) 2810 { 2811 ltb += 1; 2812 if (ltb >= 0x7fff) 2813 { 2814 eclear (c); 2815 if (ai[0] != 0) 2816 eneg (c); 2817 einfin (c); 2818 return; 2819 } 2820 break; 2821 } 2822 } 2823 bi[E] = (unsigned EMUSHORT) ltb; 2824 goto done; 2825 } 2826 if (i > 0) 2827 { /* put the larger number in bi */ 2828 emovz (bi, ci); 2829 emovz (ai, bi); 2830 emovz (ci, ai); 2831 } 2832 } 2833 if (ai[0] == bi[0]) 2834 { 2835 eaddm (ai, bi); 2836 subflg = 0; 2837 } 2838 else 2839 { 2840 esubm (ai, bi); 2841 subflg = 1; 2842 } 2843 emdnorm (bi, lost, subflg, ltb, 64); 2844 2845 done: 2846 emovo (bi, c); 2847} 2848 2849/* Divide: C = B/A, all e type. */ 2850 2851static void 2852ediv (a, b, c) 2853 unsigned EMUSHORT *a, *b, *c; 2854{ 2855 unsigned EMUSHORT ai[NI], bi[NI]; 2856 int i, sign; 2857 EMULONG lt, lta, ltb; 2858 2859/* IEEE says if result is not a NaN, the sign is "-" if and only if 2860 operands have opposite signs -- but flush -0 to 0 later if not IEEE. */ 2861 sign = eisneg(a) ^ eisneg(b); 2862 2863#ifdef NANS 2864/* Return any NaN input. */ 2865 if (eisnan (a)) 2866 { 2867 emov (a, c); 2868 return; 2869 } 2870 if (eisnan (b)) 2871 { 2872 emov (b, c); 2873 return; 2874 } 2875/* Zero over zero, or infinity over infinity, is a NaN. */ 2876 if (((ecmp (a, ezero) == 0) && (ecmp (b, ezero) == 0)) 2877 || (eisinf (a) && eisinf (b))) 2878 { 2879 mtherr ("ediv", INVALID); 2880 enan (c, sign); 2881 return; 2882 } 2883#endif 2884/* Infinity over anything else is infinity. */ 2885#ifdef INFINITY 2886 if (eisinf (b)) 2887 { 2888 einfin (c); 2889 goto divsign; 2890 } 2891/* Anything else over infinity is zero. */ 2892 if (eisinf (a)) 2893 { 2894 eclear (c); 2895 goto divsign; 2896 } 2897#endif 2898 emovi (a, ai); 2899 emovi (b, bi); 2900 lta = ai[E]; 2901 ltb = bi[E]; 2902 if (bi[E] == 0) 2903 { /* See if numerator is zero. */ 2904 for (i = 1; i < NI - 1; i++) 2905 { 2906 if (bi[i] != 0) 2907 { 2908 ltb -= enormlz (bi); 2909 goto dnzro1; 2910 } 2911 } 2912 eclear (c); 2913 goto divsign; 2914 } 2915 dnzro1: 2916 2917 if (ai[E] == 0) 2918 { /* possible divide by zero */ 2919 for (i = 1; i < NI - 1; i++) 2920 { 2921 if (ai[i] != 0) 2922 { 2923 lta -= enormlz (ai); 2924 goto dnzro2; 2925 } 2926 } 2927/* Divide by zero is not an invalid operation. 2928 It is a divide-by-zero operation! */ 2929 einfin (c); 2930 mtherr ("ediv", SING); 2931 goto divsign; 2932 } 2933 dnzro2: 2934 2935 i = edivm (ai, bi); 2936 /* calculate exponent */ 2937 lt = ltb - lta + EXONE; 2938 emdnorm (bi, i, 0, lt, 64); 2939 emovo (bi, c); 2940 2941 divsign: 2942 2943 if (sign 2944#ifndef IEEE 2945 && (ecmp (c, ezero) != 0) 2946#endif 2947 ) 2948 *(c+(NE-1)) |= 0x8000; 2949 else 2950 *(c+(NE-1)) &= ~0x8000; 2951} 2952 2953/* Multiply e-types A and B, return e-type product C. */ 2954 2955static void 2956emul (a, b, c) 2957 unsigned EMUSHORT *a, *b, *c; 2958{ 2959 unsigned EMUSHORT ai[NI], bi[NI]; 2960 int i, j, sign; 2961 EMULONG lt, lta, ltb; 2962 2963/* IEEE says if result is not a NaN, the sign is "-" if and only if 2964 operands have opposite signs -- but flush -0 to 0 later if not IEEE. */ 2965 sign = eisneg(a) ^ eisneg(b); 2966 2967#ifdef NANS 2968/* NaN times anything is the same NaN. */ 2969 if (eisnan (a)) 2970 { 2971 emov (a, c); 2972 return; 2973 } 2974 if (eisnan (b)) 2975 { 2976 emov (b, c); 2977 return; 2978 } 2979/* Zero times infinity is a NaN. */ 2980 if ((eisinf (a) && (ecmp (b, ezero) == 0)) 2981 || (eisinf (b) && (ecmp (a, ezero) == 0))) 2982 { 2983 mtherr ("emul", INVALID); 2984 enan (c, sign); 2985 return; 2986 } 2987#endif 2988/* Infinity times anything else is infinity. */ 2989#ifdef INFINITY 2990 if (eisinf (a) || eisinf (b)) 2991 { 2992 einfin (c); 2993 goto mulsign; 2994 } 2995#endif 2996 emovi (a, ai); 2997 emovi (b, bi); 2998 lta = ai[E]; 2999 ltb = bi[E]; 3000 if (ai[E] == 0) 3001 { 3002 for (i = 1; i < NI - 1; i++) 3003 { 3004 if (ai[i] != 0) 3005 { 3006 lta -= enormlz (ai); 3007 goto mnzer1; 3008 } 3009 } 3010 eclear (c); 3011 goto mulsign; 3012 } 3013 mnzer1: 3014 3015 if (bi[E] == 0) 3016 { 3017 for (i = 1; i < NI - 1; i++) 3018 { 3019 if (bi[i] != 0) 3020 { 3021 ltb -= enormlz (bi); 3022 goto mnzer2; 3023 } 3024 } 3025 eclear (c); 3026 goto mulsign; 3027 } 3028 mnzer2: 3029 3030 /* Multiply significands */ 3031 j = emulm (ai, bi); 3032 /* calculate exponent */ 3033 lt = lta + ltb - (EXONE - 1); 3034 emdnorm (bi, j, 0, lt, 64); 3035 emovo (bi, c); 3036 3037 mulsign: 3038 3039 if (sign 3040#ifndef IEEE 3041 && (ecmp (c, ezero) != 0) 3042#endif 3043 ) 3044 *(c+(NE-1)) |= 0x8000; 3045 else 3046 *(c+(NE-1)) &= ~0x8000; 3047} 3048 3049/* Convert double precision PE to e-type Y. */ 3050 3051static void 3052e53toe (pe, y) 3053 unsigned EMUSHORT *pe, *y; 3054{ 3055#ifdef DEC 3056 3057 dectoe (pe, y); 3058 3059#else 3060#ifdef IBM 3061 3062 ibmtoe (pe, y, DFmode); 3063 3064#else 3065#ifdef C4X 3066 3067 c4xtoe (pe, y, HFmode); 3068 3069#else 3070 register unsigned EMUSHORT r; 3071 register unsigned EMUSHORT *e, *p; 3072 unsigned EMUSHORT yy[NI]; 3073 int denorm, k; 3074 3075 e = pe; 3076 denorm = 0; /* flag if denormalized number */ 3077 ecleaz (yy); 3078 if (! REAL_WORDS_BIG_ENDIAN) 3079 e += 3; 3080 r = *e; 3081 yy[0] = 0; 3082 if (r & 0x8000) 3083 yy[0] = 0xffff; 3084 yy[M] = (r & 0x0f) | 0x10; 3085 r &= ~0x800f; /* strip sign and 4 significand bits */ 3086#ifdef INFINITY 3087 if (r == 0x7ff0) 3088 { 3089#ifdef NANS 3090 if (! REAL_WORDS_BIG_ENDIAN) 3091 { 3092 if (((pe[3] & 0xf) != 0) || (pe[2] != 0) 3093 || (pe[1] != 0) || (pe[0] != 0)) 3094 { 3095 enan (y, yy[0] != 0); 3096 return; 3097 } 3098 } 3099 else 3100 { 3101 if (((pe[0] & 0xf) != 0) || (pe[1] != 0) 3102 || (pe[2] != 0) || (pe[3] != 0)) 3103 { 3104 enan (y, yy[0] != 0); 3105 return; 3106 } 3107 } 3108#endif /* NANS */ 3109 eclear (y); 3110 einfin (y); 3111 if (yy[0]) 3112 eneg (y); 3113 return; 3114 } 3115#endif /* INFINITY */ 3116 r >>= 4; 3117 /* If zero exponent, then the significand is denormalized. 3118 So take back the understood high significand bit. */ 3119 3120 if (r == 0) 3121 { 3122 denorm = 1; 3123 yy[M] &= ~0x10; 3124 } 3125 r += EXONE - 01777; 3126 yy[E] = r; 3127 p = &yy[M + 1]; 3128#ifdef IEEE 3129 if (! REAL_WORDS_BIG_ENDIAN) 3130 { 3131 *p++ = *(--e); 3132 *p++ = *(--e); 3133 *p++ = *(--e); 3134 } 3135 else 3136 { 3137 ++e; 3138 *p++ = *e++; 3139 *p++ = *e++; 3140 *p++ = *e++; 3141 } 3142#endif 3143 eshift (yy, -5); 3144 if (denorm) 3145 { 3146 /* If zero exponent, then normalize the significand. */ 3147 if ((k = enormlz (yy)) > NBITS) 3148 ecleazs (yy); 3149 else 3150 yy[E] -= (unsigned EMUSHORT) (k - 1); 3151 } 3152 emovo (yy, y); 3153#endif /* not C4X */ 3154#endif /* not IBM */ 3155#endif /* not DEC */ 3156} 3157 3158/* Convert double extended precision float PE to e type Y. */ 3159 3160static void 3161e64toe (pe, y) 3162 unsigned EMUSHORT *pe, *y; 3163{ 3164 unsigned EMUSHORT yy[NI]; 3165 unsigned EMUSHORT *e, *p, *q; 3166 int i; 3167 3168 e = pe; 3169 p = yy; 3170 for (i = 0; i < NE - 5; i++) 3171 *p++ = 0; 3172/* This precision is not ordinarily supported on DEC or IBM. */ 3173#ifdef DEC 3174 for (i = 0; i < 5; i++) 3175 *p++ = *e++; 3176#endif 3177#ifdef IBM 3178 p = &yy[0] + (NE - 1); 3179 *p-- = *e++; 3180 ++e; 3181 for (i = 0; i < 5; i++) 3182 *p-- = *e++; 3183#endif 3184#ifdef IEEE 3185 if (! REAL_WORDS_BIG_ENDIAN) 3186 { 3187 for (i = 0; i < 5; i++) 3188 *p++ = *e++; 3189 3190 /* For denormal long double Intel format, shift significand up one 3191 -- but only if the top significand bit is zero. A top bit of 1 3192 is "pseudodenormal" when the exponent is zero. */ 3193 if((yy[NE-1] & 0x7fff) == 0 && (yy[NE-2] & 0x8000) == 0) 3194 { 3195 unsigned EMUSHORT temp[NI]; 3196 3197 emovi(yy, temp); 3198 eshup1(temp); 3199 emovo(temp,y); 3200 return; 3201 } 3202 } 3203 else 3204 { 3205 p = &yy[0] + (NE - 1); 3206#ifdef ARM_EXTENDED_IEEE_FORMAT 3207 /* For ARMs, the exponent is in the lowest 15 bits of the word. */ 3208 *p-- = (e[0] & 0x8000) | (e[1] & 0x7ffff); 3209 e += 2; 3210#else 3211 *p-- = *e++; 3212 ++e; 3213#endif 3214 for (i = 0; i < 4; i++) 3215 *p-- = *e++; 3216 } 3217#endif 3218#ifdef INFINITY 3219 /* Point to the exponent field and check max exponent cases. */ 3220 p = &yy[NE - 1]; 3221 if ((*p & 0x7fff) == 0x7fff) 3222 { 3223#ifdef NANS 3224 if (! REAL_WORDS_BIG_ENDIAN) 3225 { 3226 for (i = 0; i < 4; i++) 3227 { 3228 if ((i != 3 && pe[i] != 0) 3229 /* Anything but 0x8000 here, including 0, is a NaN. */ 3230 || (i == 3 && pe[i] != 0x8000)) 3231 { 3232 enan (y, (*p & 0x8000) != 0); 3233 return; 3234 } 3235 } 3236 } 3237 else 3238 { 3239#ifdef ARM_EXTENDED_IEEE_FORMAT 3240 for (i = 2; i <= 5; i++) 3241 { 3242 if (pe[i] != 0) 3243 { 3244 enan (y, (*p & 0x8000) != 0); 3245 return; 3246 } 3247 } 3248#else /* not ARM */ 3249 /* In Motorola extended precision format, the most significant 3250 bit of an infinity mantissa could be either 1 or 0. It is 3251 the lower order bits that tell whether the value is a NaN. */ 3252 if ((pe[2] & 0x7fff) != 0) 3253 goto bigend_nan; 3254 3255 for (i = 3; i <= 5; i++) 3256 { 3257 if (pe[i] != 0) 3258 { 3259bigend_nan: 3260 enan (y, (*p & 0x8000) != 0); 3261 return; 3262 } 3263 } 3264#endif /* not ARM */ 3265 } 3266#endif /* NANS */ 3267 eclear (y); 3268 einfin (y); 3269 if (*p & 0x8000) 3270 eneg (y); 3271 return; 3272 } 3273#endif /* INFINITY */ 3274 p = yy; 3275 q = y; 3276 for (i = 0; i < NE; i++) 3277 *q++ = *p++; 3278} 3279 3280/* Convert 128-bit long double precision float PE to e type Y. */ 3281 3282static void 3283e113toe (pe, y) 3284 unsigned EMUSHORT *pe, *y; 3285{ 3286 register unsigned EMUSHORT r; 3287 unsigned EMUSHORT *e, *p; 3288 unsigned EMUSHORT yy[NI]; 3289 int denorm, i; 3290 3291 e = pe; 3292 denorm = 0; 3293 ecleaz (yy); 3294#ifdef IEEE 3295 if (! REAL_WORDS_BIG_ENDIAN) 3296 e += 7; 3297#endif 3298 r = *e; 3299 yy[0] = 0; 3300 if (r & 0x8000) 3301 yy[0] = 0xffff; 3302 r &= 0x7fff; 3303#ifdef INFINITY 3304 if (r == 0x7fff) 3305 { 3306#ifdef NANS 3307 if (! REAL_WORDS_BIG_ENDIAN) 3308 { 3309 for (i = 0; i < 7; i++) 3310 { 3311 if (pe[i] != 0) 3312 { 3313 enan (y, yy[0] != 0); 3314 return; 3315 } 3316 } 3317 } 3318 else 3319 { 3320 for (i = 1; i < 8; i++) 3321 { 3322 if (pe[i] != 0) 3323 { 3324 enan (y, yy[0] != 0); 3325 return; 3326 } 3327 } 3328 } 3329#endif /* NANS */ 3330 eclear (y); 3331 einfin (y); 3332 if (yy[0]) 3333 eneg (y); 3334 return; 3335 } 3336#endif /* INFINITY */ 3337 yy[E] = r; 3338 p = &yy[M + 1]; 3339#ifdef IEEE 3340 if (! REAL_WORDS_BIG_ENDIAN) 3341 { 3342 for (i = 0; i < 7; i++) 3343 *p++ = *(--e); 3344 } 3345 else 3346 { 3347 ++e; 3348 for (i = 0; i < 7; i++) 3349 *p++ = *e++; 3350 } 3351#endif 3352/* If denormal, remove the implied bit; else shift down 1. */ 3353 if (r == 0) 3354 { 3355 yy[M] = 0; 3356 } 3357 else 3358 { 3359 yy[M] = 1; 3360 eshift (yy, -1); 3361 } 3362 emovo (yy, y); 3363} 3364 3365/* Convert single precision float PE to e type Y. */ 3366 3367static void 3368e24toe (pe, y) 3369 unsigned EMUSHORT *pe, *y; 3370{ 3371#ifdef IBM 3372 3373 ibmtoe (pe, y, SFmode); 3374 3375#else 3376 3377#ifdef C4X 3378 3379 c4xtoe (pe, y, QFmode); 3380 3381#else 3382 3383 register unsigned EMUSHORT r; 3384 register unsigned EMUSHORT *e, *p; 3385 unsigned EMUSHORT yy[NI]; 3386 int denorm, k; 3387 3388 e = pe; 3389 denorm = 0; /* flag if denormalized number */ 3390 ecleaz (yy); 3391#ifdef IEEE 3392 if (! REAL_WORDS_BIG_ENDIAN) 3393 e += 1; 3394#endif 3395#ifdef DEC 3396 e += 1; 3397#endif 3398 r = *e; 3399 yy[0] = 0; 3400 if (r & 0x8000) 3401 yy[0] = 0xffff; 3402 yy[M] = (r & 0x7f) | 0200; 3403 r &= ~0x807f; /* strip sign and 7 significand bits */ 3404#ifdef INFINITY 3405 if (r == 0x7f80) 3406 { 3407#ifdef NANS 3408 if (REAL_WORDS_BIG_ENDIAN) 3409 { 3410 if (((pe[0] & 0x7f) != 0) || (pe[1] != 0)) 3411 { 3412 enan (y, yy[0] != 0); 3413 return; 3414 } 3415 } 3416 else 3417 { 3418 if (((pe[1] & 0x7f) != 0) || (pe[0] != 0)) 3419 { 3420 enan (y, yy[0] != 0); 3421 return; 3422 } 3423 } 3424#endif /* NANS */ 3425 eclear (y); 3426 einfin (y); 3427 if (yy[0]) 3428 eneg (y); 3429 return; 3430 } 3431#endif /* INFINITY */ 3432 r >>= 7; 3433 /* If zero exponent, then the significand is denormalized. 3434 So take back the understood high significand bit. */ 3435 if (r == 0) 3436 { 3437 denorm = 1; 3438 yy[M] &= ~0200; 3439 } 3440 r += EXONE - 0177; 3441 yy[E] = r; 3442 p = &yy[M + 1]; 3443#ifdef DEC 3444 *p++ = *(--e); 3445#endif 3446#ifdef IEEE 3447 if (! REAL_WORDS_BIG_ENDIAN) 3448 *p++ = *(--e); 3449 else 3450 { 3451 ++e; 3452 *p++ = *e++; 3453 } 3454#endif 3455 eshift (yy, -8); 3456 if (denorm) 3457 { /* if zero exponent, then normalize the significand */ 3458 if ((k = enormlz (yy)) > NBITS) 3459 ecleazs (yy); 3460 else 3461 yy[E] -= (unsigned EMUSHORT) (k - 1); 3462 } 3463 emovo (yy, y); 3464#endif /* not C4X */ 3465#endif /* not IBM */ 3466} 3467 3468/* Convert e-type X to IEEE 128-bit long double format E. */ 3469 3470static void 3471etoe113 (x, e) 3472 unsigned EMUSHORT *x, *e; 3473{ 3474 unsigned EMUSHORT xi[NI]; 3475 EMULONG exp; 3476 int rndsav; 3477 3478#ifdef NANS 3479 if (eisnan (x)) 3480 { 3481 make_nan (e, eisneg (x), TFmode); 3482 return; 3483 } 3484#endif 3485 emovi (x, xi); 3486 exp = (EMULONG) xi[E]; 3487#ifdef INFINITY 3488 if (eisinf (x)) 3489 goto nonorm; 3490#endif 3491 /* round off to nearest or even */ 3492 rndsav = rndprc; 3493 rndprc = 113; 3494 emdnorm (xi, 0, 0, exp, 64); 3495 rndprc = rndsav; 3496 nonorm: 3497 toe113 (xi, e); 3498} 3499 3500/* Convert exploded e-type X, that has already been rounded to 3501 113-bit precision, to IEEE 128-bit long double format Y. */ 3502 3503static void 3504toe113 (a, b) 3505 unsigned EMUSHORT *a, *b; 3506{ 3507 register unsigned EMUSHORT *p, *q; 3508 unsigned EMUSHORT i; 3509 3510#ifdef NANS 3511 if (eiisnan (a)) 3512 { 3513 make_nan (b, eiisneg (a), TFmode); 3514 return; 3515 } 3516#endif 3517 p = a; 3518 if (REAL_WORDS_BIG_ENDIAN) 3519 q = b; 3520 else 3521 q = b + 7; /* point to output exponent */ 3522 3523 /* If not denormal, delete the implied bit. */ 3524 if (a[E] != 0) 3525 { 3526 eshup1 (a); 3527 } 3528 /* combine sign and exponent */ 3529 i = *p++; 3530 if (REAL_WORDS_BIG_ENDIAN) 3531 { 3532 if (i) 3533 *q++ = *p++ | 0x8000; 3534 else 3535 *q++ = *p++; 3536 } 3537 else 3538 { 3539 if (i) 3540 *q-- = *p++ | 0x8000; 3541 else 3542 *q-- = *p++; 3543 } 3544 /* skip over guard word */ 3545 ++p; 3546 /* move the significand */ 3547 if (REAL_WORDS_BIG_ENDIAN) 3548 { 3549 for (i = 0; i < 7; i++) 3550 *q++ = *p++; 3551 } 3552 else 3553 { 3554 for (i = 0; i < 7; i++) 3555 *q-- = *p++; 3556 } 3557} 3558 3559/* Convert e-type X to IEEE double extended format E. */ 3560 3561static void 3562etoe64 (x, e) 3563 unsigned EMUSHORT *x, *e; 3564{ 3565 unsigned EMUSHORT xi[NI]; 3566 EMULONG exp; 3567 int rndsav; 3568 3569#ifdef NANS 3570 if (eisnan (x)) 3571 { 3572 make_nan (e, eisneg (x), XFmode); 3573 return; 3574 } 3575#endif 3576 emovi (x, xi); 3577 /* adjust exponent for offset */ 3578 exp = (EMULONG) xi[E]; 3579#ifdef INFINITY 3580 if (eisinf (x)) 3581 goto nonorm; 3582#endif 3583 /* round off to nearest or even */ 3584 rndsav = rndprc; 3585 rndprc = 64; 3586 emdnorm (xi, 0, 0, exp, 64); 3587 rndprc = rndsav; 3588 nonorm: 3589 toe64 (xi, e); 3590} 3591 3592/* Convert exploded e-type X, that has already been rounded to 3593 64-bit precision, to IEEE double extended format Y. */ 3594 3595static void 3596toe64 (a, b) 3597 unsigned EMUSHORT *a, *b; 3598{ 3599 register unsigned EMUSHORT *p, *q; 3600 unsigned EMUSHORT i; 3601 3602#ifdef NANS 3603 if (eiisnan (a)) 3604 { 3605 make_nan (b, eiisneg (a), XFmode); 3606 return; 3607 } 3608#endif 3609 /* Shift denormal long double Intel format significand down one bit. */ 3610 if ((a[E] == 0) && ! REAL_WORDS_BIG_ENDIAN) 3611 eshdn1 (a); 3612 p = a; 3613#ifdef IBM 3614 q = b; 3615#endif 3616#ifdef DEC 3617 q = b + 4; 3618#endif 3619#ifdef IEEE 3620 if (REAL_WORDS_BIG_ENDIAN) 3621 q = b; 3622 else 3623 { 3624 q = b + 4; /* point to output exponent */ 3625#if LONG_DOUBLE_TYPE_SIZE == 96 3626 /* Clear the last two bytes of 12-byte Intel format */ 3627 *(q+1) = 0; 3628#endif 3629 } 3630#endif 3631 3632 /* combine sign and exponent */ 3633 i = *p++; 3634#ifdef IBM 3635 if (i) 3636 *q++ = *p++ | 0x8000; 3637 else 3638 *q++ = *p++; 3639 *q++ = 0; 3640#endif 3641#ifdef DEC 3642 if (i) 3643 *q-- = *p++ | 0x8000; 3644 else 3645 *q-- = *p++; 3646#endif 3647#ifdef IEEE 3648 if (REAL_WORDS_BIG_ENDIAN) 3649 { 3650#ifdef ARM_EXTENDED_IEEE_FORMAT 3651 /* The exponent is in the lowest 15 bits of the first word. */ 3652 *q++ = i ? 0x8000 : 0; 3653 *q++ = *p++; 3654#else 3655 if (i) 3656 *q++ = *p++ | 0x8000; 3657 else 3658 *q++ = *p++; 3659 *q++ = 0; 3660#endif 3661 } 3662 else 3663 { 3664 if (i) 3665 *q-- = *p++ | 0x8000; 3666 else 3667 *q-- = *p++; 3668 } 3669#endif 3670 /* skip over guard word */ 3671 ++p; 3672 /* move the significand */ 3673#ifdef IBM 3674 for (i = 0; i < 4; i++) 3675 *q++ = *p++; 3676#endif 3677#ifdef DEC 3678 for (i = 0; i < 4; i++) 3679 *q-- = *p++; 3680#endif 3681#ifdef IEEE 3682 if (REAL_WORDS_BIG_ENDIAN) 3683 { 3684 for (i = 0; i < 4; i++) 3685 *q++ = *p++; 3686 } 3687 else 3688 { 3689#ifdef INFINITY 3690 if (eiisinf (a)) 3691 { 3692 /* Intel long double infinity significand. */ 3693 *q-- = 0x8000; 3694 *q-- = 0; 3695 *q-- = 0; 3696 *q = 0; 3697 return; 3698 } 3699#endif 3700 for (i = 0; i < 4; i++) 3701 *q-- = *p++; 3702 } 3703#endif 3704} 3705 3706/* e type to double precision. */ 3707 3708#ifdef DEC 3709/* Convert e-type X to DEC-format double E. */ 3710 3711static void 3712etoe53 (x, e) 3713 unsigned EMUSHORT *x, *e; 3714{ 3715 etodec (x, e); /* see etodec.c */ 3716} 3717 3718/* Convert exploded e-type X, that has already been rounded to 3719 56-bit double precision, to DEC double Y. */ 3720 3721static void 3722toe53 (x, y) 3723 unsigned EMUSHORT *x, *y; 3724{ 3725 todec (x, y); 3726} 3727 3728#else 3729#ifdef IBM 3730/* Convert e-type X to IBM 370-format double E. */ 3731 3732static void 3733etoe53 (x, e) 3734 unsigned EMUSHORT *x, *e; 3735{ 3736 etoibm (x, e, DFmode); 3737} 3738 3739/* Convert exploded e-type X, that has already been rounded to 3740 56-bit precision, to IBM 370 double Y. */ 3741 3742static void 3743toe53 (x, y) 3744 unsigned EMUSHORT *x, *y; 3745{ 3746 toibm (x, y, DFmode); 3747} 3748 3749#else /* it's neither DEC nor IBM */ 3750#ifdef C4X 3751/* Convert e-type X to C4X-format long double E. */ 3752 3753static void 3754etoe53 (x, e) 3755 unsigned EMUSHORT *x, *e; 3756{ 3757 etoc4x (x, e, HFmode); 3758} 3759 3760/* Convert exploded e-type X, that has already been rounded to 3761 56-bit precision, to IBM 370 double Y. */ 3762 3763static void 3764toe53 (x, y) 3765 unsigned EMUSHORT *x, *y; 3766{ 3767 toc4x (x, y, HFmode); 3768} 3769 3770#else /* it's neither DEC nor IBM nor C4X */ 3771 3772/* Convert e-type X to IEEE double E. */ 3773 3774static void 3775etoe53 (x, e) 3776 unsigned EMUSHORT *x, *e; 3777{ 3778 unsigned EMUSHORT xi[NI]; 3779 EMULONG exp; 3780 int rndsav; 3781 3782#ifdef NANS 3783 if (eisnan (x)) 3784 { 3785 make_nan (e, eisneg (x), DFmode); 3786 return; 3787 } 3788#endif 3789 emovi (x, xi); 3790 /* adjust exponent for offsets */ 3791 exp = (EMULONG) xi[E] - (EXONE - 0x3ff); 3792#ifdef INFINITY 3793 if (eisinf (x)) 3794 goto nonorm; 3795#endif 3796 /* round off to nearest or even */ 3797 rndsav = rndprc; 3798 rndprc = 53; 3799 emdnorm (xi, 0, 0, exp, 64); 3800 rndprc = rndsav; 3801 nonorm: 3802 toe53 (xi, e); 3803} 3804 3805/* Convert exploded e-type X, that has already been rounded to 3806 53-bit precision, to IEEE double Y. */ 3807 3808static void 3809toe53 (x, y) 3810 unsigned EMUSHORT *x, *y; 3811{ 3812 unsigned EMUSHORT i; 3813 unsigned EMUSHORT *p; 3814 3815#ifdef NANS 3816 if (eiisnan (x)) 3817 { 3818 make_nan (y, eiisneg (x), DFmode); 3819 return; 3820 } 3821#endif 3822 p = &x[0]; 3823#ifdef IEEE 3824 if (! REAL_WORDS_BIG_ENDIAN) 3825 y += 3; 3826#endif 3827 *y = 0; /* output high order */ 3828 if (*p++) 3829 *y = 0x8000; /* output sign bit */ 3830 3831 i = *p++; 3832 if (i >= (unsigned int) 2047) 3833 { 3834 /* Saturate at largest number less than infinity. */ 3835#ifdef INFINITY 3836 *y |= 0x7ff0; 3837 if (! REAL_WORDS_BIG_ENDIAN) 3838 { 3839 *(--y) = 0; 3840 *(--y) = 0; 3841 *(--y) = 0; 3842 } 3843 else 3844 { 3845 ++y; 3846 *y++ = 0; 3847 *y++ = 0; 3848 *y++ = 0; 3849 } 3850#else 3851 *y |= (unsigned EMUSHORT) 0x7fef; 3852 if (! REAL_WORDS_BIG_ENDIAN) 3853 { 3854 *(--y) = 0xffff; 3855 *(--y) = 0xffff; 3856 *(--y) = 0xffff; 3857 } 3858 else 3859 { 3860 ++y; 3861 *y++ = 0xffff; 3862 *y++ = 0xffff; 3863 *y++ = 0xffff; 3864 } 3865#endif 3866 return; 3867 } 3868 if (i == 0) 3869 { 3870 eshift (x, 4); 3871 } 3872 else 3873 { 3874 i <<= 4; 3875 eshift (x, 5); 3876 } 3877 i |= *p++ & (unsigned EMUSHORT) 0x0f; /* *p = xi[M] */ 3878 *y |= (unsigned EMUSHORT) i; /* high order output already has sign bit set */ 3879 if (! REAL_WORDS_BIG_ENDIAN) 3880 { 3881 *(--y) = *p++; 3882 *(--y) = *p++; 3883 *(--y) = *p; 3884 } 3885 else 3886 { 3887 ++y; 3888 *y++ = *p++; 3889 *y++ = *p++; 3890 *y++ = *p++; 3891 } 3892} 3893 3894#endif /* not C4X */ 3895#endif /* not IBM */ 3896#endif /* not DEC */ 3897 3898 3899 3900/* e type to single precision. */ 3901 3902#ifdef IBM 3903/* Convert e-type X to IBM 370 float E. */ 3904 3905static void 3906etoe24 (x, e) 3907 unsigned EMUSHORT *x, *e; 3908{ 3909 etoibm (x, e, SFmode); 3910} 3911 3912/* Convert exploded e-type X, that has already been rounded to 3913 float precision, to IBM 370 float Y. */ 3914 3915static void 3916toe24 (x, y) 3917 unsigned EMUSHORT *x, *y; 3918{ 3919 toibm (x, y, SFmode); 3920} 3921 3922#else 3923 3924#ifdef C4X 3925/* Convert e-type X to C4X float E. */ 3926 3927static void 3928etoe24 (x, e) 3929 unsigned EMUSHORT *x, *e; 3930{ 3931 etoc4x (x, e, QFmode); 3932} 3933 3934/* Convert exploded e-type X, that has already been rounded to 3935 float precision, to IBM 370 float Y. */ 3936 3937static void 3938toe24 (x, y) 3939 unsigned EMUSHORT *x, *y; 3940{ 3941 toc4x (x, y, QFmode); 3942} 3943 3944#else 3945 3946/* Convert e-type X to IEEE float E. DEC float is the same as IEEE float. */ 3947 3948static void 3949etoe24 (x, e) 3950 unsigned EMUSHORT *x, *e; 3951{ 3952 EMULONG exp; 3953 unsigned EMUSHORT xi[NI]; 3954 int rndsav; 3955 3956#ifdef NANS 3957 if (eisnan (x)) 3958 { 3959 make_nan (e, eisneg (x), SFmode); 3960 return; 3961 } 3962#endif 3963 emovi (x, xi); 3964 /* adjust exponent for offsets */ 3965 exp = (EMULONG) xi[E] - (EXONE - 0177); 3966#ifdef INFINITY 3967 if (eisinf (x)) 3968 goto nonorm; 3969#endif 3970 /* round off to nearest or even */ 3971 rndsav = rndprc; 3972 rndprc = 24; 3973 emdnorm (xi, 0, 0, exp, 64); 3974 rndprc = rndsav; 3975 nonorm: 3976 toe24 (xi, e); 3977} 3978 3979/* Convert exploded e-type X, that has already been rounded to 3980 float precision, to IEEE float Y. */ 3981 3982static void 3983toe24 (x, y) 3984 unsigned EMUSHORT *x, *y; 3985{ 3986 unsigned EMUSHORT i; 3987 unsigned EMUSHORT *p; 3988 3989#ifdef NANS 3990 if (eiisnan (x)) 3991 { 3992 make_nan (y, eiisneg (x), SFmode); 3993 return; 3994 } 3995#endif 3996 p = &x[0]; 3997#ifdef IEEE 3998 if (! REAL_WORDS_BIG_ENDIAN) 3999 y += 1; 4000#endif 4001#ifdef DEC 4002 y += 1; 4003#endif 4004 *y = 0; /* output high order */ 4005 if (*p++) 4006 *y = 0x8000; /* output sign bit */ 4007 4008 i = *p++; 4009/* Handle overflow cases. */ 4010 if (i >= 255) 4011 { 4012#ifdef INFINITY 4013 *y |= (unsigned EMUSHORT) 0x7f80; 4014#ifdef DEC 4015 *(--y) = 0; 4016#endif 4017#ifdef IEEE 4018 if (! REAL_WORDS_BIG_ENDIAN) 4019 *(--y) = 0; 4020 else 4021 { 4022 ++y; 4023 *y = 0; 4024 } 4025#endif 4026#else /* no INFINITY */ 4027 *y |= (unsigned EMUSHORT) 0x7f7f; 4028#ifdef DEC 4029 *(--y) = 0xffff; 4030#endif 4031#ifdef IEEE 4032 if (! REAL_WORDS_BIG_ENDIAN) 4033 *(--y) = 0xffff; 4034 else 4035 { 4036 ++y; 4037 *y = 0xffff; 4038 } 4039#endif 4040#ifdef ERANGE 4041 errno = ERANGE; 4042#endif 4043#endif /* no INFINITY */ 4044 return; 4045 } 4046 if (i == 0) 4047 { 4048 eshift (x, 7); 4049 } 4050 else 4051 { 4052 i <<= 7; 4053 eshift (x, 8); 4054 } 4055 i |= *p++ & (unsigned EMUSHORT) 0x7f; /* *p = xi[M] */ 4056 /* High order output already has sign bit set. */ 4057 *y |= i; 4058#ifdef DEC 4059 *(--y) = *p; 4060#endif 4061#ifdef IEEE 4062 if (! REAL_WORDS_BIG_ENDIAN) 4063 *(--y) = *p; 4064 else 4065 { 4066 ++y; 4067 *y = *p; 4068 } 4069#endif 4070} 4071#endif /* not C4X */ 4072#endif /* not IBM */ 4073 4074/* Compare two e type numbers. 4075 Return +1 if a > b 4076 0 if a == b 4077 -1 if a < b 4078 -2 if either a or b is a NaN. */ 4079 4080static int 4081ecmp (a, b) 4082 unsigned EMUSHORT *a, *b; 4083{ 4084 unsigned EMUSHORT ai[NI], bi[NI]; 4085 register unsigned EMUSHORT *p, *q; 4086 register int i; 4087 int msign; 4088 4089#ifdef NANS 4090 if (eisnan (a) || eisnan (b)) 4091 return (-2); 4092#endif 4093 emovi (a, ai); 4094 p = ai; 4095 emovi (b, bi); 4096 q = bi; 4097 4098 if (*p != *q) 4099 { /* the signs are different */ 4100 /* -0 equals + 0 */ 4101 for (i = 1; i < NI - 1; i++) 4102 { 4103 if (ai[i] != 0) 4104 goto nzro; 4105 if (bi[i] != 0) 4106 goto nzro; 4107 } 4108 return (0); 4109 nzro: 4110 if (*p == 0) 4111 return (1); 4112 else 4113 return (-1); 4114 } 4115 /* both are the same sign */ 4116 if (*p == 0) 4117 msign = 1; 4118 else 4119 msign = -1; 4120 i = NI - 1; 4121 do 4122 { 4123 if (*p++ != *q++) 4124 { 4125 goto diff; 4126 } 4127 } 4128 while (--i > 0); 4129 4130 return (0); /* equality */ 4131 4132 diff: 4133 4134 if (*(--p) > *(--q)) 4135 return (msign); /* p is bigger */ 4136 else 4137 return (-msign); /* p is littler */ 4138} 4139 4140#if 0 4141/* Find e-type nearest integer to X, as floor (X + 0.5). */ 4142 4143static void 4144eround (x, y) 4145 unsigned EMUSHORT *x, *y; 4146{ 4147 eadd (ehalf, x, y); 4148 efloor (y, y); 4149} 4150#endif /* 0 */ 4151 4152/* Convert HOST_WIDE_INT LP to e type Y. */ 4153 4154static void 4155ltoe (lp, y) 4156 HOST_WIDE_INT *lp; 4157 unsigned EMUSHORT *y; 4158{ 4159 unsigned EMUSHORT yi[NI]; 4160 unsigned HOST_WIDE_INT ll; 4161 int k; 4162 4163 ecleaz (yi); 4164 if (*lp < 0) 4165 { 4166 /* make it positive */ 4167 ll = (unsigned HOST_WIDE_INT) (-(*lp)); 4168 yi[0] = 0xffff; /* put correct sign in the e type number */ 4169 } 4170 else 4171 { 4172 ll = (unsigned HOST_WIDE_INT) (*lp); 4173 } 4174 /* move the long integer to yi significand area */ 4175#if HOST_BITS_PER_WIDE_INT == 64 4176 yi[M] = (unsigned EMUSHORT) (ll >> 48); 4177 yi[M + 1] = (unsigned EMUSHORT) (ll >> 32); 4178 yi[M + 2] = (unsigned EMUSHORT) (ll >> 16); 4179 yi[M + 3] = (unsigned EMUSHORT) ll; 4180 yi[E] = EXONE + 47; /* exponent if normalize shift count were 0 */ 4181#else 4182 yi[M] = (unsigned EMUSHORT) (ll >> 16); 4183 yi[M + 1] = (unsigned EMUSHORT) ll; 4184 yi[E] = EXONE + 15; /* exponent if normalize shift count were 0 */ 4185#endif 4186 4187 if ((k = enormlz (yi)) > NBITS)/* normalize the significand */ 4188 ecleaz (yi); /* it was zero */ 4189 else 4190 yi[E] -= (unsigned EMUSHORT) k;/* subtract shift count from exponent */ 4191 emovo (yi, y); /* output the answer */ 4192} 4193 4194/* Convert unsigned HOST_WIDE_INT LP to e type Y. */ 4195 4196static void 4197ultoe (lp, y) 4198 unsigned HOST_WIDE_INT *lp; 4199 unsigned EMUSHORT *y; 4200{ 4201 unsigned EMUSHORT yi[NI]; 4202 unsigned HOST_WIDE_INT ll; 4203 int k; 4204 4205 ecleaz (yi); 4206 ll = *lp; 4207 4208 /* move the long integer to ayi significand area */ 4209#if HOST_BITS_PER_WIDE_INT == 64 4210 yi[M] = (unsigned EMUSHORT) (ll >> 48); 4211 yi[M + 1] = (unsigned EMUSHORT) (ll >> 32); 4212 yi[M + 2] = (unsigned EMUSHORT) (ll >> 16); 4213 yi[M + 3] = (unsigned EMUSHORT) ll; 4214 yi[E] = EXONE + 47; /* exponent if normalize shift count were 0 */ 4215#else 4216 yi[M] = (unsigned EMUSHORT) (ll >> 16); 4217 yi[M + 1] = (unsigned EMUSHORT) ll; 4218 yi[E] = EXONE + 15; /* exponent if normalize shift count were 0 */ 4219#endif 4220 4221 if ((k = enormlz (yi)) > NBITS)/* normalize the significand */ 4222 ecleaz (yi); /* it was zero */ 4223 else 4224 yi[E] -= (unsigned EMUSHORT) k; /* subtract shift count from exponent */ 4225 emovo (yi, y); /* output the answer */ 4226} 4227 4228 4229/* Find signed HOST_WIDE_INT integer I and floating point fractional 4230 part FRAC of e-type (packed internal format) floating point input X. 4231 The integer output I has the sign of the input, except that 4232 positive overflow is permitted if FIXUNS_TRUNC_LIKE_FIX_TRUNC. 4233 The output e-type fraction FRAC is the positive fractional 4234 part of abs (X). */ 4235 4236static void 4237eifrac (x, i, frac) 4238 unsigned EMUSHORT *x; 4239 HOST_WIDE_INT *i; 4240 unsigned EMUSHORT *frac; 4241{ 4242 unsigned EMUSHORT xi[NI]; 4243 int j, k; 4244 unsigned HOST_WIDE_INT ll; 4245 4246 emovi (x, xi); 4247 k = (int) xi[E] - (EXONE - 1); 4248 if (k <= 0) 4249 { 4250 /* if exponent <= 0, integer = 0 and real output is fraction */ 4251 *i = 0L; 4252 emovo (xi, frac); 4253 return; 4254 } 4255 if (k > (HOST_BITS_PER_WIDE_INT - 1)) 4256 { 4257 /* long integer overflow: output large integer 4258 and correct fraction */ 4259 if (xi[0]) 4260 *i = ((unsigned HOST_WIDE_INT) 1) << (HOST_BITS_PER_WIDE_INT - 1); 4261 else 4262 { 4263#ifdef FIXUNS_TRUNC_LIKE_FIX_TRUNC 4264 /* In this case, let it overflow and convert as if unsigned. */ 4265 euifrac (x, &ll, frac); 4266 *i = (HOST_WIDE_INT) ll; 4267 return; 4268#else 4269 /* In other cases, return the largest positive integer. */ 4270 *i = (((unsigned HOST_WIDE_INT) 1) << (HOST_BITS_PER_WIDE_INT - 1)) - 1; 4271#endif 4272 } 4273 eshift (xi, k); 4274 if (extra_warnings) 4275 warning ("overflow on truncation to integer"); 4276 } 4277 else if (k > 16) 4278 { 4279 /* Shift more than 16 bits: first shift up k-16 mod 16, 4280 then shift up by 16's. */ 4281 j = k - ((k >> 4) << 4); 4282 eshift (xi, j); 4283 ll = xi[M]; 4284 k -= j; 4285 do 4286 { 4287 eshup6 (xi); 4288 ll = (ll << 16) | xi[M]; 4289 } 4290 while ((k -= 16) > 0); 4291 *i = ll; 4292 if (xi[0]) 4293 *i = -(*i); 4294 } 4295 else 4296 { 4297 /* shift not more than 16 bits */ 4298 eshift (xi, k); 4299 *i = (HOST_WIDE_INT) xi[M] & 0xffff; 4300 if (xi[0]) 4301 *i = -(*i); 4302 } 4303 xi[0] = 0; 4304 xi[E] = EXONE - 1; 4305 xi[M] = 0; 4306 if ((k = enormlz (xi)) > NBITS) 4307 ecleaz (xi); 4308 else 4309 xi[E] -= (unsigned EMUSHORT) k; 4310 4311 emovo (xi, frac); 4312} 4313 4314 4315/* Find unsigned HOST_WIDE_INT integer I and floating point fractional part 4316 FRAC of e-type X. A negative input yields integer output = 0 but 4317 correct fraction. */ 4318 4319static void 4320euifrac (x, i, frac) 4321 unsigned EMUSHORT *x; 4322 unsigned HOST_WIDE_INT *i; 4323 unsigned EMUSHORT *frac; 4324{ 4325 unsigned HOST_WIDE_INT ll; 4326 unsigned EMUSHORT xi[NI]; 4327 int j, k; 4328 4329 emovi (x, xi); 4330 k = (int) xi[E] - (EXONE - 1); 4331 if (k <= 0) 4332 { 4333 /* if exponent <= 0, integer = 0 and argument is fraction */ 4334 *i = 0L; 4335 emovo (xi, frac); 4336 return; 4337 } 4338 if (k > HOST_BITS_PER_WIDE_INT) 4339 { 4340 /* Long integer overflow: output large integer 4341 and correct fraction. 4342 Note, the BSD microvax compiler says that ~(0UL) 4343 is a syntax error. */ 4344 *i = ~(0L); 4345 eshift (xi, k); 4346 if (extra_warnings) 4347 warning ("overflow on truncation to unsigned integer"); 4348 } 4349 else if (k > 16) 4350 { 4351 /* Shift more than 16 bits: first shift up k-16 mod 16, 4352 then shift up by 16's. */ 4353 j = k - ((k >> 4) << 4); 4354 eshift (xi, j); 4355 ll = xi[M]; 4356 k -= j; 4357 do 4358 { 4359 eshup6 (xi); 4360 ll = (ll << 16) | xi[M]; 4361 } 4362 while ((k -= 16) > 0); 4363 *i = ll; 4364 } 4365 else 4366 { 4367 /* shift not more than 16 bits */ 4368 eshift (xi, k); 4369 *i = (HOST_WIDE_INT) xi[M] & 0xffff; 4370 } 4371 4372 if (xi[0]) /* A negative value yields unsigned integer 0. */ 4373 *i = 0L; 4374 4375 xi[0] = 0; 4376 xi[E] = EXONE - 1; 4377 xi[M] = 0; 4378 if ((k = enormlz (xi)) > NBITS) 4379 ecleaz (xi); 4380 else 4381 xi[E] -= (unsigned EMUSHORT) k; 4382 4383 emovo (xi, frac); 4384} 4385 4386/* Shift the significand of exploded e-type X up or down by SC bits. */ 4387 4388static int 4389eshift (x, sc) 4390 unsigned EMUSHORT *x; 4391 int sc; 4392{ 4393 unsigned EMUSHORT lost; 4394 unsigned EMUSHORT *p; 4395 4396 if (sc == 0) 4397 return (0); 4398 4399 lost = 0; 4400 p = x + NI - 1; 4401 4402 if (sc < 0) 4403 { 4404 sc = -sc; 4405 while (sc >= 16) 4406 { 4407 lost |= *p; /* remember lost bits */ 4408 eshdn6 (x); 4409 sc -= 16; 4410 } 4411 4412 while (sc >= 8) 4413 { 4414 lost |= *p & 0xff; 4415 eshdn8 (x); 4416 sc -= 8; 4417 } 4418 4419 while (sc > 0) 4420 { 4421 lost |= *p & 1; 4422 eshdn1 (x); 4423 sc -= 1; 4424 } 4425 } 4426 else 4427 { 4428 while (sc >= 16) 4429 { 4430 eshup6 (x); 4431 sc -= 16; 4432 } 4433 4434 while (sc >= 8) 4435 { 4436 eshup8 (x); 4437 sc -= 8; 4438 } 4439 4440 while (sc > 0) 4441 { 4442 eshup1 (x); 4443 sc -= 1; 4444 } 4445 } 4446 if (lost) 4447 lost = 1; 4448 return ((int) lost); 4449} 4450 4451/* Shift normalize the significand area of exploded e-type X. 4452 Return the shift count (up = positive). */ 4453 4454static int 4455enormlz (x) 4456 unsigned EMUSHORT x[]; 4457{ 4458 register unsigned EMUSHORT *p; 4459 int sc; 4460 4461 sc = 0; 4462 p = &x[M]; 4463 if (*p != 0) 4464 goto normdn; 4465 ++p; 4466 if (*p & 0x8000) 4467 return (0); /* already normalized */ 4468 while (*p == 0) 4469 { 4470 eshup6 (x); 4471 sc += 16; 4472 4473 /* With guard word, there are NBITS+16 bits available. 4474 Return true if all are zero. */ 4475 if (sc > NBITS) 4476 return (sc); 4477 } 4478 /* see if high byte is zero */ 4479 while ((*p & 0xff00) == 0) 4480 { 4481 eshup8 (x); 4482 sc += 8; 4483 } 4484 /* now shift 1 bit at a time */ 4485 while ((*p & 0x8000) == 0) 4486 { 4487 eshup1 (x); 4488 sc += 1; 4489 if (sc > NBITS) 4490 { 4491 mtherr ("enormlz", UNDERFLOW); 4492 return (sc); 4493 } 4494 } 4495 return (sc); 4496 4497 /* Normalize by shifting down out of the high guard word 4498 of the significand */ 4499 normdn: 4500 4501 if (*p & 0xff00) 4502 { 4503 eshdn8 (x); 4504 sc -= 8; 4505 } 4506 while (*p != 0) 4507 { 4508 eshdn1 (x); 4509 sc -= 1; 4510 4511 if (sc < -NBITS) 4512 { 4513 mtherr ("enormlz", OVERFLOW); 4514 return (sc); 4515 } 4516 } 4517 return (sc); 4518} 4519 4520/* Powers of ten used in decimal <-> binary conversions. */ 4521 4522#define NTEN 12 4523#define MAXP 4096 4524 4525#if LONG_DOUBLE_TYPE_SIZE == 128 4526static unsigned EMUSHORT etens[NTEN + 1][NE] = 4527{ 4528 {0x6576, 0x4a92, 0x804a, 0x153f, 4529 0xc94c, 0x979a, 0x8a20, 0x5202, 0xc460, 0x7525,}, /* 10**4096 */ 4530 {0x6a32, 0xce52, 0x329a, 0x28ce, 4531 0xa74d, 0x5de4, 0xc53d, 0x3b5d, 0x9e8b, 0x5a92,}, /* 10**2048 */ 4532 {0x526c, 0x50ce, 0xf18b, 0x3d28, 4533 0x650d, 0x0c17, 0x8175, 0x7586, 0xc976, 0x4d48,}, 4534 {0x9c66, 0x58f8, 0xbc50, 0x5c54, 4535 0xcc65, 0x91c6, 0xa60e, 0xa0ae, 0xe319, 0x46a3,}, 4536 {0x851e, 0xeab7, 0x98fe, 0x901b, 4537 0xddbb, 0xde8d, 0x9df9, 0xebfb, 0xaa7e, 0x4351,}, 4538 {0x0235, 0x0137, 0x36b1, 0x336c, 4539 0xc66f, 0x8cdf, 0x80e9, 0x47c9, 0x93ba, 0x41a8,}, 4540 {0x50f8, 0x25fb, 0xc76b, 0x6b71, 4541 0x3cbf, 0xa6d5, 0xffcf, 0x1f49, 0xc278, 0x40d3,}, 4542 {0x0000, 0x0000, 0x0000, 0x0000, 4543 0xf020, 0xb59d, 0x2b70, 0xada8, 0x9dc5, 0x4069,}, 4544 {0x0000, 0x0000, 0x0000, 0x0000, 4545 0x0000, 0x0000, 0x0400, 0xc9bf, 0x8e1b, 0x4034,}, 4546 {0x0000, 0x0000, 0x0000, 0x0000, 4547 0x0000, 0x0000, 0x0000, 0x2000, 0xbebc, 0x4019,}, 4548 {0x0000, 0x0000, 0x0000, 0x0000, 4549 0x0000, 0x0000, 0x0000, 0x0000, 0x9c40, 0x400c,}, 4550 {0x0000, 0x0000, 0x0000, 0x0000, 4551 0x0000, 0x0000, 0x0000, 0x0000, 0xc800, 0x4005,}, 4552 {0x0000, 0x0000, 0x0000, 0x0000, 4553 0x0000, 0x0000, 0x0000, 0x0000, 0xa000, 0x4002,}, /* 10**1 */ 4554}; 4555 4556static unsigned EMUSHORT emtens[NTEN + 1][NE] = 4557{ 4558 {0x2030, 0xcffc, 0xa1c3, 0x8123, 4559 0x2de3, 0x9fde, 0xd2ce, 0x04c8, 0xa6dd, 0x0ad8,}, /* 10**-4096 */ 4560 {0x8264, 0xd2cb, 0xf2ea, 0x12d4, 4561 0x4925, 0x2de4, 0x3436, 0x534f, 0xceae, 0x256b,}, /* 10**-2048 */ 4562 {0xf53f, 0xf698, 0x6bd3, 0x0158, 4563 0x87a6, 0xc0bd, 0xda57, 0x82a5, 0xa2a6, 0x32b5,}, 4564 {0xe731, 0x04d4, 0xe3f2, 0xd332, 4565 0x7132, 0xd21c, 0xdb23, 0xee32, 0x9049, 0x395a,}, 4566 {0xa23e, 0x5308, 0xfefb, 0x1155, 4567 0xfa91, 0x1939, 0x637a, 0x4325, 0xc031, 0x3cac,}, 4568 {0xe26d, 0xdbde, 0xd05d, 0xb3f6, 4569 0xac7c, 0xe4a0, 0x64bc, 0x467c, 0xddd0, 0x3e55,}, 4570 {0x2a20, 0x6224, 0x47b3, 0x98d7, 4571 0x3f23, 0xe9a5, 0xa539, 0xea27, 0xa87f, 0x3f2a,}, 4572 {0x0b5b, 0x4af2, 0xa581, 0x18ed, 4573 0x67de, 0x94ba, 0x4539, 0x1ead, 0xcfb1, 0x3f94,}, 4574 {0xbf71, 0xa9b3, 0x7989, 0xbe68, 4575 0x4c2e, 0xe15b, 0xc44d, 0x94be, 0xe695, 0x3fc9,}, 4576 {0x3d4d, 0x7c3d, 0x36ba, 0x0d2b, 4577 0xfdc2, 0xcefc, 0x8461, 0x7711, 0xabcc, 0x3fe4,}, 4578 {0xc155, 0xa4a8, 0x404e, 0x6113, 4579 0xd3c3, 0x652b, 0xe219, 0x1758, 0xd1b7, 0x3ff1,}, 4580 {0xd70a, 0x70a3, 0x0a3d, 0xa3d7, 4581 0x3d70, 0xd70a, 0x70a3, 0x0a3d, 0xa3d7, 0x3ff8,}, 4582 {0xcccd, 0xcccc, 0xcccc, 0xcccc, 4583 0xcccc, 0xcccc, 0xcccc, 0xcccc, 0xcccc, 0x3ffb,}, /* 10**-1 */ 4584}; 4585#else 4586/* LONG_DOUBLE_TYPE_SIZE is other than 128 */ 4587static unsigned EMUSHORT etens[NTEN + 1][NE] = 4588{ 4589 {0xc94c, 0x979a, 0x8a20, 0x5202, 0xc460, 0x7525,}, /* 10**4096 */ 4590 {0xa74d, 0x5de4, 0xc53d, 0x3b5d, 0x9e8b, 0x5a92,}, /* 10**2048 */ 4591 {0x650d, 0x0c17, 0x8175, 0x7586, 0xc976, 0x4d48,}, 4592 {0xcc65, 0x91c6, 0xa60e, 0xa0ae, 0xe319, 0x46a3,}, 4593 {0xddbc, 0xde8d, 0x9df9, 0xebfb, 0xaa7e, 0x4351,}, 4594 {0xc66f, 0x8cdf, 0x80e9, 0x47c9, 0x93ba, 0x41a8,}, 4595 {0x3cbf, 0xa6d5, 0xffcf, 0x1f49, 0xc278, 0x40d3,}, 4596 {0xf020, 0xb59d, 0x2b70, 0xada8, 0x9dc5, 0x4069,}, 4597 {0x0000, 0x0000, 0x0400, 0xc9bf, 0x8e1b, 0x4034,}, 4598 {0x0000, 0x0000, 0x0000, 0x2000, 0xbebc, 0x4019,}, 4599 {0x0000, 0x0000, 0x0000, 0x0000, 0x9c40, 0x400c,}, 4600 {0x0000, 0x0000, 0x0000, 0x0000, 0xc800, 0x4005,}, 4601 {0x0000, 0x0000, 0x0000, 0x0000, 0xa000, 0x4002,}, /* 10**1 */ 4602}; 4603 4604static unsigned EMUSHORT emtens[NTEN + 1][NE] = 4605{ 4606 {0x2de4, 0x9fde, 0xd2ce, 0x04c8, 0xa6dd, 0x0ad8,}, /* 10**-4096 */ 4607 {0x4925, 0x2de4, 0x3436, 0x534f, 0xceae, 0x256b,}, /* 10**-2048 */ 4608 {0x87a6, 0xc0bd, 0xda57, 0x82a5, 0xa2a6, 0x32b5,}, 4609 {0x7133, 0xd21c, 0xdb23, 0xee32, 0x9049, 0x395a,}, 4610 {0xfa91, 0x1939, 0x637a, 0x4325, 0xc031, 0x3cac,}, 4611 {0xac7d, 0xe4a0, 0x64bc, 0x467c, 0xddd0, 0x3e55,}, 4612 {0x3f24, 0xe9a5, 0xa539, 0xea27, 0xa87f, 0x3f2a,}, 4613 {0x67de, 0x94ba, 0x4539, 0x1ead, 0xcfb1, 0x3f94,}, 4614 {0x4c2f, 0xe15b, 0xc44d, 0x94be, 0xe695, 0x3fc9,}, 4615 {0xfdc2, 0xcefc, 0x8461, 0x7711, 0xabcc, 0x3fe4,}, 4616 {0xd3c3, 0x652b, 0xe219, 0x1758, 0xd1b7, 0x3ff1,}, 4617 {0x3d71, 0xd70a, 0x70a3, 0x0a3d, 0xa3d7, 0x3ff8,}, 4618 {0xcccd, 0xcccc, 0xcccc, 0xcccc, 0xcccc, 0x3ffb,}, /* 10**-1 */ 4619}; 4620#endif 4621 4622#if 0 4623/* Convert float value X to ASCII string STRING with NDIG digits after 4624 the decimal point. */ 4625 4626static void 4627e24toasc (x, string, ndigs) 4628 unsigned EMUSHORT x[]; 4629 char *string; 4630 int ndigs; 4631{ 4632 unsigned EMUSHORT w[NI]; 4633 4634 e24toe (x, w); 4635 etoasc (w, string, ndigs); 4636} 4637 4638/* Convert double value X to ASCII string STRING with NDIG digits after 4639 the decimal point. */ 4640 4641static void 4642e53toasc (x, string, ndigs) 4643 unsigned EMUSHORT x[]; 4644 char *string; 4645 int ndigs; 4646{ 4647 unsigned EMUSHORT w[NI]; 4648 4649 e53toe (x, w); 4650 etoasc (w, string, ndigs); 4651} 4652 4653/* Convert double extended value X to ASCII string STRING with NDIG digits 4654 after the decimal point. */ 4655 4656static void 4657e64toasc (x, string, ndigs) 4658 unsigned EMUSHORT x[]; 4659 char *string; 4660 int ndigs; 4661{ 4662 unsigned EMUSHORT w[NI]; 4663 4664 e64toe (x, w); 4665 etoasc (w, string, ndigs); 4666} 4667 4668/* Convert 128-bit long double value X to ASCII string STRING with NDIG digits 4669 after the decimal point. */ 4670 4671static void 4672e113toasc (x, string, ndigs) 4673 unsigned EMUSHORT x[]; 4674 char *string; 4675 int ndigs; 4676{ 4677 unsigned EMUSHORT w[NI]; 4678 4679 e113toe (x, w); 4680 etoasc (w, string, ndigs); 4681} 4682#endif /* 0 */ 4683 4684/* Convert e-type X to ASCII string STRING with NDIGS digits after 4685 the decimal point. */ 4686 4687static char wstring[80]; /* working storage for ASCII output */ 4688 4689static void 4690etoasc (x, string, ndigs) 4691 unsigned EMUSHORT x[]; 4692 char *string; 4693 int ndigs; 4694{ 4695 EMUSHORT digit; 4696 unsigned EMUSHORT y[NI], t[NI], u[NI], w[NI]; 4697 unsigned EMUSHORT *p, *r, *ten; 4698 unsigned EMUSHORT sign; 4699 int i, j, k, expon, rndsav; 4700 char *s, *ss; 4701 unsigned EMUSHORT m; 4702 4703 4704 rndsav = rndprc; 4705 ss = string; 4706 s = wstring; 4707 *ss = '\0'; 4708 *s = '\0'; 4709#ifdef NANS 4710 if (eisnan (x)) 4711 { 4712 sprintf (wstring, " NaN "); 4713 goto bxit; 4714 } 4715#endif 4716 rndprc = NBITS; /* set to full precision */ 4717 emov (x, y); /* retain external format */ 4718 if (y[NE - 1] & 0x8000) 4719 { 4720 sign = 0xffff; 4721 y[NE - 1] &= 0x7fff; 4722 } 4723 else 4724 { 4725 sign = 0; 4726 } 4727 expon = 0; 4728 ten = &etens[NTEN][0]; 4729 emov (eone, t); 4730 /* Test for zero exponent */ 4731 if (y[NE - 1] == 0) 4732 { 4733 for (k = 0; k < NE - 1; k++) 4734 { 4735 if (y[k] != 0) 4736 goto tnzro; /* denormalized number */ 4737 } 4738 goto isone; /* valid all zeros */ 4739 } 4740 tnzro: 4741 4742 /* Test for infinity. */ 4743 if (y[NE - 1] == 0x7fff) 4744 { 4745 if (sign) 4746 sprintf (wstring, " -Infinity "); 4747 else 4748 sprintf (wstring, " Infinity "); 4749 goto bxit; 4750 } 4751 4752 /* Test for exponent nonzero but significand denormalized. 4753 * This is an error condition. 4754 */ 4755 if ((y[NE - 1] != 0) && ((y[NE - 2] & 0x8000) == 0)) 4756 { 4757 mtherr ("etoasc", DOMAIN); 4758 sprintf (wstring, "NaN"); 4759 goto bxit; 4760 } 4761 4762 /* Compare to 1.0 */ 4763 i = ecmp (eone, y); 4764 if (i == 0) 4765 goto isone; 4766 4767 if (i == -2) 4768 abort (); 4769 4770 if (i < 0) 4771 { /* Number is greater than 1 */ 4772 /* Convert significand to an integer and strip trailing decimal zeros. */ 4773 emov (y, u); 4774 u[NE - 1] = EXONE + NBITS - 1; 4775 4776 p = &etens[NTEN - 4][0]; 4777 m = 16; 4778 do 4779 { 4780 ediv (p, u, t); 4781 efloor (t, w); 4782 for (j = 0; j < NE - 1; j++) 4783 { 4784 if (t[j] != w[j]) 4785 goto noint; 4786 } 4787 emov (t, u); 4788 expon += (int) m; 4789 noint: 4790 p += NE; 4791 m >>= 1; 4792 } 4793 while (m != 0); 4794 4795 /* Rescale from integer significand */ 4796 u[NE - 1] += y[NE - 1] - (unsigned int) (EXONE + NBITS - 1); 4797 emov (u, y); 4798 /* Find power of 10 */ 4799 emov (eone, t); 4800 m = MAXP; 4801 p = &etens[0][0]; 4802 /* An unordered compare result shouldn't happen here. */ 4803 while (ecmp (ten, u) <= 0) 4804 { 4805 if (ecmp (p, u) <= 0) 4806 { 4807 ediv (p, u, u); 4808 emul (p, t, t); 4809 expon += (int) m; 4810 } 4811 m >>= 1; 4812 if (m == 0) 4813 break; 4814 p += NE; 4815 } 4816 } 4817 else 4818 { /* Number is less than 1.0 */ 4819 /* Pad significand with trailing decimal zeros. */ 4820 if (y[NE - 1] == 0) 4821 { 4822 while ((y[NE - 2] & 0x8000) == 0) 4823 { 4824 emul (ten, y, y); 4825 expon -= 1; 4826 } 4827 } 4828 else 4829 { 4830 emovi (y, w); 4831 for (i = 0; i < NDEC + 1; i++) 4832 { 4833 if ((w[NI - 1] & 0x7) != 0) 4834 break; 4835 /* multiply by 10 */ 4836 emovz (w, u); 4837 eshdn1 (u); 4838 eshdn1 (u); 4839 eaddm (w, u); 4840 u[1] += 3; 4841 while (u[2] != 0) 4842 { 4843 eshdn1 (u); 4844 u[1] += 1; 4845 } 4846 if (u[NI - 1] != 0) 4847 break; 4848 if (eone[NE - 1] <= u[1]) 4849 break; 4850 emovz (u, w); 4851 expon -= 1; 4852 } 4853 emovo (w, y); 4854 } 4855 k = -MAXP; 4856 p = &emtens[0][0]; 4857 r = &etens[0][0]; 4858 emov (y, w); 4859 emov (eone, t); 4860 while (ecmp (eone, w) > 0) 4861 { 4862 if (ecmp (p, w) >= 0) 4863 { 4864 emul (r, w, w); 4865 emul (r, t, t); 4866 expon += k; 4867 } 4868 k /= 2; 4869 if (k == 0) 4870 break; 4871 p += NE; 4872 r += NE; 4873 } 4874 ediv (t, eone, t); 4875 } 4876 isone: 4877 /* Find the first (leading) digit. */ 4878 emovi (t, w); 4879 emovz (w, t); 4880 emovi (y, w); 4881 emovz (w, y); 4882 eiremain (t, y); 4883 digit = equot[NI - 1]; 4884 while ((digit == 0) && (ecmp (y, ezero) != 0)) 4885 { 4886 eshup1 (y); 4887 emovz (y, u); 4888 eshup1 (u); 4889 eshup1 (u); 4890 eaddm (u, y); 4891 eiremain (t, y); 4892 digit = equot[NI - 1]; 4893 expon -= 1; 4894 } 4895 s = wstring; 4896 if (sign) 4897 *s++ = '-'; 4898 else 4899 *s++ = ' '; 4900 /* Examine number of digits requested by caller. */ 4901 if (ndigs < 0) 4902 ndigs = 0; 4903 if (ndigs > NDEC) 4904 ndigs = NDEC; 4905 if (digit == 10) 4906 { 4907 *s++ = '1'; 4908 *s++ = '.'; 4909 if (ndigs > 0) 4910 { 4911 *s++ = '0'; 4912 ndigs -= 1; 4913 } 4914 expon += 1; 4915 } 4916 else 4917 { 4918 *s++ = (char)digit + '0'; 4919 *s++ = '.'; 4920 } 4921 /* Generate digits after the decimal point. */ 4922 for (k = 0; k <= ndigs; k++) 4923 { 4924 /* multiply current number by 10, without normalizing */ 4925 eshup1 (y); 4926 emovz (y, u); 4927 eshup1 (u); 4928 eshup1 (u); 4929 eaddm (u, y); 4930 eiremain (t, y); 4931 *s++ = (char) equot[NI - 1] + '0'; 4932 } 4933 digit = equot[NI - 1]; 4934 --s; 4935 ss = s; 4936 /* round off the ASCII string */ 4937 if (digit > 4) 4938 { 4939 /* Test for critical rounding case in ASCII output. */ 4940 if (digit == 5) 4941 { 4942 emovo (y, t); 4943 if (ecmp (t, ezero) != 0) 4944 goto roun; /* round to nearest */ 4945#ifndef C4X 4946 if ((*(s - 1) & 1) == 0) 4947 goto doexp; /* round to even */ 4948#endif 4949 } 4950 /* Round up and propagate carry-outs */ 4951 roun: 4952 --s; 4953 k = *s & 0x7f; 4954 /* Carry out to most significant digit? */ 4955 if (k == '.') 4956 { 4957 --s; 4958 k = *s; 4959 k += 1; 4960 *s = (char) k; 4961 /* Most significant digit carries to 10? */ 4962 if (k > '9') 4963 { 4964 expon += 1; 4965 *s = '1'; 4966 } 4967 goto doexp; 4968 } 4969 /* Round up and carry out from less significant digits */ 4970 k += 1; 4971 *s = (char) k; 4972 if (k > '9') 4973 { 4974 *s = '0'; 4975 goto roun; 4976 } 4977 } 4978 doexp: 4979 /* 4980 if (expon >= 0) 4981 sprintf (ss, "e+%d", expon); 4982 else 4983 sprintf (ss, "e%d", expon); 4984 */ 4985 sprintf (ss, "e%d", expon); 4986 bxit: 4987 rndprc = rndsav; 4988 /* copy out the working string */ 4989 s = string; 4990 ss = wstring; 4991 while (*ss == ' ') /* strip possible leading space */ 4992 ++ss; 4993 while ((*s++ = *ss++) != '\0') 4994 ; 4995} 4996 4997 4998/* Convert ASCII string to floating point. 4999 5000 Numeric input is a free format decimal number of any length, with 5001 or without decimal point. Entering E after the number followed by an 5002 integer number causes the second number to be interpreted as a power of 5003 10 to be multiplied by the first number (i.e., "scientific" notation). */ 5004 5005/* Convert ASCII string S to single precision float value Y. */ 5006 5007static void 5008asctoe24 (s, y) 5009 const char *s; 5010 unsigned EMUSHORT *y; 5011{ 5012 asctoeg (s, y, 24); 5013} 5014 5015 5016/* Convert ASCII string S to double precision value Y. */ 5017 5018static void 5019asctoe53 (s, y) 5020 const char *s; 5021 unsigned EMUSHORT *y; 5022{ 5023#if defined(DEC) || defined(IBM) 5024 asctoeg (s, y, 56); 5025#else 5026#if defined(C4X) 5027 asctoeg (s, y, 32); 5028#else 5029 asctoeg (s, y, 53); 5030#endif 5031#endif 5032} 5033 5034 5035/* Convert ASCII string S to double extended value Y. */ 5036 5037static void 5038asctoe64 (s, y) 5039 const char *s; 5040 unsigned EMUSHORT *y; 5041{ 5042 asctoeg (s, y, 64); 5043} 5044 5045/* Convert ASCII string S to 128-bit long double Y. */ 5046 5047static void 5048asctoe113 (s, y) 5049 const char *s; 5050 unsigned EMUSHORT *y; 5051{ 5052 asctoeg (s, y, 113); 5053} 5054 5055/* Convert ASCII string S to e type Y. */ 5056 5057static void 5058asctoe (s, y) 5059 const char *s; 5060 unsigned EMUSHORT *y; 5061{ 5062 asctoeg (s, y, NBITS); 5063} 5064 5065/* Convert ASCII string SS to e type Y, with a specified rounding precision 5066 of OPREC bits. BASE is 16 for C9X hexadecimal floating constants. */ 5067 5068static void 5069asctoeg (ss, y, oprec) 5070 const char *ss; 5071 unsigned EMUSHORT *y; 5072 int oprec; 5073{ 5074 unsigned EMUSHORT yy[NI], xt[NI], tt[NI]; 5075 int esign, decflg, sgnflg, nexp, exp, prec, lost; 5076 int k, trail, c, rndsav; 5077 EMULONG lexp; 5078 unsigned EMUSHORT nsign, *p; 5079 char *sp, *s, *lstr; 5080 int base = 10; 5081 5082 /* Copy the input string. */ 5083 lstr = (char *) alloca (strlen (ss) + 1); 5084 5085 while (*ss == ' ') /* skip leading spaces */ 5086 ++ss; 5087 5088 sp = lstr; 5089 while ((*sp++ = *ss++) != '\0') 5090 ; 5091 s = lstr; 5092 5093 if (s[0] == '0' && (s[1] == 'x' || s[1] == 'X')) 5094 { 5095 base = 16; 5096 s += 2; 5097 } 5098 5099 rndsav = rndprc; 5100 rndprc = NBITS; /* Set to full precision */ 5101 lost = 0; 5102 nsign = 0; 5103 decflg = 0; 5104 sgnflg = 0; 5105 nexp = 0; 5106 exp = 0; 5107 prec = 0; 5108 ecleaz (yy); 5109 trail = 0; 5110 5111 nxtcom: 5112 if (*s >= '0' && *s <= '9') 5113 k = *s - '0'; 5114 else if (*s >= 'a') 5115 k = 10 + *s - 'a'; 5116 else 5117 k = 10 + *s - 'A'; 5118 if ((k >= 0) && (k < base)) 5119 { 5120 /* Ignore leading zeros */ 5121 if ((prec == 0) && (decflg == 0) && (k == 0)) 5122 goto donchr; 5123 /* Identify and strip trailing zeros after the decimal point. */ 5124 if ((trail == 0) && (decflg != 0)) 5125 { 5126 sp = s; 5127 while ((*sp >= '0' && *sp <= '9') 5128 || (base == 16 && ((*sp >= 'a' && *sp <= 'f') 5129 || (*sp >= 'A' && *sp <= 'F')))) 5130 ++sp; 5131 /* Check for syntax error */ 5132 c = *sp & 0x7f; 5133 if ((base != 10 || ((c != 'e') && (c != 'E'))) 5134 && (base != 16 || ((c != 'p') && (c != 'P'))) 5135 && (c != '\0') 5136 && (c != '\n') && (c != '\r') && (c != ' ') 5137 && (c != ',')) 5138 goto error; 5139 --sp; 5140 while (*sp == '0') 5141 *sp-- = 'z'; 5142 trail = 1; 5143 if (*s == 'z') 5144 goto donchr; 5145 } 5146 5147 /* If enough digits were given to more than fill up the yy register, 5148 continuing until overflow into the high guard word yy[2] 5149 guarantees that there will be a roundoff bit at the top 5150 of the low guard word after normalization. */ 5151 5152 if (yy[2] == 0) 5153 { 5154 if (base == 16) 5155 { 5156 if (decflg) 5157 nexp += 4; /* count digits after decimal point */ 5158 5159 eshup1 (yy); /* multiply current number by 16 */ 5160 eshup1 (yy); 5161 eshup1 (yy); 5162 eshup1 (yy); 5163 } 5164 else 5165 { 5166 if (decflg) 5167 nexp += 1; /* count digits after decimal point */ 5168 5169 eshup1 (yy); /* multiply current number by 10 */ 5170 emovz (yy, xt); 5171 eshup1 (xt); 5172 eshup1 (xt); 5173 eaddm (xt, yy); 5174 } 5175 /* Insert the current digit. */ 5176 ecleaz (xt); 5177 xt[NI - 2] = (unsigned EMUSHORT) k; 5178 eaddm (xt, yy); 5179 } 5180 else 5181 { 5182 /* Mark any lost non-zero digit. */ 5183 lost |= k; 5184 /* Count lost digits before the decimal point. */ 5185 if (decflg == 0) 5186 { 5187 if (base == 10) 5188 nexp -= 1; 5189 else 5190 nexp -= 4; 5191 } 5192 } 5193 prec += 1; 5194 goto donchr; 5195 } 5196 5197 switch (*s) 5198 { 5199 case 'z': 5200 break; 5201 case 'E': 5202 case 'e': 5203 case 'P': 5204 case 'p': 5205 goto expnt; 5206 case '.': /* decimal point */ 5207 if (decflg) 5208 goto error; 5209 ++decflg; 5210 break; 5211 case '-': 5212 nsign = 0xffff; 5213 if (sgnflg) 5214 goto error; 5215 ++sgnflg; 5216 break; 5217 case '+': 5218 if (sgnflg) 5219 goto error; 5220 ++sgnflg; 5221 break; 5222 case ',': 5223 case ' ': 5224 case '\0': 5225 case '\n': 5226 case '\r': 5227 goto daldone; 5228 case 'i': 5229 case 'I': 5230 goto infinite; 5231 default: 5232 error: 5233#ifdef NANS 5234 einan (yy); 5235#else 5236 mtherr ("asctoe", DOMAIN); 5237 eclear (yy); 5238#endif 5239 goto aexit; 5240 } 5241 donchr: 5242 ++s; 5243 goto nxtcom; 5244 5245 /* Exponent interpretation */ 5246 expnt: 5247 /* 0.0eXXX is zero, regardless of XXX. Check for the 0.0. */ 5248 for (k = 0; k < NI; k++) 5249 { 5250 if (yy[k] != 0) 5251 goto read_expnt; 5252 } 5253 goto aexit; 5254 5255read_expnt: 5256 esign = 1; 5257 exp = 0; 5258 ++s; 5259 /* check for + or - */ 5260 if (*s == '-') 5261 { 5262 esign = -1; 5263 ++s; 5264 } 5265 if (*s == '+') 5266 ++s; 5267 while ((*s >= '0') && (*s <= '9')) 5268 { 5269 exp *= 10; 5270 exp += *s++ - '0'; 5271 if (exp > 999999) 5272 break; 5273 } 5274 if (esign < 0) 5275 exp = -exp; 5276 if ((exp > MAXDECEXP) && (base == 10)) 5277 { 5278 infinite: 5279 ecleaz (yy); 5280 yy[E] = 0x7fff; /* infinity */ 5281 goto aexit; 5282 } 5283 if ((exp < MINDECEXP) && (base == 10)) 5284 { 5285 zero: 5286 ecleaz (yy); 5287 goto aexit; 5288 } 5289 5290 daldone: 5291 if (base == 16) 5292 { 5293 /* Base 16 hexadecimal floating constant. */ 5294 if ((k = enormlz (yy)) > NBITS) 5295 { 5296 ecleaz (yy); 5297 goto aexit; 5298 } 5299 /* Adjust the exponent. NEXP is the number of hex digits, 5300 EXP is a power of 2. */ 5301 lexp = (EXONE - 1 + NBITS) - k + yy[E] + exp - nexp; 5302 if (lexp > 0x7fff) 5303 goto infinite; 5304 if (lexp < 0) 5305 goto zero; 5306 yy[E] = lexp; 5307 goto expdon; 5308 } 5309 5310 nexp = exp - nexp; 5311 /* Pad trailing zeros to minimize power of 10, per IEEE spec. */ 5312 while ((nexp > 0) && (yy[2] == 0)) 5313 { 5314 emovz (yy, xt); 5315 eshup1 (xt); 5316 eshup1 (xt); 5317 eaddm (yy, xt); 5318 eshup1 (xt); 5319 if (xt[2] != 0) 5320 break; 5321 nexp -= 1; 5322 emovz (xt, yy); 5323 } 5324 if ((k = enormlz (yy)) > NBITS) 5325 { 5326 ecleaz (yy); 5327 goto aexit; 5328 } 5329 lexp = (EXONE - 1 + NBITS) - k; 5330 emdnorm (yy, lost, 0, lexp, 64); 5331 lost = 0; 5332 5333 /* Convert to external format: 5334 5335 Multiply by 10**nexp. If precision is 64 bits, 5336 the maximum relative error incurred in forming 10**n 5337 for 0 <= n <= 324 is 8.2e-20, at 10**180. 5338 For 0 <= n <= 999, the peak relative error is 1.4e-19 at 10**947. 5339 For 0 >= n >= -999, it is -1.55e-19 at 10**-435. */ 5340 5341 lexp = yy[E]; 5342 if (nexp == 0) 5343 { 5344 k = 0; 5345 goto expdon; 5346 } 5347 esign = 1; 5348 if (nexp < 0) 5349 { 5350 nexp = -nexp; 5351 esign = -1; 5352 if (nexp > 4096) 5353 { 5354 /* Punt. Can't handle this without 2 divides. */ 5355 emovi (etens[0], tt); 5356 lexp -= tt[E]; 5357 k = edivm (tt, yy); 5358 lexp += EXONE; 5359 nexp -= 4096; 5360 } 5361 } 5362 p = &etens[NTEN][0]; 5363 emov (eone, xt); 5364 exp = 1; 5365 do 5366 { 5367 if (exp & nexp) 5368 emul (p, xt, xt); 5369 p -= NE; 5370 exp = exp + exp; 5371 } 5372 while (exp <= MAXP); 5373 5374 emovi (xt, tt); 5375 if (esign < 0) 5376 { 5377 lexp -= tt[E]; 5378 k = edivm (tt, yy); 5379 lexp += EXONE; 5380 } 5381 else 5382 { 5383 lexp += tt[E]; 5384 k = emulm (tt, yy); 5385 lexp -= EXONE - 1; 5386 } 5387 lost = k; 5388 5389 expdon: 5390 5391 /* Round and convert directly to the destination type */ 5392 if (oprec == 53) 5393 lexp -= EXONE - 0x3ff; 5394#ifdef C4X 5395 else if (oprec == 24 || oprec == 32) 5396 lexp -= (EXONE - 0x7f); 5397#else 5398#ifdef IBM 5399 else if (oprec == 24 || oprec == 56) 5400 lexp -= EXONE - (0x41 << 2); 5401#else 5402 else if (oprec == 24) 5403 lexp -= EXONE - 0177; 5404#endif /* IBM */ 5405#endif /* C4X */ 5406#ifdef DEC 5407 else if (oprec == 56) 5408 lexp -= EXONE - 0201; 5409#endif 5410 rndprc = oprec; 5411 emdnorm (yy, lost, 0, lexp, 64); 5412 5413 aexit: 5414 5415 rndprc = rndsav; 5416 yy[0] = nsign; 5417 switch (oprec) 5418 { 5419#ifdef DEC 5420 case 56: 5421 todec (yy, y); /* see etodec.c */ 5422 break; 5423#endif 5424#ifdef IBM 5425 case 56: 5426 toibm (yy, y, DFmode); 5427 break; 5428#endif 5429#ifdef C4X 5430 case 32: 5431 toc4x (yy, y, HFmode); 5432 break; 5433#endif 5434 5435 case 53: 5436 toe53 (yy, y); 5437 break; 5438 case 24: 5439 toe24 (yy, y); 5440 break; 5441 case 64: 5442 toe64 (yy, y); 5443 break; 5444 case 113: 5445 toe113 (yy, y); 5446 break; 5447 case NBITS: 5448 emovo (yy, y); 5449 break; 5450 } 5451} 5452 5453 5454 5455/* Return Y = largest integer not greater than X (truncated toward minus 5456 infinity). */ 5457 5458static unsigned EMUSHORT bmask[] = 5459{ 5460 0xffff, 5461 0xfffe, 5462 0xfffc, 5463 0xfff8, 5464 0xfff0, 5465 0xffe0, 5466 0xffc0, 5467 0xff80, 5468 0xff00, 5469 0xfe00, 5470 0xfc00, 5471 0xf800, 5472 0xf000, 5473 0xe000, 5474 0xc000, 5475 0x8000, 5476 0x0000, 5477}; 5478 5479static void 5480efloor (x, y) 5481 unsigned EMUSHORT x[], y[]; 5482{ 5483 register unsigned EMUSHORT *p; 5484 int e, expon, i; 5485 unsigned EMUSHORT f[NE]; 5486 5487 emov (x, f); /* leave in external format */ 5488 expon = (int) f[NE - 1]; 5489 e = (expon & 0x7fff) - (EXONE - 1); 5490 if (e <= 0) 5491 { 5492 eclear (y); 5493 goto isitneg; 5494 } 5495 /* number of bits to clear out */ 5496 e = NBITS - e; 5497 emov (f, y); 5498 if (e <= 0) 5499 return; 5500 5501 p = &y[0]; 5502 while (e >= 16) 5503 { 5504 *p++ = 0; 5505 e -= 16; 5506 } 5507 /* clear the remaining bits */ 5508 *p &= bmask[e]; 5509 /* truncate negatives toward minus infinity */ 5510 isitneg: 5511 5512 if ((unsigned EMUSHORT) expon & (unsigned EMUSHORT) 0x8000) 5513 { 5514 for (i = 0; i < NE - 1; i++) 5515 { 5516 if (f[i] != y[i]) 5517 { 5518 esub (eone, y, y); 5519 break; 5520 } 5521 } 5522 } 5523} 5524 5525 5526#if 0 5527/* Return S and EXP such that S * 2^EXP = X and .5 <= S < 1. 5528 For example, 1.1 = 0.55 * 2^1. */ 5529 5530static void 5531efrexp (x, exp, s) 5532 unsigned EMUSHORT x[]; 5533 int *exp; 5534 unsigned EMUSHORT s[]; 5535{ 5536 unsigned EMUSHORT xi[NI]; 5537 EMULONG li; 5538 5539 emovi (x, xi); 5540 /* Handle denormalized numbers properly using long integer exponent. */ 5541 li = (EMULONG) ((EMUSHORT) xi[1]); 5542 5543 if (li == 0) 5544 { 5545 li -= enormlz (xi); 5546 } 5547 xi[1] = 0x3ffe; 5548 emovo (xi, s); 5549 *exp = (int) (li - 0x3ffe); 5550} 5551#endif 5552 5553/* Return e type Y = X * 2^PWR2. */ 5554 5555static void 5556eldexp (x, pwr2, y) 5557 unsigned EMUSHORT x[]; 5558 int pwr2; 5559 unsigned EMUSHORT y[]; 5560{ 5561 unsigned EMUSHORT xi[NI]; 5562 EMULONG li; 5563 int i; 5564 5565 emovi (x, xi); 5566 li = xi[1]; 5567 li += pwr2; 5568 i = 0; 5569 emdnorm (xi, i, i, li, 64); 5570 emovo (xi, y); 5571} 5572 5573 5574#if 0 5575/* C = remainder after dividing B by A, all e type values. 5576 Least significant integer quotient bits left in EQUOT. */ 5577 5578static void 5579eremain (a, b, c) 5580 unsigned EMUSHORT a[], b[], c[]; 5581{ 5582 unsigned EMUSHORT den[NI], num[NI]; 5583 5584#ifdef NANS 5585 if (eisinf (b) 5586 || (ecmp (a, ezero) == 0) 5587 || eisnan (a) 5588 || eisnan (b)) 5589 { 5590 enan (c, 0); 5591 return; 5592 } 5593#endif 5594 if (ecmp (a, ezero) == 0) 5595 { 5596 mtherr ("eremain", SING); 5597 eclear (c); 5598 return; 5599 } 5600 emovi (a, den); 5601 emovi (b, num); 5602 eiremain (den, num); 5603 /* Sign of remainder = sign of quotient */ 5604 if (a[0] == b[0]) 5605 num[0] = 0; 5606 else 5607 num[0] = 0xffff; 5608 emovo (num, c); 5609} 5610#endif 5611 5612/* Return quotient of exploded e-types NUM / DEN in EQUOT, 5613 remainder in NUM. */ 5614 5615static void 5616eiremain (den, num) 5617 unsigned EMUSHORT den[], num[]; 5618{ 5619 EMULONG ld, ln; 5620 unsigned EMUSHORT j; 5621 5622 ld = den[E]; 5623 ld -= enormlz (den); 5624 ln = num[E]; 5625 ln -= enormlz (num); 5626 ecleaz (equot); 5627 while (ln >= ld) 5628 { 5629 if (ecmpm (den, num) <= 0) 5630 { 5631 esubm (den, num); 5632 j = 1; 5633 } 5634 else 5635 j = 0; 5636 eshup1 (equot); 5637 equot[NI - 1] |= j; 5638 eshup1 (num); 5639 ln -= 1; 5640 } 5641 emdnorm (num, 0, 0, ln, 0); 5642} 5643 5644/* Report an error condition CODE encountered in function NAME. 5645 5646 Mnemonic Value Significance 5647 5648 DOMAIN 1 argument domain error 5649 SING 2 function singularity 5650 OVERFLOW 3 overflow range error 5651 UNDERFLOW 4 underflow range error 5652 TLOSS 5 total loss of precision 5653 PLOSS 6 partial loss of precision 5654 INVALID 7 NaN - producing operation 5655 EDOM 33 Unix domain error code 5656 ERANGE 34 Unix range error code 5657 5658 The order of appearance of the following messages is bound to the 5659 error codes defined above. */ 5660 5661int merror = 0; 5662extern int merror; 5663 5664static void 5665mtherr (name, code) 5666 const char *name; 5667 int code; 5668{ 5669 /* The string passed by the calling program is supposed to be the 5670 name of the function in which the error occurred. 5671 The code argument selects which error message string will be printed. */ 5672 5673 if (strcmp (name, "esub") == 0) 5674 name = "subtraction"; 5675 else if (strcmp (name, "ediv") == 0) 5676 name = "division"; 5677 else if (strcmp (name, "emul") == 0) 5678 name = "multiplication"; 5679 else if (strcmp (name, "enormlz") == 0) 5680 name = "normalization"; 5681 else if (strcmp (name, "etoasc") == 0) 5682 name = "conversion to text"; 5683 else if (strcmp (name, "asctoe") == 0) 5684 name = "parsing"; 5685 else if (strcmp (name, "eremain") == 0) 5686 name = "modulus"; 5687 else if (strcmp (name, "esqrt") == 0) 5688 name = "square root"; 5689 if (extra_warnings) 5690 { 5691 switch (code) 5692 { 5693 case DOMAIN: warning ("%s: argument domain error" , name); break; 5694 case SING: warning ("%s: function singularity" , name); break; 5695 case OVERFLOW: warning ("%s: overflow range error" , name); break; 5696 case UNDERFLOW: warning ("%s: underflow range error" , name); break; 5697 case TLOSS: warning ("%s: total loss of precision" , name); break; 5698 case PLOSS: warning ("%s: partial loss of precision", name); break; 5699 case INVALID: warning ("%s: NaN - producing operation", name); break; 5700 default: abort (); 5701 } 5702 } 5703 5704 /* Set global error message word */ 5705 merror = code + 1; 5706} 5707 5708#ifdef DEC 5709/* Convert DEC double precision D to e type E. */ 5710 5711static void 5712dectoe (d, e) 5713 unsigned EMUSHORT *d; 5714 unsigned EMUSHORT *e; 5715{ 5716 unsigned EMUSHORT y[NI]; 5717 register unsigned EMUSHORT r, *p; 5718 5719 ecleaz (y); /* start with a zero */ 5720 p = y; /* point to our number */ 5721 r = *d; /* get DEC exponent word */ 5722 if (*d & (unsigned int) 0x8000) 5723 *p = 0xffff; /* fill in our sign */ 5724 ++p; /* bump pointer to our exponent word */ 5725 r &= 0x7fff; /* strip the sign bit */ 5726 if (r == 0) /* answer = 0 if high order DEC word = 0 */ 5727 goto done; 5728 5729 5730 r >>= 7; /* shift exponent word down 7 bits */ 5731 r += EXONE - 0201; /* subtract DEC exponent offset */ 5732 /* add our e type exponent offset */ 5733 *p++ = r; /* to form our exponent */ 5734 5735 r = *d++; /* now do the high order mantissa */ 5736 r &= 0177; /* strip off the DEC exponent and sign bits */ 5737 r |= 0200; /* the DEC understood high order mantissa bit */ 5738 *p++ = r; /* put result in our high guard word */ 5739 5740 *p++ = *d++; /* fill in the rest of our mantissa */ 5741 *p++ = *d++; 5742 *p = *d; 5743 5744 eshdn8 (y); /* shift our mantissa down 8 bits */ 5745 done: 5746 emovo (y, e); 5747} 5748 5749/* Convert e type X to DEC double precision D. */ 5750 5751static void 5752etodec (x, d) 5753 unsigned EMUSHORT *x, *d; 5754{ 5755 unsigned EMUSHORT xi[NI]; 5756 EMULONG exp; 5757 int rndsav; 5758 5759 emovi (x, xi); 5760 /* Adjust exponent for offsets. */ 5761 exp = (EMULONG) xi[E] - (EXONE - 0201); 5762 /* Round off to nearest or even. */ 5763 rndsav = rndprc; 5764 rndprc = 56; 5765 emdnorm (xi, 0, 0, exp, 64); 5766 rndprc = rndsav; 5767 todec (xi, d); 5768} 5769 5770/* Convert exploded e-type X, that has already been rounded to 5771 56-bit precision, to DEC format double Y. */ 5772 5773static void 5774todec (x, y) 5775 unsigned EMUSHORT *x, *y; 5776{ 5777 unsigned EMUSHORT i; 5778 unsigned EMUSHORT *p; 5779 5780 p = x; 5781 *y = 0; 5782 if (*p++) 5783 *y = 0100000; 5784 i = *p++; 5785 if (i == 0) 5786 { 5787 *y++ = 0; 5788 *y++ = 0; 5789 *y++ = 0; 5790 *y++ = 0; 5791 return; 5792 } 5793 if (i > 0377) 5794 { 5795 *y++ |= 077777; 5796 *y++ = 0xffff; 5797 *y++ = 0xffff; 5798 *y++ = 0xffff; 5799#ifdef ERANGE 5800 errno = ERANGE; 5801#endif 5802 return; 5803 } 5804 i &= 0377; 5805 i <<= 7; 5806 eshup8 (x); 5807 x[M] &= 0177; 5808 i |= x[M]; 5809 *y++ |= i; 5810 *y++ = x[M + 1]; 5811 *y++ = x[M + 2]; 5812 *y++ = x[M + 3]; 5813} 5814#endif /* DEC */ 5815 5816#ifdef IBM 5817/* Convert IBM single/double precision to e type. */ 5818 5819static void 5820ibmtoe (d, e, mode) 5821 unsigned EMUSHORT *d; 5822 unsigned EMUSHORT *e; 5823 enum machine_mode mode; 5824{ 5825 unsigned EMUSHORT y[NI]; 5826 register unsigned EMUSHORT r, *p; 5827 int rndsav; 5828 5829 ecleaz (y); /* start with a zero */ 5830 p = y; /* point to our number */ 5831 r = *d; /* get IBM exponent word */ 5832 if (*d & (unsigned int) 0x8000) 5833 *p = 0xffff; /* fill in our sign */ 5834 ++p; /* bump pointer to our exponent word */ 5835 r &= 0x7f00; /* strip the sign bit */ 5836 r >>= 6; /* shift exponent word down 6 bits */ 5837 /* in fact shift by 8 right and 2 left */ 5838 r += EXONE - (0x41 << 2); /* subtract IBM exponent offset */ 5839 /* add our e type exponent offset */ 5840 *p++ = r; /* to form our exponent */ 5841 5842 *p++ = *d++ & 0xff; /* now do the high order mantissa */ 5843 /* strip off the IBM exponent and sign bits */ 5844 if (mode != SFmode) /* there are only 2 words in SFmode */ 5845 { 5846 *p++ = *d++; /* fill in the rest of our mantissa */ 5847 *p++ = *d++; 5848 } 5849 *p = *d; 5850 5851 if (y[M] == 0 && y[M+1] == 0 && y[M+2] == 0 && y[M+3] == 0) 5852 y[0] = y[E] = 0; 5853 else 5854 y[E] -= 5 + enormlz (y); /* now normalise the mantissa */ 5855 /* handle change in RADIX */ 5856 emovo (y, e); 5857} 5858 5859 5860 5861/* Convert e type to IBM single/double precision. */ 5862 5863static void 5864etoibm (x, d, mode) 5865 unsigned EMUSHORT *x, *d; 5866 enum machine_mode mode; 5867{ 5868 unsigned EMUSHORT xi[NI]; 5869 EMULONG exp; 5870 int rndsav; 5871 5872 emovi (x, xi); 5873 exp = (EMULONG) xi[E] - (EXONE - (0x41 << 2)); /* adjust exponent for offsets */ 5874 /* round off to nearest or even */ 5875 rndsav = rndprc; 5876 rndprc = 56; 5877 emdnorm (xi, 0, 0, exp, 64); 5878 rndprc = rndsav; 5879 toibm (xi, d, mode); 5880} 5881 5882static void 5883toibm (x, y, mode) 5884 unsigned EMUSHORT *x, *y; 5885 enum machine_mode mode; 5886{ 5887 unsigned EMUSHORT i; 5888 unsigned EMUSHORT *p; 5889 int r; 5890 5891 p = x; 5892 *y = 0; 5893 if (*p++) 5894 *y = 0x8000; 5895 i = *p++; 5896 if (i == 0) 5897 { 5898 *y++ = 0; 5899 *y++ = 0; 5900 if (mode != SFmode) 5901 { 5902 *y++ = 0; 5903 *y++ = 0; 5904 } 5905 return; 5906 } 5907 r = i & 0x3; 5908 i >>= 2; 5909 if (i > 0x7f) 5910 { 5911 *y++ |= 0x7fff; 5912 *y++ = 0xffff; 5913 if (mode != SFmode) 5914 { 5915 *y++ = 0xffff; 5916 *y++ = 0xffff; 5917 } 5918#ifdef ERANGE 5919 errno = ERANGE; 5920#endif 5921 return; 5922 } 5923 i &= 0x7f; 5924 *y |= (i << 8); 5925 eshift (x, r + 5); 5926 *y++ |= x[M]; 5927 *y++ = x[M + 1]; 5928 if (mode != SFmode) 5929 { 5930 *y++ = x[M + 2]; 5931 *y++ = x[M + 3]; 5932 } 5933} 5934#endif /* IBM */ 5935 5936 5937#ifdef C4X 5938/* Convert C4X single/double precision to e type. */ 5939 5940static void 5941c4xtoe (d, e, mode) 5942 unsigned EMUSHORT *d; 5943 unsigned EMUSHORT *e; 5944 enum machine_mode mode; 5945{ 5946 unsigned EMUSHORT y[NI]; 5947 int r; 5948 int isnegative; 5949 int size; 5950 int i; 5951 int carry; 5952 5953 /* Short-circuit the zero case. */ 5954 if ((d[0] == 0x8000) 5955 && (d[1] == 0x0000) 5956 && ((mode == QFmode) || ((d[2] == 0x0000) && (d[3] == 0x0000)))) 5957 { 5958 e[0] = 0; 5959 e[1] = 0; 5960 e[2] = 0; 5961 e[3] = 0; 5962 e[4] = 0; 5963 e[5] = 0; 5964 return; 5965 } 5966 5967 ecleaz (y); /* start with a zero */ 5968 r = d[0]; /* get sign/exponent part */ 5969 if (r & (unsigned int) 0x0080) 5970 { 5971 y[0] = 0xffff; /* fill in our sign */ 5972 isnegative = TRUE; 5973 } 5974 else 5975 { 5976 isnegative = FALSE; 5977 } 5978 5979 r >>= 8; /* Shift exponent word down 8 bits. */ 5980 if (r & 0x80) /* Make the exponent negative if it is. */ 5981 { 5982 r = r | (~0 & ~0xff); 5983 } 5984 5985 if (isnegative) 5986 { 5987 /* Now do the high order mantissa. We don't "or" on the high bit 5988 because it is 2 (not 1) and is handled a little differently 5989 below. */ 5990 y[M] = d[0] & 0x7f; 5991 5992 y[M+1] = d[1]; 5993 if (mode != QFmode) /* There are only 2 words in QFmode. */ 5994 { 5995 y[M+2] = d[2]; /* Fill in the rest of our mantissa. */ 5996 y[M+3] = d[3]; 5997 size = 4; 5998 } 5999 else 6000 { 6001 size = 2; 6002 } 6003 eshift(y, -8); 6004 6005 /* Now do the two's complement on the data. */ 6006 6007 carry = 1; /* Initially add 1 for the two's complement. */ 6008 for (i=size + M; i > M; i--) 6009 { 6010 if (carry && (y[i] == 0x0000)) 6011 { 6012 /* We overflowed into the next word, carry is the same. */ 6013 y[i] = carry ? 0x0000 : 0xffff; 6014 } 6015 else 6016 { 6017 /* No overflow, just invert and add carry. */ 6018 y[i] = ((~y[i]) + carry) & 0xffff; 6019 carry = 0; 6020 } 6021 } 6022 6023 if (carry) 6024 { 6025 eshift(y, -1); 6026 y[M+1] |= 0x8000; 6027 r++; 6028 } 6029 y[1] = r + EXONE; 6030 } 6031 else 6032 { 6033 /* Add our e type exponent offset to form our exponent. */ 6034 r += EXONE; 6035 y[1] = r; 6036 6037 /* Now do the high order mantissa strip off the exponent and sign 6038 bits and add the high 1 bit. */ 6039 y[M] = (d[0] & 0x7f) | 0x80; 6040 6041 y[M+1] = d[1]; 6042 if (mode != QFmode) /* There are only 2 words in QFmode. */ 6043 { 6044 y[M+2] = d[2]; /* Fill in the rest of our mantissa. */ 6045 y[M+3] = d[3]; 6046 } 6047 eshift(y, -8); 6048 } 6049 6050 emovo (y, e); 6051} 6052 6053 6054/* Convert e type to C4X single/double precision. */ 6055 6056static void 6057etoc4x (x, d, mode) 6058 unsigned EMUSHORT *x, *d; 6059 enum machine_mode mode; 6060{ 6061 unsigned EMUSHORT xi[NI]; 6062 EMULONG exp; 6063 int rndsav; 6064 6065 emovi (x, xi); 6066 6067 /* Adjust exponent for offsets. */ 6068 exp = (EMULONG) xi[E] - (EXONE - 0x7f); 6069 6070 /* Round off to nearest or even. */ 6071 rndsav = rndprc; 6072 rndprc = mode == QFmode ? 24 : 32; 6073 emdnorm (xi, 0, 0, exp, 64); 6074 rndprc = rndsav; 6075 toc4x (xi, d, mode); 6076} 6077 6078static void 6079toc4x (x, y, mode) 6080 unsigned EMUSHORT *x, *y; 6081 enum machine_mode mode; 6082{ 6083 int i; 6084 int v; 6085 int carry; 6086 6087 /* Short-circuit the zero case */ 6088 if ((x[0] == 0) /* Zero exponent and sign */ 6089 && (x[1] == 0) 6090 && (x[M] == 0) /* The rest is for zero mantissa */ 6091 && (x[M+1] == 0) 6092 /* Only check for double if necessary */ 6093 && ((mode == QFmode) || ((x[M+2] == 0) && (x[M+3] == 0)))) 6094 { 6095 /* We have a zero. Put it into the output and return. */ 6096 *y++ = 0x8000; 6097 *y++ = 0x0000; 6098 if (mode != QFmode) 6099 { 6100 *y++ = 0x0000; 6101 *y++ = 0x0000; 6102 } 6103 return; 6104 } 6105 6106 *y = 0; 6107 6108 /* Negative number require a two's complement conversion of the 6109 mantissa. */ 6110 if (x[0]) 6111 { 6112 *y = 0x0080; 6113 6114 i = ((int) x[1]) - 0x7f; 6115 6116 /* Now add 1 to the inverted data to do the two's complement. */ 6117 if (mode != QFmode) 6118 v = 4 + M; 6119 else 6120 v = 2 + M; 6121 carry = 1; 6122 while (v > M) 6123 { 6124 if (x[v] == 0x0000) 6125 { 6126 x[v] = carry ? 0x0000 : 0xffff; 6127 } 6128 else 6129 { 6130 x[v] = ((~x[v]) + carry) & 0xffff; 6131 carry = 0; 6132 } 6133 v--; 6134 } 6135 6136 /* The following is a special case. The C4X negative float requires 6137 a zero in the high bit (because the format is (2 - x) x 2^m), so 6138 if a one is in that bit, we have to shift left one to get rid 6139 of it. This only occurs if the number is -1 x 2^m. */ 6140 if (x[M+1] & 0x8000) 6141 { 6142 /* This is the case of -1 x 2^m, we have to rid ourselves of the 6143 high sign bit and shift the exponent. */ 6144 eshift(x, 1); 6145 i--; 6146 } 6147 } 6148 else 6149 { 6150 i = ((int) x[1]) - 0x7f; 6151 } 6152 6153 if ((i < -128) || (i > 127)) 6154 { 6155 y[0] |= 0xff7f; 6156 y[1] = 0xffff; 6157 if (mode != QFmode) 6158 { 6159 y[2] = 0xffff; 6160 y[3] = 0xffff; 6161 } 6162#ifdef ERANGE 6163 errno = ERANGE; 6164#endif 6165 return; 6166 } 6167 6168 y[0] |= ((i & 0xff) << 8); 6169 6170 eshift (x, 8); 6171 6172 y[0] |= x[M] & 0x7f; 6173 y[1] = x[M + 1]; 6174 if (mode != QFmode) 6175 { 6176 y[2] = x[M + 2]; 6177 y[3] = x[M + 3]; 6178 } 6179} 6180#endif /* C4X */ 6181 6182/* Output a binary NaN bit pattern in the target machine's format. */ 6183 6184/* If special NaN bit patterns are required, define them in tm.h 6185 as arrays of unsigned 16-bit shorts. Otherwise, use the default 6186 patterns here. */ 6187#ifdef TFMODE_NAN 6188TFMODE_NAN; 6189#else 6190#ifdef IEEE 6191unsigned EMUSHORT TFbignan[8] = 6192 {0x7fff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff}; 6193unsigned EMUSHORT TFlittlenan[8] = {0, 0, 0, 0, 0, 0, 0x8000, 0xffff}; 6194#endif 6195#endif 6196 6197#ifdef XFMODE_NAN 6198XFMODE_NAN; 6199#else 6200#ifdef IEEE 6201unsigned EMUSHORT XFbignan[6] = 6202 {0x7fff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff}; 6203unsigned EMUSHORT XFlittlenan[6] = {0, 0, 0, 0xc000, 0xffff, 0}; 6204#endif 6205#endif 6206 6207#ifdef DFMODE_NAN 6208DFMODE_NAN; 6209#else 6210#ifdef IEEE 6211unsigned EMUSHORT DFbignan[4] = {0x7fff, 0xffff, 0xffff, 0xffff}; 6212unsigned EMUSHORT DFlittlenan[4] = {0, 0, 0, 0xfff8}; 6213#endif 6214#endif 6215 6216#ifdef SFMODE_NAN 6217SFMODE_NAN; 6218#else 6219#ifdef IEEE 6220unsigned EMUSHORT SFbignan[2] = {0x7fff, 0xffff}; 6221unsigned EMUSHORT SFlittlenan[2] = {0, 0xffc0}; 6222#endif 6223#endif 6224 6225 6226static void 6227make_nan (nan, sign, mode) 6228 unsigned EMUSHORT *nan; 6229 int sign; 6230 enum machine_mode mode; 6231{ 6232 int n; 6233 unsigned EMUSHORT *p; 6234 6235 switch (mode) 6236 { 6237/* Possibly the `reserved operand' patterns on a VAX can be 6238 used like NaN's, but probably not in the same way as IEEE. */ 6239#if !defined(DEC) && !defined(IBM) && !defined(C4X) 6240 case TFmode: 6241 n = 8; 6242 if (REAL_WORDS_BIG_ENDIAN) 6243 p = TFbignan; 6244 else 6245 p = TFlittlenan; 6246 break; 6247 6248 case XFmode: 6249 n = 6; 6250 if (REAL_WORDS_BIG_ENDIAN) 6251 p = XFbignan; 6252 else 6253 p = XFlittlenan; 6254 break; 6255 6256 case DFmode: 6257 n = 4; 6258 if (REAL_WORDS_BIG_ENDIAN) 6259 p = DFbignan; 6260 else 6261 p = DFlittlenan; 6262 break; 6263 6264 case SFmode: 6265 case HFmode: 6266 n = 2; 6267 if (REAL_WORDS_BIG_ENDIAN) 6268 p = SFbignan; 6269 else 6270 p = SFlittlenan; 6271 break; 6272#endif 6273 6274 default: 6275 abort (); 6276 } 6277 if (REAL_WORDS_BIG_ENDIAN) 6278 *nan++ = (sign << 15) | (*p++ & 0x7fff); 6279 while (--n != 0) 6280 *nan++ = *p++; 6281 if (! REAL_WORDS_BIG_ENDIAN) 6282 *nan = (sign << 15) | (*p & 0x7fff); 6283} 6284 6285/* This is the inverse of the function `etarsingle' invoked by 6286 REAL_VALUE_TO_TARGET_SINGLE. */ 6287 6288REAL_VALUE_TYPE 6289ereal_unto_float (f) 6290 long f; 6291{ 6292 REAL_VALUE_TYPE r; 6293 unsigned EMUSHORT s[2]; 6294 unsigned EMUSHORT e[NE]; 6295 6296 /* Convert 32 bit integer to array of 16 bit pieces in target machine order. 6297 This is the inverse operation to what the function `endian' does. */ 6298 if (REAL_WORDS_BIG_ENDIAN) 6299 { 6300 s[0] = (unsigned EMUSHORT) (f >> 16); 6301 s[1] = (unsigned EMUSHORT) f; 6302 } 6303 else 6304 { 6305 s[0] = (unsigned EMUSHORT) f; 6306 s[1] = (unsigned EMUSHORT) (f >> 16); 6307 } 6308 /* Convert and promote the target float to E-type. */ 6309 e24toe (s, e); 6310 /* Output E-type to REAL_VALUE_TYPE. */ 6311 PUT_REAL (e, &r); 6312 return r; 6313} 6314 6315 6316/* This is the inverse of the function `etardouble' invoked by 6317 REAL_VALUE_TO_TARGET_DOUBLE. */ 6318 6319REAL_VALUE_TYPE 6320ereal_unto_double (d) 6321 long d[]; 6322{ 6323 REAL_VALUE_TYPE r; 6324 unsigned EMUSHORT s[4]; 6325 unsigned EMUSHORT e[NE]; 6326 6327 /* Convert array of HOST_WIDE_INT to equivalent array of 16-bit pieces. */ 6328 if (REAL_WORDS_BIG_ENDIAN) 6329 { 6330 s[0] = (unsigned EMUSHORT) (d[0] >> 16); 6331 s[1] = (unsigned EMUSHORT) d[0]; 6332 s[2] = (unsigned EMUSHORT) (d[1] >> 16); 6333 s[3] = (unsigned EMUSHORT) d[1]; 6334 } 6335 else 6336 { 6337 /* Target float words are little-endian. */ 6338 s[0] = (unsigned EMUSHORT) d[0]; 6339 s[1] = (unsigned EMUSHORT) (d[0] >> 16); 6340 s[2] = (unsigned EMUSHORT) d[1]; 6341 s[3] = (unsigned EMUSHORT) (d[1] >> 16); 6342 } 6343 /* Convert target double to E-type. */ 6344 e53toe (s, e); 6345 /* Output E-type to REAL_VALUE_TYPE. */ 6346 PUT_REAL (e, &r); 6347 return r; 6348} 6349 6350 6351/* Convert an SFmode target `float' value to a REAL_VALUE_TYPE. 6352 This is somewhat like ereal_unto_float, but the input types 6353 for these are different. */ 6354 6355REAL_VALUE_TYPE 6356ereal_from_float (f) 6357 HOST_WIDE_INT f; 6358{ 6359 REAL_VALUE_TYPE r; 6360 unsigned EMUSHORT s[2]; 6361 unsigned EMUSHORT e[NE]; 6362 6363 /* Convert 32 bit integer to array of 16 bit pieces in target machine order. 6364 This is the inverse operation to what the function `endian' does. */ 6365 if (REAL_WORDS_BIG_ENDIAN) 6366 { 6367 s[0] = (unsigned EMUSHORT) (f >> 16); 6368 s[1] = (unsigned EMUSHORT) f; 6369 } 6370 else 6371 { 6372 s[0] = (unsigned EMUSHORT) f; 6373 s[1] = (unsigned EMUSHORT) (f >> 16); 6374 } 6375 /* Convert and promote the target float to E-type. */ 6376 e24toe (s, e); 6377 /* Output E-type to REAL_VALUE_TYPE. */ 6378 PUT_REAL (e, &r); 6379 return r; 6380} 6381 6382 6383/* Convert a DFmode target `double' value to a REAL_VALUE_TYPE. 6384 This is somewhat like ereal_unto_double, but the input types 6385 for these are different. 6386 6387 The DFmode is stored as an array of HOST_WIDE_INT in the target's 6388 data format, with no holes in the bit packing. The first element 6389 of the input array holds the bits that would come first in the 6390 target computer's memory. */ 6391 6392REAL_VALUE_TYPE 6393ereal_from_double (d) 6394 HOST_WIDE_INT d[]; 6395{ 6396 REAL_VALUE_TYPE r; 6397 unsigned EMUSHORT s[4]; 6398 unsigned EMUSHORT e[NE]; 6399 6400 /* Convert array of HOST_WIDE_INT to equivalent array of 16-bit pieces. */ 6401 if (REAL_WORDS_BIG_ENDIAN) 6402 { 6403#if HOST_BITS_PER_WIDE_INT == 32 6404 s[0] = (unsigned EMUSHORT) (d[0] >> 16); 6405 s[1] = (unsigned EMUSHORT) d[0]; 6406 s[2] = (unsigned EMUSHORT) (d[1] >> 16); 6407 s[3] = (unsigned EMUSHORT) d[1]; 6408#else 6409 /* In this case the entire target double is contained in the 6410 first array element. The second element of the input is 6411 ignored. */ 6412 s[0] = (unsigned EMUSHORT) (d[0] >> 48); 6413 s[1] = (unsigned EMUSHORT) (d[0] >> 32); 6414 s[2] = (unsigned EMUSHORT) (d[0] >> 16); 6415 s[3] = (unsigned EMUSHORT) d[0]; 6416#endif 6417 } 6418 else 6419 { 6420 /* Target float words are little-endian. */ 6421 s[0] = (unsigned EMUSHORT) d[0]; 6422 s[1] = (unsigned EMUSHORT) (d[0] >> 16); 6423#if HOST_BITS_PER_WIDE_INT == 32 6424 s[2] = (unsigned EMUSHORT) d[1]; 6425 s[3] = (unsigned EMUSHORT) (d[1] >> 16); 6426#else 6427 s[2] = (unsigned EMUSHORT) (d[0] >> 32); 6428 s[3] = (unsigned EMUSHORT) (d[0] >> 48); 6429#endif 6430 } 6431 /* Convert target double to E-type. */ 6432 e53toe (s, e); 6433 /* Output E-type to REAL_VALUE_TYPE. */ 6434 PUT_REAL (e, &r); 6435 return r; 6436} 6437 6438 6439#if 0 6440/* Convert target computer unsigned 64-bit integer to e-type. 6441 The endian-ness of DImode follows the convention for integers, 6442 so we use WORDS_BIG_ENDIAN here, not REAL_WORDS_BIG_ENDIAN. */ 6443 6444static void 6445uditoe (di, e) 6446 unsigned EMUSHORT *di; /* Address of the 64-bit int. */ 6447 unsigned EMUSHORT *e; 6448{ 6449 unsigned EMUSHORT yi[NI]; 6450 int k; 6451 6452 ecleaz (yi); 6453 if (WORDS_BIG_ENDIAN) 6454 { 6455 for (k = M; k < M + 4; k++) 6456 yi[k] = *di++; 6457 } 6458 else 6459 { 6460 for (k = M + 3; k >= M; k--) 6461 yi[k] = *di++; 6462 } 6463 yi[E] = EXONE + 47; /* exponent if normalize shift count were 0 */ 6464 if ((k = enormlz (yi)) > NBITS)/* normalize the significand */ 6465 ecleaz (yi); /* it was zero */ 6466 else 6467 yi[E] -= (unsigned EMUSHORT) k;/* subtract shift count from exponent */ 6468 emovo (yi, e); 6469} 6470 6471/* Convert target computer signed 64-bit integer to e-type. */ 6472 6473static void 6474ditoe (di, e) 6475 unsigned EMUSHORT *di; /* Address of the 64-bit int. */ 6476 unsigned EMUSHORT *e; 6477{ 6478 unsigned EMULONG acc; 6479 unsigned EMUSHORT yi[NI]; 6480 unsigned EMUSHORT carry; 6481 int k, sign; 6482 6483 ecleaz (yi); 6484 if (WORDS_BIG_ENDIAN) 6485 { 6486 for (k = M; k < M + 4; k++) 6487 yi[k] = *di++; 6488 } 6489 else 6490 { 6491 for (k = M + 3; k >= M; k--) 6492 yi[k] = *di++; 6493 } 6494 /* Take absolute value */ 6495 sign = 0; 6496 if (yi[M] & 0x8000) 6497 { 6498 sign = 1; 6499 carry = 0; 6500 for (k = M + 3; k >= M; k--) 6501 { 6502 acc = (unsigned EMULONG) (~yi[k] & 0xffff) + carry; 6503 yi[k] = acc; 6504 carry = 0; 6505 if (acc & 0x10000) 6506 carry = 1; 6507 } 6508 } 6509 yi[E] = EXONE + 47; /* exponent if normalize shift count were 0 */ 6510 if ((k = enormlz (yi)) > NBITS)/* normalize the significand */ 6511 ecleaz (yi); /* it was zero */ 6512 else 6513 yi[E] -= (unsigned EMUSHORT) k;/* subtract shift count from exponent */ 6514 emovo (yi, e); 6515 if (sign) 6516 eneg (e); 6517} 6518 6519 6520/* Convert e-type to unsigned 64-bit int. */ 6521 6522static void 6523etoudi (x, i) 6524 unsigned EMUSHORT *x; 6525 unsigned EMUSHORT *i; 6526{ 6527 unsigned EMUSHORT xi[NI]; 6528 int j, k; 6529 6530 emovi (x, xi); 6531 if (xi[0]) 6532 { 6533 xi[M] = 0; 6534 goto noshift; 6535 } 6536 k = (int) xi[E] - (EXONE - 1); 6537 if (k <= 0) 6538 { 6539 for (j = 0; j < 4; j++) 6540 *i++ = 0; 6541 return; 6542 } 6543 if (k > 64) 6544 { 6545 for (j = 0; j < 4; j++) 6546 *i++ = 0xffff; 6547 if (extra_warnings) 6548 warning ("overflow on truncation to integer"); 6549 return; 6550 } 6551 if (k > 16) 6552 { 6553 /* Shift more than 16 bits: first shift up k-16 mod 16, 6554 then shift up by 16's. */ 6555 j = k - ((k >> 4) << 4); 6556 if (j == 0) 6557 j = 16; 6558 eshift (xi, j); 6559 if (WORDS_BIG_ENDIAN) 6560 *i++ = xi[M]; 6561 else 6562 { 6563 i += 3; 6564 *i-- = xi[M]; 6565 } 6566 k -= j; 6567 do 6568 { 6569 eshup6 (xi); 6570 if (WORDS_BIG_ENDIAN) 6571 *i++ = xi[M]; 6572 else 6573 *i-- = xi[M]; 6574 } 6575 while ((k -= 16) > 0); 6576 } 6577 else 6578 { 6579 /* shift not more than 16 bits */ 6580 eshift (xi, k); 6581 6582noshift: 6583 6584 if (WORDS_BIG_ENDIAN) 6585 { 6586 i += 3; 6587 *i-- = xi[M]; 6588 *i-- = 0; 6589 *i-- = 0; 6590 *i = 0; 6591 } 6592 else 6593 { 6594 *i++ = xi[M]; 6595 *i++ = 0; 6596 *i++ = 0; 6597 *i = 0; 6598 } 6599 } 6600} 6601 6602 6603/* Convert e-type to signed 64-bit int. */ 6604 6605static void 6606etodi (x, i) 6607 unsigned EMUSHORT *x; 6608 unsigned EMUSHORT *i; 6609{ 6610 unsigned EMULONG acc; 6611 unsigned EMUSHORT xi[NI]; 6612 unsigned EMUSHORT carry; 6613 unsigned EMUSHORT *isave; 6614 int j, k; 6615 6616 emovi (x, xi); 6617 k = (int) xi[E] - (EXONE - 1); 6618 if (k <= 0) 6619 { 6620 for (j = 0; j < 4; j++) 6621 *i++ = 0; 6622 return; 6623 } 6624 if (k > 64) 6625 { 6626 for (j = 0; j < 4; j++) 6627 *i++ = 0xffff; 6628 if (extra_warnings) 6629 warning ("overflow on truncation to integer"); 6630 return; 6631 } 6632 isave = i; 6633 if (k > 16) 6634 { 6635 /* Shift more than 16 bits: first shift up k-16 mod 16, 6636 then shift up by 16's. */ 6637 j = k - ((k >> 4) << 4); 6638 if (j == 0) 6639 j = 16; 6640 eshift (xi, j); 6641 if (WORDS_BIG_ENDIAN) 6642 *i++ = xi[M]; 6643 else 6644 { 6645 i += 3; 6646 *i-- = xi[M]; 6647 } 6648 k -= j; 6649 do 6650 { 6651 eshup6 (xi); 6652 if (WORDS_BIG_ENDIAN) 6653 *i++ = xi[M]; 6654 else 6655 *i-- = xi[M]; 6656 } 6657 while ((k -= 16) > 0); 6658 } 6659 else 6660 { 6661 /* shift not more than 16 bits */ 6662 eshift (xi, k); 6663 6664 if (WORDS_BIG_ENDIAN) 6665 { 6666 i += 3; 6667 *i = xi[M]; 6668 *i-- = 0; 6669 *i-- = 0; 6670 *i = 0; 6671 } 6672 else 6673 { 6674 *i++ = xi[M]; 6675 *i++ = 0; 6676 *i++ = 0; 6677 *i = 0; 6678 } 6679 } 6680 /* Negate if negative */ 6681 if (xi[0]) 6682 { 6683 carry = 0; 6684 if (WORDS_BIG_ENDIAN) 6685 isave += 3; 6686 for (k = 0; k < 4; k++) 6687 { 6688 acc = (unsigned EMULONG) (~(*isave) & 0xffff) + carry; 6689 if (WORDS_BIG_ENDIAN) 6690 *isave-- = acc; 6691 else 6692 *isave++ = acc; 6693 carry = 0; 6694 if (acc & 0x10000) 6695 carry = 1; 6696 } 6697 } 6698} 6699 6700 6701/* Longhand square root routine. */ 6702 6703 6704static int esqinited = 0; 6705static unsigned short sqrndbit[NI]; 6706 6707static void 6708esqrt (x, y) 6709 unsigned EMUSHORT *x, *y; 6710{ 6711 unsigned EMUSHORT temp[NI], num[NI], sq[NI], xx[NI]; 6712 EMULONG m, exp; 6713 int i, j, k, n, nlups; 6714 6715 if (esqinited == 0) 6716 { 6717 ecleaz (sqrndbit); 6718 sqrndbit[NI - 2] = 1; 6719 esqinited = 1; 6720 } 6721 /* Check for arg <= 0 */ 6722 i = ecmp (x, ezero); 6723 if (i <= 0) 6724 { 6725 if (i == -1) 6726 { 6727 mtherr ("esqrt", DOMAIN); 6728 eclear (y); 6729 } 6730 else 6731 emov (x, y); 6732 return; 6733 } 6734 6735#ifdef INFINITY 6736 if (eisinf (x)) 6737 { 6738 eclear (y); 6739 einfin (y); 6740 return; 6741 } 6742#endif 6743 /* Bring in the arg and renormalize if it is denormal. */ 6744 emovi (x, xx); 6745 m = (EMULONG) xx[1]; /* local long word exponent */ 6746 if (m == 0) 6747 m -= enormlz (xx); 6748 6749 /* Divide exponent by 2 */ 6750 m -= 0x3ffe; 6751 exp = (unsigned short) ((m / 2) + 0x3ffe); 6752 6753 /* Adjust if exponent odd */ 6754 if ((m & 1) != 0) 6755 { 6756 if (m > 0) 6757 exp += 1; 6758 eshdn1 (xx); 6759 } 6760 6761 ecleaz (sq); 6762 ecleaz (num); 6763 n = 8; /* get 8 bits of result per inner loop */ 6764 nlups = rndprc; 6765 j = 0; 6766 6767 while (nlups > 0) 6768 { 6769 /* bring in next word of arg */ 6770 if (j < NE) 6771 num[NI - 1] = xx[j + 3]; 6772 /* Do additional bit on last outer loop, for roundoff. */ 6773 if (nlups <= 8) 6774 n = nlups + 1; 6775 for (i = 0; i < n; i++) 6776 { 6777 /* Next 2 bits of arg */ 6778 eshup1 (num); 6779 eshup1 (num); 6780 /* Shift up answer */ 6781 eshup1 (sq); 6782 /* Make trial divisor */ 6783 for (k = 0; k < NI; k++) 6784 temp[k] = sq[k]; 6785 eshup1 (temp); 6786 eaddm (sqrndbit, temp); 6787 /* Subtract and insert answer bit if it goes in */ 6788 if (ecmpm (temp, num) <= 0) 6789 { 6790 esubm (temp, num); 6791 sq[NI - 2] |= 1; 6792 } 6793 } 6794 nlups -= n; 6795 j += 1; 6796 } 6797 6798 /* Adjust for extra, roundoff loop done. */ 6799 exp += (NBITS - 1) - rndprc; 6800 6801 /* Sticky bit = 1 if the remainder is nonzero. */ 6802 k = 0; 6803 for (i = 3; i < NI; i++) 6804 k |= (int) num[i]; 6805 6806 /* Renormalize and round off. */ 6807 emdnorm (sq, k, 0, exp, 64); 6808 emovo (sq, y); 6809} 6810#endif 6811#endif /* EMU_NON_COMPILE not defined */ 6812 6813/* Return the binary precision of the significand for a given 6814 floating point mode. The mode can hold an integer value 6815 that many bits wide, without losing any bits. */ 6816 6817int 6818significand_size (mode) 6819 enum machine_mode mode; 6820{ 6821 6822/* Don't test the modes, but their sizes, lest this 6823 code won't work for BITS_PER_UNIT != 8 . */ 6824 6825switch (GET_MODE_BITSIZE (mode)) 6826 { 6827 case 32: 6828 6829#if TARGET_FLOAT_FORMAT == C4X_FLOAT_FORMAT 6830 return 56; 6831#endif 6832 6833 return 24; 6834 6835 case 64: 6836#if TARGET_FLOAT_FORMAT == IEEE_FLOAT_FORMAT 6837 return 53; 6838#else 6839#if TARGET_FLOAT_FORMAT == IBM_FLOAT_FORMAT 6840 return 56; 6841#else 6842#if TARGET_FLOAT_FORMAT == VAX_FLOAT_FORMAT 6843 return 56; 6844#else 6845#if TARGET_FLOAT_FORMAT == C4X_FLOAT_FORMAT 6846 return 56; 6847#else 6848 abort (); 6849#endif 6850#endif 6851#endif 6852#endif 6853 6854 case 96: 6855 return 64; 6856 case 128: 6857 return 113; 6858 6859 default: 6860 abort (); 6861 } 6862} 6863