real.c revision 90075
11722Sjkh/* real.c - implementation of REAL_ARITHMETIC, REAL_VALUE_ATOF, 21722Sjkh and support for XFmode IEEE extended real floating point arithmetic. 31722Sjkh Copyright (C) 1993, 1994, 1995, 1996, 1997, 1998, 41722Sjkh 1999, 2000, 2002 Free Software Foundation, Inc. 51722Sjkh Contributed by Stephen L. Moshier (moshier@world.std.com). 61722Sjkh 71722SjkhThis file is part of GCC. 81722Sjkh 91722SjkhGCC is free software; you can redistribute it and/or modify it under 101722Sjkhthe terms of the GNU General Public License as published by the Free 111722SjkhSoftware Foundation; either version 2, or (at your option) any later 121722Sjkhversion. 131722Sjkh 141722SjkhGCC is distributed in the hope that it will be useful, but WITHOUT ANY 151722SjkhWARRANTY; without even the implied warranty of MERCHANTABILITY or 161722SjkhFITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License 171722Sjkhfor more details. 181722Sjkh 191722SjkhYou should have received a copy of the GNU General Public License 201722Sjkhalong with GCC; see the file COPYING. If not, write to the Free 211722SjkhSoftware Foundation, 59 Temple Place - Suite 330, Boston, MA 221722Sjkh02111-1307, USA. */ 231722Sjkh 241722Sjkh#include "config.h" 251722Sjkh#include "system.h" 261722Sjkh#include "tree.h" 278857Srgrimes#include "toplev.h" 281722Sjkh#include "tm_p.h" 291722Sjkh 301722Sjkh/* To enable support of XFmode extended real floating point, define 311722SjkhLONG_DOUBLE_TYPE_SIZE 96 in the tm.h file (m68k.h or i386.h). 321722Sjkh 331722SjkhTo support cross compilation between IEEE, VAX and IBM floating 34237625Sobrienpoint formats, define REAL_ARITHMETIC in the tm.h file. 35237625Sobrien 36237625SobrienIn either case the machine files (tm.h) must not contain any code 37237625Sobrienthat tries to use host floating point arithmetic to convert 381722SjkhREAL_VALUE_TYPEs from `double' to `float', pass them to fprintf, 39173067Syaretc. In cross-compile situations a REAL_VALUE_TYPE may not 401722Sjkhbe intelligible to the host computer's native arithmetic. 411722Sjkh 421722SjkhThe emulator defaults to the host's floating point format so that 431722Sjkhits decimal conversion functions can be used if desired (see 441722Sjkhreal.h). 451722Sjkh 461722SjkhThe first part of this file interfaces gcc to a floating point 47173065Syararithmetic suite that was not written with gcc in mind. Avoid 481722Sjkhchanging the low-level arithmetic routines unless you have suitable 491722Sjkhtest programs available. A special version of the PARANOIA floating 50237625Sobrienpoint arithmetic tester, modified for this purpose, can be found on 51237625Sobrienusc.edu: /pub/C-numanal/ieeetest.zoo. Other tests, and libraries of 521722SjkhXFmode and TFmode transcendental functions, can be obtained by ftp from 531722Sjkhnetlib.att.com: netlib/cephes. */ 541722Sjkh 551722Sjkh/* Type of computer arithmetic. 561722Sjkh Only one of DEC, IBM, IEEE, C4X, or UNK should get defined. 571722Sjkh 581722Sjkh `IEEE', when REAL_WORDS_BIG_ENDIAN is non-zero, refers generically 591722Sjkh to big-endian IEEE floating-point data structure. This definition 601722Sjkh should work in SFmode `float' type and DFmode `double' type on 611722Sjkh virtually all big-endian IEEE machines. If LONG_DOUBLE_TYPE_SIZE 621722Sjkh has been defined to be 96, then IEEE also invokes the particular 631722Sjkh XFmode (`long double' type) data structure used by the Motorola 641722Sjkh 680x0 series processors. 651722Sjkh 66101341Sgshapiro `IEEE', when REAL_WORDS_BIG_ENDIAN is zero, refers generally to 671722Sjkh little-endian IEEE machines. In this case, if LONG_DOUBLE_TYPE_SIZE 681722Sjkh has been defined to be 96, then IEEE also invokes the particular 691722Sjkh XFmode `long double' data structure used by the Intel 80x86 series 701722Sjkh processors. 711722Sjkh 721722Sjkh `DEC' refers specifically to the Digital Equipment Corp PDP-11 731722Sjkh and VAX floating point data structure. This model currently 74237625Sobrien supports no type wider than DFmode. 75237625Sobrien 768664Sphk `IBM' refers specifically to the IBM System/370 and compatible 778664Sphk floating point data structure. This model currently supports 788664Sphk no type wider than DFmode. The IBM conversions were contributed by 798664Sphk frank@atom.ansto.gov.au (Frank Crawford). 808664Sphk 818664Sphk `C4X' refers specifically to the floating point format used on 828664Sphk Texas Instruments TMS320C3x and TMS320C4x digital signal 838664Sphk processors. This supports QFmode (32-bit float, double) and HFmode 848664Sphk (40-bit long double) where BITS_PER_BYTE is 32. Unlike IEEE 858664Sphk floats, C4x floats are not rounded to be even. The C4x conversions 868664Sphk were contributed by m.hayes@elec.canterbury.ac.nz (Michael Hayes) and 878664Sphk Haj.Ten.Brugge@net.HCC.nl (Herman ten Brugge). 888664Sphk 898664Sphk If LONG_DOUBLE_TYPE_SIZE = 64 (the default, unless tm.h defines it) 90237625Sobrien then `long double' and `double' are both implemented, but they 91237625Sobrien both mean DFmode. In this case, the software floating-point 921722Sjkh support available here is activated by writing 93173065Syar #define REAL_ARITHMETIC 941722Sjkh in tm.h. 951722Sjkh 961722Sjkh The case LONG_DOUBLE_TYPE_SIZE = 128 activates TFmode support 978857Srgrimes and may deactivate XFmode since `long double' is used to refer 981722Sjkh to both modes. Defining INTEL_EXTENDED_IEEE_FORMAT to non-zero 991722Sjkh at the same time enables 80387-style 80-bit floats in a 128-bit 100173065Syar padded image, as seen on IA-64. 101173065Syar 102173065Syar The macros FLOAT_WORDS_BIG_ENDIAN, HOST_FLOAT_WORDS_BIG_ENDIAN, 103101550Sgshapiro contributed by Richard Earnshaw <Richard.Earnshaw@cl.cam.ac.uk>, 1041722Sjkh separate the floating point unit's endian-ness from that of 1051722Sjkh the integer addressing. This permits one to define a big-endian 1061722Sjkh FPU on a little-endian machine (e.g., ARM). An extension to 107237625Sobrien BYTES_BIG_ENDIAN may be required for some machines in the future. 108237625Sobrien These optional macros may be defined in tm.h. In real.h, they 1091722Sjkh default to WORDS_BIG_ENDIAN, etc., so there is no need to define 1101722Sjkh them for any normal host or target machine on which the floats 1111722Sjkh and the integers have the same endian-ness. */ 1121722Sjkh 11329453Scharnier 1141722Sjkh/* The following converts gcc macros into the ones used by this file. */ 1151722Sjkh 1161722Sjkh/* REAL_ARITHMETIC defined means that macros in real.h are 1171722Sjkh defined to call emulator functions. */ 1181722Sjkh#ifdef REAL_ARITHMETIC 1191722Sjkh 1201722Sjkh#if TARGET_FLOAT_FORMAT == VAX_FLOAT_FORMAT 1211722Sjkh/* PDP-11, Pro350, VAX: */ 1221722Sjkh#define DEC 1 1231722Sjkh#else /* it's not VAX */ 1241722Sjkh#if TARGET_FLOAT_FORMAT == IBM_FLOAT_FORMAT 1251722Sjkh/* IBM System/370 style */ 1261722Sjkh#define IBM 1 1271722Sjkh#else /* it's also not an IBM */ 1281722Sjkh#if TARGET_FLOAT_FORMAT == C4X_FLOAT_FORMAT 1291722Sjkh/* TMS320C3x/C4x style */ 1301722Sjkh#define C4X 1 131#else /* it's also not a C4X */ 132#if TARGET_FLOAT_FORMAT == IEEE_FLOAT_FORMAT 133#define IEEE 134#else /* it's not IEEE either */ 135/* UNKnown arithmetic. We don't support this and can't go on. */ 136unknown arithmetic type 137#define UNK 1 138#endif /* not IEEE */ 139#endif /* not C4X */ 140#endif /* not IBM */ 141#endif /* not VAX */ 142 143#define REAL_WORDS_BIG_ENDIAN FLOAT_WORDS_BIG_ENDIAN 144 145#else 146/* REAL_ARITHMETIC not defined means that the *host's* data 147 structure will be used. It may differ by endian-ness from the 148 target machine's structure and will get its ends swapped 149 accordingly (but not here). Probably only the decimal <-> binary 150 functions in this file will actually be used in this case. */ 151 152#if HOST_FLOAT_FORMAT == VAX_FLOAT_FORMAT 153#define DEC 1 154#else /* it's not VAX */ 155#if HOST_FLOAT_FORMAT == IBM_FLOAT_FORMAT 156/* IBM System/370 style */ 157#define IBM 1 158#else /* it's also not an IBM */ 159#if HOST_FLOAT_FORMAT == IEEE_FLOAT_FORMAT 160#define IEEE 161#else /* it's not IEEE either */ 162unknown arithmetic type 163#define UNK 1 164#endif /* not IEEE */ 165#endif /* not IBM */ 166#endif /* not VAX */ 167 168#define REAL_WORDS_BIG_ENDIAN HOST_FLOAT_WORDS_BIG_ENDIAN 169 170#endif /* REAL_ARITHMETIC not defined */ 171 172/* Define INFINITY for support of infinity. 173 Define NANS for support of Not-a-Number's (NaN's). */ 174#if !defined(DEC) && !defined(IBM) && !defined(C4X) 175#define INFINITY 176#define NANS 177#endif 178 179/* Support of NaNs requires support of infinity. */ 180#ifdef NANS 181#ifndef INFINITY 182#define INFINITY 183#endif 184#endif 185 186/* Find a host integer type that is at least 16 bits wide, 187 and another type at least twice whatever that size is. */ 188 189#if HOST_BITS_PER_CHAR >= 16 190#define EMUSHORT char 191#define EMUSHORT_SIZE HOST_BITS_PER_CHAR 192#define EMULONG_SIZE (2 * HOST_BITS_PER_CHAR) 193#else 194#if HOST_BITS_PER_SHORT >= 16 195#define EMUSHORT short 196#define EMUSHORT_SIZE HOST_BITS_PER_SHORT 197#define EMULONG_SIZE (2 * HOST_BITS_PER_SHORT) 198#else 199#if HOST_BITS_PER_INT >= 16 200#define EMUSHORT int 201#define EMUSHORT_SIZE HOST_BITS_PER_INT 202#define EMULONG_SIZE (2 * HOST_BITS_PER_INT) 203#else 204#if HOST_BITS_PER_LONG >= 16 205#define EMUSHORT long 206#define EMUSHORT_SIZE HOST_BITS_PER_LONG 207#define EMULONG_SIZE (2 * HOST_BITS_PER_LONG) 208#else 209/* You will have to modify this program to have a smaller unit size. */ 210#define EMU_NON_COMPILE 211#endif 212#endif 213#endif 214#endif 215 216/* If no 16-bit type has been found and the compiler is GCC, try HImode. */ 217#if defined(__GNUC__) && EMUSHORT_SIZE != 16 218typedef int HItype __attribute__ ((mode (HI))); 219typedef unsigned int UHItype __attribute__ ((mode (HI))); 220#undef EMUSHORT 221#undef EMUSHORT_SIZE 222#undef EMULONG_SIZE 223#define EMUSHORT HItype 224#define UEMUSHORT UHItype 225#define EMUSHORT_SIZE 16 226#define EMULONG_SIZE 32 227#else 228#define UEMUSHORT unsigned EMUSHORT 229#endif 230 231#if HOST_BITS_PER_SHORT >= EMULONG_SIZE 232#define EMULONG short 233#else 234#if HOST_BITS_PER_INT >= EMULONG_SIZE 235#define EMULONG int 236#else 237#if HOST_BITS_PER_LONG >= EMULONG_SIZE 238#define EMULONG long 239#else 240#if HOST_BITS_PER_LONGLONG >= EMULONG_SIZE 241#define EMULONG long long int 242#else 243/* You will have to modify this program to have a smaller unit size. */ 244#define EMU_NON_COMPILE 245#endif 246#endif 247#endif 248#endif 249 250 251/* The host interface doesn't work if no 16-bit size exists. */ 252#if EMUSHORT_SIZE != 16 253#define EMU_NON_COMPILE 254#endif 255 256/* OK to continue compilation. */ 257#ifndef EMU_NON_COMPILE 258 259/* Construct macros to translate between REAL_VALUE_TYPE and e type. 260 In GET_REAL and PUT_REAL, r and e are pointers. 261 A REAL_VALUE_TYPE is guaranteed to occupy contiguous locations 262 in memory, with no holes. */ 263 264#if MAX_LONG_DOUBLE_TYPE_SIZE == 96 || \ 265 ((INTEL_EXTENDED_IEEE_FORMAT != 0) && MAX_LONG_DOUBLE_TYPE_SIZE == 128) 266/* Number of 16 bit words in external e type format */ 267# define NE 6 268# define MAXDECEXP 4932 269# define MINDECEXP -4956 270# define GET_REAL(r,e) memcpy ((e), (r), 2*NE) 271# define PUT_REAL(e,r) \ 272 do { \ 273 memcpy ((r), (e), 2*NE); \ 274 if (2*NE < sizeof (*r)) \ 275 memset ((char *) (r) + 2*NE, 0, sizeof (*r) - 2*NE); \ 276 } while (0) 277# else /* no XFmode */ 278# if MAX_LONG_DOUBLE_TYPE_SIZE == 128 279# define NE 10 280# define MAXDECEXP 4932 281# define MINDECEXP -4977 282# define GET_REAL(r,e) memcpy ((e), (r), 2*NE) 283# define PUT_REAL(e,r) \ 284 do { \ 285 memcpy ((r), (e), 2*NE); \ 286 if (2*NE < sizeof (*r)) \ 287 memset ((char *) (r) + 2*NE, 0, sizeof (*r) - 2*NE); \ 288 } while (0) 289#else 290#define NE 6 291#define MAXDECEXP 4932 292#define MINDECEXP -4956 293#ifdef REAL_ARITHMETIC 294/* Emulator uses target format internally 295 but host stores it in host endian-ness. */ 296 297#define GET_REAL(r,e) \ 298do { \ 299 if (HOST_FLOAT_WORDS_BIG_ENDIAN == REAL_WORDS_BIG_ENDIAN) \ 300 e53toe ((const UEMUSHORT *) (r), (e)); \ 301 else \ 302 { \ 303 UEMUSHORT w[4]; \ 304 memcpy (&w[3], ((const EMUSHORT *) r), sizeof (EMUSHORT)); \ 305 memcpy (&w[2], ((const EMUSHORT *) r) + 1, sizeof (EMUSHORT)); \ 306 memcpy (&w[1], ((const EMUSHORT *) r) + 2, sizeof (EMUSHORT)); \ 307 memcpy (&w[0], ((const EMUSHORT *) r) + 3, sizeof (EMUSHORT)); \ 308 e53toe (w, (e)); \ 309 } \ 310 } while (0) 311 312#define PUT_REAL(e,r) \ 313do { \ 314 if (HOST_FLOAT_WORDS_BIG_ENDIAN == REAL_WORDS_BIG_ENDIAN) \ 315 etoe53 ((e), (UEMUSHORT *) (r)); \ 316 else \ 317 { \ 318 UEMUSHORT w[4]; \ 319 etoe53 ((e), w); \ 320 memcpy (((EMUSHORT *) r), &w[3], sizeof (EMUSHORT)); \ 321 memcpy (((EMUSHORT *) r) + 1, &w[2], sizeof (EMUSHORT)); \ 322 memcpy (((EMUSHORT *) r) + 2, &w[1], sizeof (EMUSHORT)); \ 323 memcpy (((EMUSHORT *) r) + 3, &w[0], sizeof (EMUSHORT)); \ 324 } \ 325 } while (0) 326 327#else /* not REAL_ARITHMETIC */ 328 329/* emulator uses host format */ 330#define GET_REAL(r,e) e53toe ((const UEMUSHORT *) (r), (e)) 331#define PUT_REAL(e,r) etoe53 ((e), (UEMUSHORT *) (r)) 332 333#endif /* not REAL_ARITHMETIC */ 334#endif /* not TFmode */ 335#endif /* not XFmode */ 336 337 338/* Number of 16 bit words in internal format */ 339#define NI (NE+3) 340 341/* Array offset to exponent */ 342#define E 1 343 344/* Array offset to high guard word */ 345#define M 2 346 347/* Number of bits of precision */ 348#define NBITS ((NI-4)*16) 349 350/* Maximum number of decimal digits in ASCII conversion 351 * = NBITS*log10(2) 352 */ 353#define NDEC (NBITS*8/27) 354 355/* The exponent of 1.0 */ 356#define EXONE (0x3fff) 357 358#if defined(HOST_EBCDIC) 359/* bit 8 is significant in EBCDIC */ 360#define CHARMASK 0xff 361#else 362#define CHARMASK 0x7f 363#endif 364 365extern int extra_warnings; 366extern const UEMUSHORT ezero[NE], ehalf[NE], eone[NE], etwo[NE]; 367extern const UEMUSHORT elog2[NE], esqrt2[NE]; 368 369static void endian PARAMS ((const UEMUSHORT *, long *, 370 enum machine_mode)); 371static void eclear PARAMS ((UEMUSHORT *)); 372static void emov PARAMS ((const UEMUSHORT *, UEMUSHORT *)); 373#if 0 374static void eabs PARAMS ((UEMUSHORT *)); 375#endif 376static void eneg PARAMS ((UEMUSHORT *)); 377static int eisneg PARAMS ((const UEMUSHORT *)); 378static int eisinf PARAMS ((const UEMUSHORT *)); 379static int eisnan PARAMS ((const UEMUSHORT *)); 380static void einfin PARAMS ((UEMUSHORT *)); 381#ifdef NANS 382static void enan PARAMS ((UEMUSHORT *, int)); 383static void einan PARAMS ((UEMUSHORT *)); 384static int eiisnan PARAMS ((const UEMUSHORT *)); 385static int eiisneg PARAMS ((const UEMUSHORT *)); 386static void make_nan PARAMS ((UEMUSHORT *, int, enum machine_mode)); 387#endif 388static void emovi PARAMS ((const UEMUSHORT *, UEMUSHORT *)); 389static void emovo PARAMS ((const UEMUSHORT *, UEMUSHORT *)); 390static void ecleaz PARAMS ((UEMUSHORT *)); 391static void ecleazs PARAMS ((UEMUSHORT *)); 392static void emovz PARAMS ((const UEMUSHORT *, UEMUSHORT *)); 393#if 0 394static void eiinfin PARAMS ((UEMUSHORT *)); 395#endif 396#ifdef INFINITY 397static int eiisinf PARAMS ((const UEMUSHORT *)); 398#endif 399static int ecmpm PARAMS ((const UEMUSHORT *, const UEMUSHORT *)); 400static void eshdn1 PARAMS ((UEMUSHORT *)); 401static void eshup1 PARAMS ((UEMUSHORT *)); 402static void eshdn8 PARAMS ((UEMUSHORT *)); 403static void eshup8 PARAMS ((UEMUSHORT *)); 404static void eshup6 PARAMS ((UEMUSHORT *)); 405static void eshdn6 PARAMS ((UEMUSHORT *)); 406static void eaddm PARAMS ((const UEMUSHORT *, UEMUSHORT *)); 407static void esubm PARAMS ((const UEMUSHORT *, UEMUSHORT *)); 408static void m16m PARAMS ((unsigned int, const UEMUSHORT *, UEMUSHORT *)); 409static int edivm PARAMS ((const UEMUSHORT *, UEMUSHORT *)); 410static int emulm PARAMS ((const UEMUSHORT *, UEMUSHORT *)); 411static void emdnorm PARAMS ((UEMUSHORT *, int, int, EMULONG, int)); 412static void esub PARAMS ((const UEMUSHORT *, const UEMUSHORT *, 413 UEMUSHORT *)); 414static void eadd PARAMS ((const UEMUSHORT *, const UEMUSHORT *, 415 UEMUSHORT *)); 416static void eadd1 PARAMS ((const UEMUSHORT *, const UEMUSHORT *, 417 UEMUSHORT *)); 418static void ediv PARAMS ((const UEMUSHORT *, const UEMUSHORT *, 419 UEMUSHORT *)); 420static void emul PARAMS ((const UEMUSHORT *, const UEMUSHORT *, 421 UEMUSHORT *)); 422static void e53toe PARAMS ((const UEMUSHORT *, UEMUSHORT *)); 423static void e64toe PARAMS ((const UEMUSHORT *, UEMUSHORT *)); 424#if (INTEL_EXTENDED_IEEE_FORMAT == 0) 425static void e113toe PARAMS ((const UEMUSHORT *, UEMUSHORT *)); 426#endif 427static void e24toe PARAMS ((const UEMUSHORT *, UEMUSHORT *)); 428#if (INTEL_EXTENDED_IEEE_FORMAT == 0) 429static void etoe113 PARAMS ((const UEMUSHORT *, UEMUSHORT *)); 430static void toe113 PARAMS ((UEMUSHORT *, UEMUSHORT *)); 431#endif 432static void etoe64 PARAMS ((const UEMUSHORT *, UEMUSHORT *)); 433static void toe64 PARAMS ((UEMUSHORT *, UEMUSHORT *)); 434static void etoe53 PARAMS ((const UEMUSHORT *, UEMUSHORT *)); 435static void toe53 PARAMS ((UEMUSHORT *, UEMUSHORT *)); 436static void etoe24 PARAMS ((const UEMUSHORT *, UEMUSHORT *)); 437static void toe24 PARAMS ((UEMUSHORT *, UEMUSHORT *)); 438static int ecmp PARAMS ((const UEMUSHORT *, const UEMUSHORT *)); 439#if 0 440static void eround PARAMS ((const UEMUSHORT *, UEMUSHORT *)); 441#endif 442static void ltoe PARAMS ((const HOST_WIDE_INT *, UEMUSHORT *)); 443static void ultoe PARAMS ((const unsigned HOST_WIDE_INT *, UEMUSHORT *)); 444static void eifrac PARAMS ((const UEMUSHORT *, HOST_WIDE_INT *, 445 UEMUSHORT *)); 446static void euifrac PARAMS ((const UEMUSHORT *, unsigned HOST_WIDE_INT *, 447 UEMUSHORT *)); 448static int eshift PARAMS ((UEMUSHORT *, int)); 449static int enormlz PARAMS ((UEMUSHORT *)); 450#if 0 451static void e24toasc PARAMS ((const UEMUSHORT *, char *, int)); 452static void e53toasc PARAMS ((const UEMUSHORT *, char *, int)); 453static void e64toasc PARAMS ((const UEMUSHORT *, char *, int)); 454static void e113toasc PARAMS ((const UEMUSHORT *, char *, int)); 455#endif /* 0 */ 456static void etoasc PARAMS ((const UEMUSHORT *, char *, int)); 457static void asctoe24 PARAMS ((const char *, UEMUSHORT *)); 458static void asctoe53 PARAMS ((const char *, UEMUSHORT *)); 459static void asctoe64 PARAMS ((const char *, UEMUSHORT *)); 460#if (INTEL_EXTENDED_IEEE_FORMAT == 0) 461static void asctoe113 PARAMS ((const char *, UEMUSHORT *)); 462#endif 463static void asctoe PARAMS ((const char *, UEMUSHORT *)); 464static void asctoeg PARAMS ((const char *, UEMUSHORT *, int)); 465static void efloor PARAMS ((const UEMUSHORT *, UEMUSHORT *)); 466#if 0 467static void efrexp PARAMS ((const UEMUSHORT *, int *, 468 UEMUSHORT *)); 469#endif 470static void eldexp PARAMS ((const UEMUSHORT *, int, UEMUSHORT *)); 471#if 0 472static void eremain PARAMS ((const UEMUSHORT *, const UEMUSHORT *, 473 UEMUSHORT *)); 474#endif 475static void eiremain PARAMS ((UEMUSHORT *, UEMUSHORT *)); 476static void mtherr PARAMS ((const char *, int)); 477#ifdef DEC 478static void dectoe PARAMS ((const UEMUSHORT *, UEMUSHORT *)); 479static void etodec PARAMS ((const UEMUSHORT *, UEMUSHORT *)); 480static void todec PARAMS ((UEMUSHORT *, UEMUSHORT *)); 481#endif 482#ifdef IBM 483static void ibmtoe PARAMS ((const UEMUSHORT *, UEMUSHORT *, 484 enum machine_mode)); 485static void etoibm PARAMS ((const UEMUSHORT *, UEMUSHORT *, 486 enum machine_mode)); 487static void toibm PARAMS ((UEMUSHORT *, UEMUSHORT *, 488 enum machine_mode)); 489#endif 490#ifdef C4X 491static void c4xtoe PARAMS ((const UEMUSHORT *, UEMUSHORT *, 492 enum machine_mode)); 493static void etoc4x PARAMS ((const UEMUSHORT *, UEMUSHORT *, 494 enum machine_mode)); 495static void toc4x PARAMS ((UEMUSHORT *, UEMUSHORT *, 496 enum machine_mode)); 497#endif 498#if 0 499static void uditoe PARAMS ((const UEMUSHORT *, UEMUSHORT *)); 500static void ditoe PARAMS ((const UEMUSHORT *, UEMUSHORT *)); 501static void etoudi PARAMS ((const UEMUSHORT *, UEMUSHORT *)); 502static void etodi PARAMS ((const UEMUSHORT *, UEMUSHORT *)); 503static void esqrt PARAMS ((const UEMUSHORT *, UEMUSHORT *)); 504#endif 505 506/* Copy 32-bit numbers obtained from array containing 16-bit numbers, 507 swapping ends if required, into output array of longs. The 508 result is normally passed to fprintf by the ASM_OUTPUT_ macros. */ 509 510static void 511endian (e, x, mode) 512 const UEMUSHORT e[]; 513 long x[]; 514 enum machine_mode mode; 515{ 516 unsigned long th, t; 517 518 if (REAL_WORDS_BIG_ENDIAN) 519 { 520 switch (mode) 521 { 522 case TFmode: 523#if (INTEL_EXTENDED_IEEE_FORMAT == 0) 524 /* Swap halfwords in the fourth long. */ 525 th = (unsigned long) e[6] & 0xffff; 526 t = (unsigned long) e[7] & 0xffff; 527 t |= th << 16; 528 x[3] = (long) t; 529#else 530 x[3] = 0; 531#endif 532 /* FALLTHRU */ 533 534 case XFmode: 535 /* Swap halfwords in the third long. */ 536 th = (unsigned long) e[4] & 0xffff; 537 t = (unsigned long) e[5] & 0xffff; 538 t |= th << 16; 539 x[2] = (long) t; 540 /* FALLTHRU */ 541 542 case DFmode: 543 /* Swap halfwords in the second word. */ 544 th = (unsigned long) e[2] & 0xffff; 545 t = (unsigned long) e[3] & 0xffff; 546 t |= th << 16; 547 x[1] = (long) t; 548 /* FALLTHRU */ 549 550 case SFmode: 551 case HFmode: 552 /* Swap halfwords in the first word. */ 553 th = (unsigned long) e[0] & 0xffff; 554 t = (unsigned long) e[1] & 0xffff; 555 t |= th << 16; 556 x[0] = (long) t; 557 break; 558 559 default: 560 abort (); 561 } 562 } 563 else 564 { 565 /* Pack the output array without swapping. */ 566 567 switch (mode) 568 { 569 case TFmode: 570#if (INTEL_EXTENDED_IEEE_FORMAT == 0) 571 /* Pack the fourth long. */ 572 th = (unsigned long) e[7] & 0xffff; 573 t = (unsigned long) e[6] & 0xffff; 574 t |= th << 16; 575 x[3] = (long) t; 576#else 577 x[3] = 0; 578#endif 579 /* FALLTHRU */ 580 581 case XFmode: 582 /* Pack the third long. 583 Each element of the input REAL_VALUE_TYPE array has 16 useful bits 584 in it. */ 585 th = (unsigned long) e[5] & 0xffff; 586 t = (unsigned long) e[4] & 0xffff; 587 t |= th << 16; 588 x[2] = (long) t; 589 /* FALLTHRU */ 590 591 case DFmode: 592 /* Pack the second long */ 593 th = (unsigned long) e[3] & 0xffff; 594 t = (unsigned long) e[2] & 0xffff; 595 t |= th << 16; 596 x[1] = (long) t; 597 /* FALLTHRU */ 598 599 case SFmode: 600 case HFmode: 601 /* Pack the first long */ 602 th = (unsigned long) e[1] & 0xffff; 603 t = (unsigned long) e[0] & 0xffff; 604 t |= th << 16; 605 x[0] = (long) t; 606 break; 607 608 default: 609 abort (); 610 } 611 } 612} 613 614 615/* This is the implementation of the REAL_ARITHMETIC macro. */ 616 617void 618earith (value, icode, r1, r2) 619 REAL_VALUE_TYPE *value; 620 int icode; 621 REAL_VALUE_TYPE *r1; 622 REAL_VALUE_TYPE *r2; 623{ 624 UEMUSHORT d1[NE], d2[NE], v[NE]; 625 enum tree_code code; 626 627 GET_REAL (r1, d1); 628 GET_REAL (r2, d2); 629#ifdef NANS 630/* Return NaN input back to the caller. */ 631 if (eisnan (d1)) 632 { 633 PUT_REAL (d1, value); 634 return; 635 } 636 if (eisnan (d2)) 637 { 638 PUT_REAL (d2, value); 639 return; 640 } 641#endif 642 code = (enum tree_code) icode; 643 switch (code) 644 { 645 case PLUS_EXPR: 646 eadd (d2, d1, v); 647 break; 648 649 case MINUS_EXPR: 650 esub (d2, d1, v); /* d1 - d2 */ 651 break; 652 653 case MULT_EXPR: 654 emul (d2, d1, v); 655 break; 656 657 case RDIV_EXPR: 658#ifndef REAL_INFINITY 659 if (ecmp (d2, ezero) == 0) 660 { 661#ifdef NANS 662 enan (v, eisneg (d1) ^ eisneg (d2)); 663 break; 664#else 665 abort (); 666#endif 667 } 668#endif 669 ediv (d2, d1, v); /* d1/d2 */ 670 break; 671 672 case MIN_EXPR: /* min (d1,d2) */ 673 if (ecmp (d1, d2) < 0) 674 emov (d1, v); 675 else 676 emov (d2, v); 677 break; 678 679 case MAX_EXPR: /* max (d1,d2) */ 680 if (ecmp (d1, d2) > 0) 681 emov (d1, v); 682 else 683 emov (d2, v); 684 break; 685 default: 686 emov (ezero, v); 687 break; 688 } 689PUT_REAL (v, value); 690} 691 692 693/* Truncate REAL_VALUE_TYPE toward zero to signed HOST_WIDE_INT. 694 implements REAL_VALUE_RNDZINT (x) (etrunci (x)). */ 695 696REAL_VALUE_TYPE 697etrunci (x) 698 REAL_VALUE_TYPE x; 699{ 700 UEMUSHORT f[NE], g[NE]; 701 REAL_VALUE_TYPE r; 702 HOST_WIDE_INT l; 703 704 GET_REAL (&x, g); 705#ifdef NANS 706 if (eisnan (g)) 707 return (x); 708#endif 709 eifrac (g, &l, f); 710 ltoe (&l, g); 711 PUT_REAL (g, &r); 712 return (r); 713} 714 715 716/* Truncate REAL_VALUE_TYPE toward zero to unsigned HOST_WIDE_INT; 717 implements REAL_VALUE_UNSIGNED_RNDZINT (x) (etruncui (x)). */ 718 719REAL_VALUE_TYPE 720etruncui (x) 721 REAL_VALUE_TYPE x; 722{ 723 UEMUSHORT f[NE], g[NE]; 724 REAL_VALUE_TYPE r; 725 unsigned HOST_WIDE_INT l; 726 727 GET_REAL (&x, g); 728#ifdef NANS 729 if (eisnan (g)) 730 return (x); 731#endif 732 euifrac (g, &l, f); 733 ultoe (&l, g); 734 PUT_REAL (g, &r); 735 return (r); 736} 737 738 739/* This is the REAL_VALUE_ATOF function. It converts a decimal or hexadecimal 740 string to binary, rounding off as indicated by the machine_mode argument. 741 Then it promotes the rounded value to REAL_VALUE_TYPE. */ 742 743REAL_VALUE_TYPE 744ereal_atof (s, t) 745 const char *s; 746 enum machine_mode t; 747{ 748 UEMUSHORT tem[NE], e[NE]; 749 REAL_VALUE_TYPE r; 750 751 switch (t) 752 { 753#ifdef C4X 754 case QFmode: 755 case HFmode: 756 asctoe53 (s, tem); 757 e53toe (tem, e); 758 break; 759#else 760 case HFmode: 761#endif 762 763 case SFmode: 764 asctoe24 (s, tem); 765 e24toe (tem, e); 766 break; 767 768 case DFmode: 769 asctoe53 (s, tem); 770 e53toe (tem, e); 771 break; 772 773 case TFmode: 774#if (INTEL_EXTENDED_IEEE_FORMAT == 0) 775 asctoe113 (s, tem); 776 e113toe (tem, e); 777 break; 778#endif 779 /* FALLTHRU */ 780 781 case XFmode: 782 asctoe64 (s, tem); 783 e64toe (tem, e); 784 break; 785 786 default: 787 asctoe (s, e); 788 } 789 PUT_REAL (e, &r); 790 return (r); 791} 792 793 794/* Expansion of REAL_NEGATE. */ 795 796REAL_VALUE_TYPE 797ereal_negate (x) 798 REAL_VALUE_TYPE x; 799{ 800 UEMUSHORT e[NE]; 801 REAL_VALUE_TYPE r; 802 803 GET_REAL (&x, e); 804 eneg (e); 805 PUT_REAL (e, &r); 806 return (r); 807} 808 809 810/* Round real toward zero to HOST_WIDE_INT; 811 implements REAL_VALUE_FIX (x). */ 812 813HOST_WIDE_INT 814efixi (x) 815 REAL_VALUE_TYPE x; 816{ 817 UEMUSHORT f[NE], g[NE]; 818 HOST_WIDE_INT l; 819 820 GET_REAL (&x, f); 821#ifdef NANS 822 if (eisnan (f)) 823 { 824 warning ("conversion from NaN to int"); 825 return (-1); 826 } 827#endif 828 eifrac (f, &l, g); 829 return l; 830} 831 832/* Round real toward zero to unsigned HOST_WIDE_INT 833 implements REAL_VALUE_UNSIGNED_FIX (x). 834 Negative input returns zero. */ 835 836unsigned HOST_WIDE_INT 837efixui (x) 838 REAL_VALUE_TYPE x; 839{ 840 UEMUSHORT f[NE], g[NE]; 841 unsigned HOST_WIDE_INT l; 842 843 GET_REAL (&x, f); 844#ifdef NANS 845 if (eisnan (f)) 846 { 847 warning ("conversion from NaN to unsigned int"); 848 return (-1); 849 } 850#endif 851 euifrac (f, &l, g); 852 return l; 853} 854 855 856/* REAL_VALUE_FROM_INT macro. */ 857 858void 859ereal_from_int (d, i, j, mode) 860 REAL_VALUE_TYPE *d; 861 HOST_WIDE_INT i, j; 862 enum machine_mode mode; 863{ 864 UEMUSHORT df[NE], dg[NE]; 865 HOST_WIDE_INT low, high; 866 int sign; 867 868 if (GET_MODE_CLASS (mode) != MODE_FLOAT) 869 abort (); 870 sign = 0; 871 low = i; 872 if ((high = j) < 0) 873 { 874 sign = 1; 875 /* complement and add 1 */ 876 high = ~high; 877 if (low) 878 low = -low; 879 else 880 high += 1; 881 } 882 eldexp (eone, HOST_BITS_PER_WIDE_INT, df); 883 ultoe ((unsigned HOST_WIDE_INT *) &high, dg); 884 emul (dg, df, dg); 885 ultoe ((unsigned HOST_WIDE_INT *) &low, df); 886 eadd (df, dg, dg); 887 if (sign) 888 eneg (dg); 889 890 /* A REAL_VALUE_TYPE may not be wide enough to hold the two HOST_WIDE_INTS. 891 Avoid double-rounding errors later by rounding off now from the 892 extra-wide internal format to the requested precision. */ 893 switch (GET_MODE_BITSIZE (mode)) 894 { 895 case 32: 896 etoe24 (dg, df); 897 e24toe (df, dg); 898 break; 899 900 case 64: 901 etoe53 (dg, df); 902 e53toe (df, dg); 903 break; 904 905 case 96: 906 etoe64 (dg, df); 907 e64toe (df, dg); 908 break; 909 910 case 128: 911#if (INTEL_EXTENDED_IEEE_FORMAT == 0) 912 etoe113 (dg, df); 913 e113toe (df, dg); 914#else 915 etoe64 (dg, df); 916 e64toe (df, dg); 917#endif 918 break; 919 920 default: 921 abort (); 922 } 923 924 PUT_REAL (dg, d); 925} 926 927 928/* REAL_VALUE_FROM_UNSIGNED_INT macro. */ 929 930void 931ereal_from_uint (d, i, j, mode) 932 REAL_VALUE_TYPE *d; 933 unsigned HOST_WIDE_INT i, j; 934 enum machine_mode mode; 935{ 936 UEMUSHORT df[NE], dg[NE]; 937 unsigned HOST_WIDE_INT low, high; 938 939 if (GET_MODE_CLASS (mode) != MODE_FLOAT) 940 abort (); 941 low = i; 942 high = j; 943 eldexp (eone, HOST_BITS_PER_WIDE_INT, df); 944 ultoe (&high, dg); 945 emul (dg, df, dg); 946 ultoe (&low, df); 947 eadd (df, dg, dg); 948 949 /* A REAL_VALUE_TYPE may not be wide enough to hold the two HOST_WIDE_INTS. 950 Avoid double-rounding errors later by rounding off now from the 951 extra-wide internal format to the requested precision. */ 952 switch (GET_MODE_BITSIZE (mode)) 953 { 954 case 32: 955 etoe24 (dg, df); 956 e24toe (df, dg); 957 break; 958 959 case 64: 960 etoe53 (dg, df); 961 e53toe (df, dg); 962 break; 963 964 case 96: 965 etoe64 (dg, df); 966 e64toe (df, dg); 967 break; 968 969 case 128: 970#if (INTEL_EXTENDED_IEEE_FORMAT == 0) 971 etoe113 (dg, df); 972 e113toe (df, dg); 973#else 974 etoe64 (dg, df); 975 e64toe (df, dg); 976#endif 977 break; 978 979 default: 980 abort (); 981 } 982 983 PUT_REAL (dg, d); 984} 985 986 987/* REAL_VALUE_TO_INT macro. */ 988 989void 990ereal_to_int (low, high, rr) 991 HOST_WIDE_INT *low, *high; 992 REAL_VALUE_TYPE rr; 993{ 994 UEMUSHORT d[NE], df[NE], dg[NE], dh[NE]; 995 int s; 996 997 GET_REAL (&rr, d); 998#ifdef NANS 999 if (eisnan (d)) 1000 { 1001 warning ("conversion from NaN to int"); 1002 *low = -1; 1003 *high = -1; 1004 return; 1005 } 1006#endif 1007 /* convert positive value */ 1008 s = 0; 1009 if (eisneg (d)) 1010 { 1011 eneg (d); 1012 s = 1; 1013 } 1014 eldexp (eone, HOST_BITS_PER_WIDE_INT, df); 1015 ediv (df, d, dg); /* dg = d / 2^32 is the high word */ 1016 euifrac (dg, (unsigned HOST_WIDE_INT *) high, dh); 1017 emul (df, dh, dg); /* fractional part is the low word */ 1018 euifrac (dg, (unsigned HOST_WIDE_INT *) low, dh); 1019 if (s) 1020 { 1021 /* complement and add 1 */ 1022 *high = ~(*high); 1023 if (*low) 1024 *low = -(*low); 1025 else 1026 *high += 1; 1027 } 1028} 1029 1030 1031/* REAL_VALUE_LDEXP macro. */ 1032 1033REAL_VALUE_TYPE 1034ereal_ldexp (x, n) 1035 REAL_VALUE_TYPE x; 1036 int n; 1037{ 1038 UEMUSHORT e[NE], y[NE]; 1039 REAL_VALUE_TYPE r; 1040 1041 GET_REAL (&x, e); 1042#ifdef NANS 1043 if (eisnan (e)) 1044 return (x); 1045#endif 1046 eldexp (e, n, y); 1047 PUT_REAL (y, &r); 1048 return (r); 1049} 1050 1051/* These routines are conditionally compiled because functions 1052 of the same names may be defined in fold-const.c. */ 1053 1054#ifdef REAL_ARITHMETIC 1055 1056/* Check for infinity in a REAL_VALUE_TYPE. */ 1057 1058int 1059target_isinf (x) 1060 REAL_VALUE_TYPE x ATTRIBUTE_UNUSED; 1061{ 1062#ifdef INFINITY 1063 UEMUSHORT e[NE]; 1064 1065 GET_REAL (&x, e); 1066 return (eisinf (e)); 1067#else 1068 return 0; 1069#endif 1070} 1071 1072/* Check whether a REAL_VALUE_TYPE item is a NaN. */ 1073 1074int 1075target_isnan (x) 1076 REAL_VALUE_TYPE x ATTRIBUTE_UNUSED; 1077{ 1078#ifdef NANS 1079 UEMUSHORT e[NE]; 1080 1081 GET_REAL (&x, e); 1082 return (eisnan (e)); 1083#else 1084 return (0); 1085#endif 1086} 1087 1088 1089/* Check for a negative REAL_VALUE_TYPE number. 1090 This just checks the sign bit, so that -0 counts as negative. */ 1091 1092int 1093target_negative (x) 1094 REAL_VALUE_TYPE x; 1095{ 1096 return ereal_isneg (x); 1097} 1098 1099/* Expansion of REAL_VALUE_TRUNCATE. 1100 The result is in floating point, rounded to nearest or even. */ 1101 1102REAL_VALUE_TYPE 1103real_value_truncate (mode, arg) 1104 enum machine_mode mode; 1105 REAL_VALUE_TYPE arg; 1106{ 1107 UEMUSHORT e[NE], t[NE]; 1108 REAL_VALUE_TYPE r; 1109 1110 GET_REAL (&arg, e); 1111#ifdef NANS 1112 if (eisnan (e)) 1113 return (arg); 1114#endif 1115 eclear (t); 1116 switch (mode) 1117 { 1118 case TFmode: 1119#if (INTEL_EXTENDED_IEEE_FORMAT == 0) 1120 etoe113 (e, t); 1121 e113toe (t, t); 1122 break; 1123#endif 1124 /* FALLTHRU */ 1125 1126 case XFmode: 1127 etoe64 (e, t); 1128 e64toe (t, t); 1129 break; 1130 1131 case DFmode: 1132 etoe53 (e, t); 1133 e53toe (t, t); 1134 break; 1135 1136 case SFmode: 1137#ifndef C4X 1138 case HFmode: 1139#endif 1140 etoe24 (e, t); 1141 e24toe (t, t); 1142 break; 1143 1144#ifdef C4X 1145 case HFmode: 1146 case QFmode: 1147 etoe53 (e, t); 1148 e53toe (t, t); 1149 break; 1150#endif 1151 1152 case SImode: 1153 r = etrunci (arg); 1154 return (r); 1155 1156 /* If an unsupported type was requested, presume that 1157 the machine files know something useful to do with 1158 the unmodified value. */ 1159 1160 default: 1161 return (arg); 1162 } 1163 PUT_REAL (t, &r); 1164 return (r); 1165} 1166 1167/* Try to change R into its exact multiplicative inverse in machine mode 1168 MODE. Return nonzero function value if successful. */ 1169 1170int 1171exact_real_inverse (mode, r) 1172 enum machine_mode mode; 1173 REAL_VALUE_TYPE *r; 1174{ 1175 UEMUSHORT e[NE], einv[NE]; 1176 REAL_VALUE_TYPE rinv; 1177 int i; 1178 1179 GET_REAL (r, e); 1180 1181 /* Test for input in range. Don't transform IEEE special values. */ 1182 if (eisinf (e) || eisnan (e) || (ecmp (e, ezero) == 0)) 1183 return 0; 1184 1185 /* Test for a power of 2: all significand bits zero except the MSB. 1186 We are assuming the target has binary (or hex) arithmetic. */ 1187 if (e[NE - 2] != 0x8000) 1188 return 0; 1189 1190 for (i = 0; i < NE - 2; i++) 1191 { 1192 if (e[i] != 0) 1193 return 0; 1194 } 1195 1196 /* Compute the inverse and truncate it to the required mode. */ 1197 ediv (e, eone, einv); 1198 PUT_REAL (einv, &rinv); 1199 rinv = real_value_truncate (mode, rinv); 1200 1201#ifdef CHECK_FLOAT_VALUE 1202 /* This check is not redundant. It may, for example, flush 1203 a supposedly IEEE denormal value to zero. */ 1204 i = 0; 1205 if (CHECK_FLOAT_VALUE (mode, rinv, i)) 1206 return 0; 1207#endif 1208 GET_REAL (&rinv, einv); 1209 1210 /* Check the bits again, because the truncation might have 1211 generated an arbitrary saturation value on overflow. */ 1212 if (einv[NE - 2] != 0x8000) 1213 return 0; 1214 1215 for (i = 0; i < NE - 2; i++) 1216 { 1217 if (einv[i] != 0) 1218 return 0; 1219 } 1220 1221 /* Fail if the computed inverse is out of range. */ 1222 if (eisinf (einv) || eisnan (einv) || (ecmp (einv, ezero) == 0)) 1223 return 0; 1224 1225 /* Output the reciprocal and return success flag. */ 1226 PUT_REAL (einv, r); 1227 return 1; 1228} 1229#endif /* REAL_ARITHMETIC defined */ 1230 1231/* Used for debugging--print the value of R in human-readable format 1232 on stderr. */ 1233 1234void 1235debug_real (r) 1236 REAL_VALUE_TYPE r; 1237{ 1238 char dstr[30]; 1239 1240 REAL_VALUE_TO_DECIMAL (r, "%.20g", dstr); 1241 fprintf (stderr, "%s", dstr); 1242} 1243 1244 1245/* The following routines convert REAL_VALUE_TYPE to the various floating 1246 point formats that are meaningful to supported computers. 1247 1248 The results are returned in 32-bit pieces, each piece stored in a `long'. 1249 This is so they can be printed by statements like 1250 1251 fprintf (file, "%lx, %lx", L[0], L[1]); 1252 1253 that will work on both narrow- and wide-word host computers. */ 1254 1255/* Convert R to a 128-bit long double precision value. The output array L 1256 contains four 32-bit pieces of the result, in the order they would appear 1257 in memory. */ 1258 1259void 1260etartdouble (r, l) 1261 REAL_VALUE_TYPE r; 1262 long l[]; 1263{ 1264 UEMUSHORT e[NE]; 1265 1266 GET_REAL (&r, e); 1267#if INTEL_EXTENDED_IEEE_FORMAT == 0 1268 etoe113 (e, e); 1269#else 1270 etoe64 (e, e); 1271#endif 1272 endian (e, l, TFmode); 1273} 1274 1275/* Convert R to a double extended precision value. The output array L 1276 contains three 32-bit pieces of the result, in the order they would 1277 appear in memory. */ 1278 1279void 1280etarldouble (r, l) 1281 REAL_VALUE_TYPE r; 1282 long l[]; 1283{ 1284 UEMUSHORT e[NE]; 1285 1286 GET_REAL (&r, e); 1287 etoe64 (e, e); 1288 endian (e, l, XFmode); 1289} 1290 1291/* Convert R to a double precision value. The output array L contains two 1292 32-bit pieces of the result, in the order they would appear in memory. */ 1293 1294void 1295etardouble (r, l) 1296 REAL_VALUE_TYPE r; 1297 long l[]; 1298{ 1299 UEMUSHORT e[NE]; 1300 1301 GET_REAL (&r, e); 1302 etoe53 (e, e); 1303 endian (e, l, DFmode); 1304} 1305 1306/* Convert R to a single precision float value stored in the least-significant 1307 bits of a `long'. */ 1308 1309long 1310etarsingle (r) 1311 REAL_VALUE_TYPE r; 1312{ 1313 UEMUSHORT e[NE]; 1314 long l; 1315 1316 GET_REAL (&r, e); 1317 etoe24 (e, e); 1318 endian (e, &l, SFmode); 1319 return ((long) l); 1320} 1321 1322/* Convert X to a decimal ASCII string S for output to an assembly 1323 language file. Note, there is no standard way to spell infinity or 1324 a NaN, so these values may require special treatment in the tm.h 1325 macros. */ 1326 1327void 1328ereal_to_decimal (x, s) 1329 REAL_VALUE_TYPE x; 1330 char *s; 1331{ 1332 UEMUSHORT e[NE]; 1333 1334 GET_REAL (&x, e); 1335 etoasc (e, s, 20); 1336} 1337 1338/* Compare X and Y. Return 1 if X > Y, 0 if X == Y, -1 if X < Y, 1339 or -2 if either is a NaN. */ 1340 1341int 1342ereal_cmp (x, y) 1343 REAL_VALUE_TYPE x, y; 1344{ 1345 UEMUSHORT ex[NE], ey[NE]; 1346 1347 GET_REAL (&x, ex); 1348 GET_REAL (&y, ey); 1349 return (ecmp (ex, ey)); 1350} 1351 1352/* Return 1 if the sign bit of X is set, else return 0. */ 1353 1354int 1355ereal_isneg (x) 1356 REAL_VALUE_TYPE x; 1357{ 1358 UEMUSHORT ex[NE]; 1359 1360 GET_REAL (&x, ex); 1361 return (eisneg (ex)); 1362} 1363 1364/* End of REAL_ARITHMETIC interface */ 1365 1366/* 1367 Extended precision IEEE binary floating point arithmetic routines 1368 1369 Numbers are stored in C language as arrays of 16-bit unsigned 1370 short integers. The arguments of the routines are pointers to 1371 the arrays. 1372 1373 External e type data structure, similar to Intel 8087 chip 1374 temporary real format but possibly with a larger significand: 1375 1376 NE-1 significand words (least significant word first, 1377 most significant bit is normally set) 1378 exponent (value = EXONE for 1.0, 1379 top bit is the sign) 1380 1381 1382 Internal exploded e-type data structure of a number (a "word" is 16 bits): 1383 1384 ei[0] sign word (0 for positive, 0xffff for negative) 1385 ei[1] biased exponent (value = EXONE for the number 1.0) 1386 ei[2] high guard word (always zero after normalization) 1387 ei[3] 1388 to ei[NI-2] significand (NI-4 significand words, 1389 most significant word first, 1390 most significant bit is set) 1391 ei[NI-1] low guard word (0x8000 bit is rounding place) 1392 1393 1394 1395 Routines for external format e-type numbers 1396 1397 asctoe (string, e) ASCII string to extended double e type 1398 asctoe64 (string, &d) ASCII string to long double 1399 asctoe53 (string, &d) ASCII string to double 1400 asctoe24 (string, &f) ASCII string to single 1401 asctoeg (string, e, prec) ASCII string to specified precision 1402 e24toe (&f, e) IEEE single precision to e type 1403 e53toe (&d, e) IEEE double precision to e type 1404 e64toe (&d, e) IEEE long double precision to e type 1405 e113toe (&d, e) 128-bit long double precision to e type 1406#if 0 1407 eabs (e) absolute value 1408#endif 1409 eadd (a, b, c) c = b + a 1410 eclear (e) e = 0 1411 ecmp (a, b) Returns 1 if a > b, 0 if a == b, 1412 -1 if a < b, -2 if either a or b is a NaN. 1413 ediv (a, b, c) c = b / a 1414 efloor (a, b) truncate to integer, toward -infinity 1415 efrexp (a, exp, s) extract exponent and significand 1416 eifrac (e, &l, frac) e to HOST_WIDE_INT and e type fraction 1417 euifrac (e, &l, frac) e to unsigned HOST_WIDE_INT and e type fraction 1418 einfin (e) set e to infinity, leaving its sign alone 1419 eldexp (a, n, b) multiply by 2**n 1420 emov (a, b) b = a 1421 emul (a, b, c) c = b * a 1422 eneg (e) e = -e 1423#if 0 1424 eround (a, b) b = nearest integer value to a 1425#endif 1426 esub (a, b, c) c = b - a 1427#if 0 1428 e24toasc (&f, str, n) single to ASCII string, n digits after decimal 1429 e53toasc (&d, str, n) double to ASCII string, n digits after decimal 1430 e64toasc (&d, str, n) 80-bit long double to ASCII string 1431 e113toasc (&d, str, n) 128-bit long double to ASCII string 1432#endif 1433 etoasc (e, str, n) e to ASCII string, n digits after decimal 1434 etoe24 (e, &f) convert e type to IEEE single precision 1435 etoe53 (e, &d) convert e type to IEEE double precision 1436 etoe64 (e, &d) convert e type to IEEE long double precision 1437 ltoe (&l, e) HOST_WIDE_INT to e type 1438 ultoe (&l, e) unsigned HOST_WIDE_INT to e type 1439 eisneg (e) 1 if sign bit of e != 0, else 0 1440 eisinf (e) 1 if e has maximum exponent (non-IEEE) 1441 or is infinite (IEEE) 1442 eisnan (e) 1 if e is a NaN 1443 1444 1445 Routines for internal format exploded e-type numbers 1446 1447 eaddm (ai, bi) add significands, bi = bi + ai 1448 ecleaz (ei) ei = 0 1449 ecleazs (ei) set ei = 0 but leave its sign alone 1450 ecmpm (ai, bi) compare significands, return 1, 0, or -1 1451 edivm (ai, bi) divide significands, bi = bi / ai 1452 emdnorm (ai,l,s,exp) normalize and round off 1453 emovi (a, ai) convert external a to internal ai 1454 emovo (ai, a) convert internal ai to external a 1455 emovz (ai, bi) bi = ai, low guard word of bi = 0 1456 emulm (ai, bi) multiply significands, bi = bi * ai 1457 enormlz (ei) left-justify the significand 1458 eshdn1 (ai) shift significand and guards down 1 bit 1459 eshdn8 (ai) shift down 8 bits 1460 eshdn6 (ai) shift down 16 bits 1461 eshift (ai, n) shift ai n bits up (or down if n < 0) 1462 eshup1 (ai) shift significand and guards up 1 bit 1463 eshup8 (ai) shift up 8 bits 1464 eshup6 (ai) shift up 16 bits 1465 esubm (ai, bi) subtract significands, bi = bi - ai 1466 eiisinf (ai) 1 if infinite 1467 eiisnan (ai) 1 if a NaN 1468 eiisneg (ai) 1 if sign bit of ai != 0, else 0 1469 einan (ai) set ai = NaN 1470#if 0 1471 eiinfin (ai) set ai = infinity 1472#endif 1473 1474 The result is always normalized and rounded to NI-4 word precision 1475 after each arithmetic operation. 1476 1477 Exception flags are NOT fully supported. 1478 1479 Signaling NaN's are NOT supported; they are treated the same 1480 as quiet NaN's. 1481 1482 Define INFINITY for support of infinity; otherwise a 1483 saturation arithmetic is implemented. 1484 1485 Define NANS for support of Not-a-Number items; otherwise the 1486 arithmetic will never produce a NaN output, and might be confused 1487 by a NaN input. 1488 If NaN's are supported, the output of `ecmp (a,b)' is -2 if 1489 either a or b is a NaN. This means asking `if (ecmp (a,b) < 0)' 1490 may not be legitimate. Use `if (ecmp (a,b) == -1)' for `less than' 1491 if in doubt. 1492 1493 Denormals are always supported here where appropriate (e.g., not 1494 for conversion to DEC numbers). */ 1495 1496/* Definitions for error codes that are passed to the common error handling 1497 routine mtherr. 1498 1499 For Digital Equipment PDP-11 and VAX computers, certain 1500 IBM systems, and others that use numbers with a 56-bit 1501 significand, the symbol DEC should be defined. In this 1502 mode, most floating point constants are given as arrays 1503 of octal integers to eliminate decimal to binary conversion 1504 errors that might be introduced by the compiler. 1505 1506 For computers, such as IBM PC, that follow the IEEE 1507 Standard for Binary Floating Point Arithmetic (ANSI/IEEE 1508 Std 754-1985), the symbol IEEE should be defined. 1509 These numbers have 53-bit significands. In this mode, constants 1510 are provided as arrays of hexadecimal 16 bit integers. 1511 The endian-ness of generated values is controlled by 1512 REAL_WORDS_BIG_ENDIAN. 1513 1514 To accommodate other types of computer arithmetic, all 1515 constants are also provided in a normal decimal radix 1516 which one can hope are correctly converted to a suitable 1517 format by the available C language compiler. To invoke 1518 this mode, the symbol UNK is defined. 1519 1520 An important difference among these modes is a predefined 1521 set of machine arithmetic constants for each. The numbers 1522 MACHEP (the machine roundoff error), MAXNUM (largest number 1523 represented), and several other parameters are preset by 1524 the configuration symbol. Check the file const.c to 1525 ensure that these values are correct for your computer. 1526 1527 For ANSI C compatibility, define ANSIC equal to 1. Currently 1528 this affects only the atan2 function and others that use it. */ 1529 1530/* Constant definitions for math error conditions. */ 1531 1532#define DOMAIN 1 /* argument domain error */ 1533#define SING 2 /* argument singularity */ 1534#define OVERFLOW 3 /* overflow range error */ 1535#define UNDERFLOW 4 /* underflow range error */ 1536#define TLOSS 5 /* total loss of precision */ 1537#define PLOSS 6 /* partial loss of precision */ 1538#define INVALID 7 /* NaN-producing operation */ 1539 1540/* e type constants used by high precision check routines */ 1541 1542#if MAX_LONG_DOUBLE_TYPE_SIZE == 128 && (INTEL_EXTENDED_IEEE_FORMAT == 0) 1543/* 0.0 */ 1544const UEMUSHORT ezero[NE] = 1545 {0x0000, 0x0000, 0x0000, 0x0000, 1546 0x0000, 0x0000, 0x0000, 0x0000, 0x0000, 0x0000,}; 1547 1548/* 5.0E-1 */ 1549const UEMUSHORT ehalf[NE] = 1550 {0x0000, 0x0000, 0x0000, 0x0000, 1551 0x0000, 0x0000, 0x0000, 0x0000, 0x8000, 0x3ffe,}; 1552 1553/* 1.0E0 */ 1554const UEMUSHORT eone[NE] = 1555 {0x0000, 0x0000, 0x0000, 0x0000, 1556 0x0000, 0x0000, 0x0000, 0x0000, 0x8000, 0x3fff,}; 1557 1558/* 2.0E0 */ 1559const UEMUSHORT etwo[NE] = 1560 {0x0000, 0x0000, 0x0000, 0x0000, 1561 0x0000, 0x0000, 0x0000, 0x0000, 0x8000, 0x4000,}; 1562 1563/* 3.2E1 */ 1564const UEMUSHORT e32[NE] = 1565 {0x0000, 0x0000, 0x0000, 0x0000, 1566 0x0000, 0x0000, 0x0000, 0x0000, 0x8000, 0x4004,}; 1567 1568/* 6.93147180559945309417232121458176568075500134360255E-1 */ 1569const UEMUSHORT elog2[NE] = 1570 {0x40f3, 0xf6af, 0x03f2, 0xb398, 1571 0xc9e3, 0x79ab, 0150717, 0013767, 0130562, 0x3ffe,}; 1572 1573/* 1.41421356237309504880168872420969807856967187537695E0 */ 1574const UEMUSHORT esqrt2[NE] = 1575 {0x1d6f, 0xbe9f, 0x754a, 0x89b3, 1576 0x597d, 0x6484, 0174736, 0171463, 0132404, 0x3fff,}; 1577 1578/* 3.14159265358979323846264338327950288419716939937511E0 */ 1579const UEMUSHORT epi[NE] = 1580 {0x2902, 0x1cd1, 0x80dc, 0x628b, 1581 0xc4c6, 0xc234, 0020550, 0155242, 0144417, 0040000,}; 1582 1583#else 1584/* LONG_DOUBLE_TYPE_SIZE is other than 128 */ 1585const UEMUSHORT ezero[NE] = 1586 {0, 0000000, 0000000, 0000000, 0000000, 0000000,}; 1587const UEMUSHORT ehalf[NE] = 1588 {0, 0000000, 0000000, 0000000, 0100000, 0x3ffe,}; 1589const UEMUSHORT eone[NE] = 1590 {0, 0000000, 0000000, 0000000, 0100000, 0x3fff,}; 1591const UEMUSHORT etwo[NE] = 1592 {0, 0000000, 0000000, 0000000, 0100000, 0040000,}; 1593const UEMUSHORT e32[NE] = 1594 {0, 0000000, 0000000, 0000000, 0100000, 0040004,}; 1595const UEMUSHORT elog2[NE] = 1596 {0xc9e4, 0x79ab, 0150717, 0013767, 0130562, 0x3ffe,}; 1597const UEMUSHORT esqrt2[NE] = 1598 {0x597e, 0x6484, 0174736, 0171463, 0132404, 0x3fff,}; 1599const UEMUSHORT epi[NE] = 1600 {0xc4c6, 0xc234, 0020550, 0155242, 0144417, 0040000,}; 1601#endif 1602 1603/* Control register for rounding precision. 1604 This can be set to 113 (if NE=10), 80 (if NE=6), 64, 56, 53, or 24 bits. */ 1605 1606int rndprc = NBITS; 1607extern int rndprc; 1608 1609/* Clear out entire e-type number X. */ 1610 1611static void 1612eclear (x) 1613 UEMUSHORT *x; 1614{ 1615 int i; 1616 1617 for (i = 0; i < NE; i++) 1618 *x++ = 0; 1619} 1620 1621/* Move e-type number from A to B. */ 1622 1623static void 1624emov (a, b) 1625 const UEMUSHORT *a; 1626 UEMUSHORT *b; 1627{ 1628 int i; 1629 1630 for (i = 0; i < NE; i++) 1631 *b++ = *a++; 1632} 1633 1634 1635#if 0 1636/* Absolute value of e-type X. */ 1637 1638static void 1639eabs (x) 1640 UEMUSHORT x[]; 1641{ 1642 /* sign is top bit of last word of external format */ 1643 x[NE - 1] &= 0x7fff; 1644} 1645#endif /* 0 */ 1646 1647/* Negate the e-type number X. */ 1648 1649static void 1650eneg (x) 1651 UEMUSHORT x[]; 1652{ 1653 1654 x[NE - 1] ^= 0x8000; /* Toggle the sign bit */ 1655} 1656 1657/* Return 1 if sign bit of e-type number X is nonzero, else zero. */ 1658 1659static int 1660eisneg (x) 1661 const UEMUSHORT x[]; 1662{ 1663 1664 if (x[NE - 1] & 0x8000) 1665 return (1); 1666 else 1667 return (0); 1668} 1669 1670/* Return 1 if e-type number X is infinity, else return zero. */ 1671 1672static int 1673eisinf (x) 1674 const UEMUSHORT x[]; 1675{ 1676 1677#ifdef NANS 1678 if (eisnan (x)) 1679 return (0); 1680#endif 1681 if ((x[NE - 1] & 0x7fff) == 0x7fff) 1682 return (1); 1683 else 1684 return (0); 1685} 1686 1687/* Check if e-type number is not a number. The bit pattern is one that we 1688 defined, so we know for sure how to detect it. */ 1689 1690static int 1691eisnan (x) 1692 const UEMUSHORT x[] ATTRIBUTE_UNUSED; 1693{ 1694#ifdef NANS 1695 int i; 1696 1697 /* NaN has maximum exponent */ 1698 if ((x[NE - 1] & 0x7fff) != 0x7fff) 1699 return (0); 1700 /* ... and non-zero significand field. */ 1701 for (i = 0; i < NE - 1; i++) 1702 { 1703 if (*x++ != 0) 1704 return (1); 1705 } 1706#endif 1707 1708 return (0); 1709} 1710 1711/* Fill e-type number X with infinity pattern (IEEE) 1712 or largest possible number (non-IEEE). */ 1713 1714static void 1715einfin (x) 1716 UEMUSHORT *x; 1717{ 1718 int i; 1719 1720#ifdef INFINITY 1721 for (i = 0; i < NE - 1; i++) 1722 *x++ = 0; 1723 *x |= 32767; 1724#else 1725 for (i = 0; i < NE - 1; i++) 1726 *x++ = 0xffff; 1727 *x |= 32766; 1728 if (rndprc < NBITS) 1729 { 1730 if (rndprc == 113) 1731 { 1732 *(x - 9) = 0; 1733 *(x - 8) = 0; 1734 } 1735 if (rndprc == 64) 1736 { 1737 *(x - 5) = 0; 1738 } 1739 if (rndprc == 53) 1740 { 1741 *(x - 4) = 0xf800; 1742 } 1743 else 1744 { 1745 *(x - 4) = 0; 1746 *(x - 3) = 0; 1747 *(x - 2) = 0xff00; 1748 } 1749 } 1750#endif 1751} 1752 1753/* Output an e-type NaN. 1754 This generates Intel's quiet NaN pattern for extended real. 1755 The exponent is 7fff, the leading mantissa word is c000. */ 1756 1757#ifdef NANS 1758static void 1759enan (x, sign) 1760 UEMUSHORT *x; 1761 int sign; 1762{ 1763 int i; 1764 1765 for (i = 0; i < NE - 2; i++) 1766 *x++ = 0; 1767 *x++ = 0xc000; 1768 *x = (sign << 15) | 0x7fff; 1769} 1770#endif /* NANS */ 1771 1772/* Move in an e-type number A, converting it to exploded e-type B. */ 1773 1774static void 1775emovi (a, b) 1776 const UEMUSHORT *a; 1777 UEMUSHORT *b; 1778{ 1779 const UEMUSHORT *p; 1780 UEMUSHORT *q; 1781 int i; 1782 1783 q = b; 1784 p = a + (NE - 1); /* point to last word of external number */ 1785 /* get the sign bit */ 1786 if (*p & 0x8000) 1787 *q++ = 0xffff; 1788 else 1789 *q++ = 0; 1790 /* get the exponent */ 1791 *q = *p--; 1792 *q++ &= 0x7fff; /* delete the sign bit */ 1793#ifdef INFINITY 1794 if ((*(q - 1) & 0x7fff) == 0x7fff) 1795 { 1796#ifdef NANS 1797 if (eisnan (a)) 1798 { 1799 *q++ = 0; 1800 for (i = 3; i < NI; i++) 1801 *q++ = *p--; 1802 return; 1803 } 1804#endif 1805 1806 for (i = 2; i < NI; i++) 1807 *q++ = 0; 1808 return; 1809 } 1810#endif 1811 1812 /* clear high guard word */ 1813 *q++ = 0; 1814 /* move in the significand */ 1815 for (i = 0; i < NE - 1; i++) 1816 *q++ = *p--; 1817 /* clear low guard word */ 1818 *q = 0; 1819} 1820 1821/* Move out exploded e-type number A, converting it to e type B. */ 1822 1823static void 1824emovo (a, b) 1825 const UEMUSHORT *a; 1826 UEMUSHORT *b; 1827{ 1828 const UEMUSHORT *p; 1829 UEMUSHORT *q; 1830 UEMUSHORT i; 1831 int j; 1832 1833 p = a; 1834 q = b + (NE - 1); /* point to output exponent */ 1835 /* combine sign and exponent */ 1836 i = *p++; 1837 if (i) 1838 *q-- = *p++ | 0x8000; 1839 else 1840 *q-- = *p++; 1841#ifdef INFINITY 1842 if (*(p - 1) == 0x7fff) 1843 { 1844#ifdef NANS 1845 if (eiisnan (a)) 1846 { 1847 enan (b, eiisneg (a)); 1848 return; 1849 } 1850#endif 1851 einfin (b); 1852 return; 1853 } 1854#endif 1855 /* skip over guard word */ 1856 ++p; 1857 /* move the significand */ 1858 for (j = 0; j < NE - 1; j++) 1859 *q-- = *p++; 1860} 1861 1862/* Clear out exploded e-type number XI. */ 1863 1864static void 1865ecleaz (xi) 1866 UEMUSHORT *xi; 1867{ 1868 int i; 1869 1870 for (i = 0; i < NI; i++) 1871 *xi++ = 0; 1872} 1873 1874/* Clear out exploded e-type XI, but don't touch the sign. */ 1875 1876static void 1877ecleazs (xi) 1878 UEMUSHORT *xi; 1879{ 1880 int i; 1881 1882 ++xi; 1883 for (i = 0; i < NI - 1; i++) 1884 *xi++ = 0; 1885} 1886 1887/* Move exploded e-type number from A to B. */ 1888 1889static void 1890emovz (a, b) 1891 const UEMUSHORT *a; 1892 UEMUSHORT *b; 1893{ 1894 int i; 1895 1896 for (i = 0; i < NI - 1; i++) 1897 *b++ = *a++; 1898 /* clear low guard word */ 1899 *b = 0; 1900} 1901 1902/* Generate exploded e-type NaN. 1903 The explicit pattern for this is maximum exponent and 1904 top two significant bits set. */ 1905 1906#ifdef NANS 1907static void 1908einan (x) 1909 UEMUSHORT x[]; 1910{ 1911 1912 ecleaz (x); 1913 x[E] = 0x7fff; 1914 x[M + 1] = 0xc000; 1915} 1916#endif /* NANS */ 1917 1918/* Return nonzero if exploded e-type X is a NaN. */ 1919 1920#ifdef NANS 1921static int 1922eiisnan (x) 1923 const UEMUSHORT x[]; 1924{ 1925 int i; 1926 1927 if ((x[E] & 0x7fff) == 0x7fff) 1928 { 1929 for (i = M + 1; i < NI; i++) 1930 { 1931 if (x[i] != 0) 1932 return (1); 1933 } 1934 } 1935 return (0); 1936} 1937#endif /* NANS */ 1938 1939/* Return nonzero if sign of exploded e-type X is nonzero. */ 1940 1941#ifdef NANS 1942static int 1943eiisneg (x) 1944 const UEMUSHORT x[]; 1945{ 1946 1947 return x[0] != 0; 1948} 1949#endif /* NANS */ 1950 1951#if 0 1952/* Fill exploded e-type X with infinity pattern. 1953 This has maximum exponent and significand all zeros. */ 1954 1955static void 1956eiinfin (x) 1957 UEMUSHORT x[]; 1958{ 1959 1960 ecleaz (x); 1961 x[E] = 0x7fff; 1962} 1963#endif /* 0 */ 1964 1965/* Return nonzero if exploded e-type X is infinite. */ 1966 1967#ifdef INFINITY 1968static int 1969eiisinf (x) 1970 const UEMUSHORT x[]; 1971{ 1972 1973#ifdef NANS 1974 if (eiisnan (x)) 1975 return (0); 1976#endif 1977 if ((x[E] & 0x7fff) == 0x7fff) 1978 return (1); 1979 return (0); 1980} 1981#endif /* INFINITY */ 1982 1983/* Compare significands of numbers in internal exploded e-type format. 1984 Guard words are included in the comparison. 1985 1986 Returns +1 if a > b 1987 0 if a == b 1988 -1 if a < b */ 1989 1990static int 1991ecmpm (a, b) 1992 const UEMUSHORT *a, *b; 1993{ 1994 int i; 1995 1996 a += M; /* skip up to significand area */ 1997 b += M; 1998 for (i = M; i < NI; i++) 1999 { 2000 if (*a++ != *b++) 2001 goto difrnt; 2002 } 2003 return (0); 2004 2005 difrnt: 2006 if (*(--a) > *(--b)) 2007 return (1); 2008 else 2009 return (-1); 2010} 2011 2012/* Shift significand of exploded e-type X down by 1 bit. */ 2013 2014static void 2015eshdn1 (x) 2016 UEMUSHORT *x; 2017{ 2018 UEMUSHORT bits; 2019 int i; 2020 2021 x += M; /* point to significand area */ 2022 2023 bits = 0; 2024 for (i = M; i < NI; i++) 2025 { 2026 if (*x & 1) 2027 bits |= 1; 2028 *x >>= 1; 2029 if (bits & 2) 2030 *x |= 0x8000; 2031 bits <<= 1; 2032 ++x; 2033 } 2034} 2035 2036/* Shift significand of exploded e-type X up by 1 bit. */ 2037 2038static void 2039eshup1 (x) 2040 UEMUSHORT *x; 2041{ 2042 UEMUSHORT bits; 2043 int i; 2044 2045 x += NI - 1; 2046 bits = 0; 2047 2048 for (i = M; i < NI; i++) 2049 { 2050 if (*x & 0x8000) 2051 bits |= 1; 2052 *x <<= 1; 2053 if (bits & 2) 2054 *x |= 1; 2055 bits <<= 1; 2056 --x; 2057 } 2058} 2059 2060 2061/* Shift significand of exploded e-type X down by 8 bits. */ 2062 2063static void 2064eshdn8 (x) 2065 UEMUSHORT *x; 2066{ 2067 UEMUSHORT newbyt, oldbyt; 2068 int i; 2069 2070 x += M; 2071 oldbyt = 0; 2072 for (i = M; i < NI; i++) 2073 { 2074 newbyt = *x << 8; 2075 *x >>= 8; 2076 *x |= oldbyt; 2077 oldbyt = newbyt; 2078 ++x; 2079 } 2080} 2081 2082/* Shift significand of exploded e-type X up by 8 bits. */ 2083 2084static void 2085eshup8 (x) 2086 UEMUSHORT *x; 2087{ 2088 int i; 2089 UEMUSHORT newbyt, oldbyt; 2090 2091 x += NI - 1; 2092 oldbyt = 0; 2093 2094 for (i = M; i < NI; i++) 2095 { 2096 newbyt = *x >> 8; 2097 *x <<= 8; 2098 *x |= oldbyt; 2099 oldbyt = newbyt; 2100 --x; 2101 } 2102} 2103 2104/* Shift significand of exploded e-type X up by 16 bits. */ 2105 2106static void 2107eshup6 (x) 2108 UEMUSHORT *x; 2109{ 2110 int i; 2111 UEMUSHORT *p; 2112 2113 p = x + M; 2114 x += M + 1; 2115 2116 for (i = M; i < NI - 1; i++) 2117 *p++ = *x++; 2118 2119 *p = 0; 2120} 2121 2122/* Shift significand of exploded e-type X down by 16 bits. */ 2123 2124static void 2125eshdn6 (x) 2126 UEMUSHORT *x; 2127{ 2128 int i; 2129 UEMUSHORT *p; 2130 2131 x += NI - 1; 2132 p = x + 1; 2133 2134 for (i = M; i < NI - 1; i++) 2135 *(--p) = *(--x); 2136 2137 *(--p) = 0; 2138} 2139 2140/* Add significands of exploded e-type X and Y. X + Y replaces Y. */ 2141 2142static void 2143eaddm (x, y) 2144 const UEMUSHORT *x; 2145 UEMUSHORT *y; 2146{ 2147 unsigned EMULONG a; 2148 int i; 2149 unsigned int carry; 2150 2151 x += NI - 1; 2152 y += NI - 1; 2153 carry = 0; 2154 for (i = M; i < NI; i++) 2155 { 2156 a = (unsigned EMULONG) (*x) + (unsigned EMULONG) (*y) + carry; 2157 if (a & 0x10000) 2158 carry = 1; 2159 else 2160 carry = 0; 2161 *y = (UEMUSHORT) a; 2162 --x; 2163 --y; 2164 } 2165} 2166 2167/* Subtract significands of exploded e-type X and Y. Y - X replaces Y. */ 2168 2169static void 2170esubm (x, y) 2171 const UEMUSHORT *x; 2172 UEMUSHORT *y; 2173{ 2174 unsigned EMULONG a; 2175 int i; 2176 unsigned int carry; 2177 2178 x += NI - 1; 2179 y += NI - 1; 2180 carry = 0; 2181 for (i = M; i < NI; i++) 2182 { 2183 a = (unsigned EMULONG) (*y) - (unsigned EMULONG) (*x) - carry; 2184 if (a & 0x10000) 2185 carry = 1; 2186 else 2187 carry = 0; 2188 *y = (UEMUSHORT) a; 2189 --x; 2190 --y; 2191 } 2192} 2193 2194 2195static UEMUSHORT equot[NI]; 2196 2197 2198#if 0 2199/* Radix 2 shift-and-add versions of multiply and divide */ 2200 2201 2202/* Divide significands */ 2203 2204int 2205edivm (den, num) 2206 UEMUSHORT den[], num[]; 2207{ 2208 int i; 2209 UEMUSHORT *p, *q; 2210 UEMUSHORT j; 2211 2212 p = &equot[0]; 2213 *p++ = num[0]; 2214 *p++ = num[1]; 2215 2216 for (i = M; i < NI; i++) 2217 { 2218 *p++ = 0; 2219 } 2220 2221 /* Use faster compare and subtraction if denominator has only 15 bits of 2222 significance. */ 2223 2224 p = &den[M + 2]; 2225 if (*p++ == 0) 2226 { 2227 for (i = M + 3; i < NI; i++) 2228 { 2229 if (*p++ != 0) 2230 goto fulldiv; 2231 } 2232 if ((den[M + 1] & 1) != 0) 2233 goto fulldiv; 2234 eshdn1 (num); 2235 eshdn1 (den); 2236 2237 p = &den[M + 1]; 2238 q = &num[M + 1]; 2239 2240 for (i = 0; i < NBITS + 2; i++) 2241 { 2242 if (*p <= *q) 2243 { 2244 *q -= *p; 2245 j = 1; 2246 } 2247 else 2248 { 2249 j = 0; 2250 } 2251 eshup1 (equot); 2252 equot[NI - 2] |= j; 2253 eshup1 (num); 2254 } 2255 goto divdon; 2256 } 2257 2258 /* The number of quotient bits to calculate is NBITS + 1 scaling guard 2259 bit + 1 roundoff bit. */ 2260 2261 fulldiv: 2262 2263 p = &equot[NI - 2]; 2264 for (i = 0; i < NBITS + 2; i++) 2265 { 2266 if (ecmpm (den, num) <= 0) 2267 { 2268 esubm (den, num); 2269 j = 1; /* quotient bit = 1 */ 2270 } 2271 else 2272 j = 0; 2273 eshup1 (equot); 2274 *p |= j; 2275 eshup1 (num); 2276 } 2277 2278 divdon: 2279 2280 eshdn1 (equot); 2281 eshdn1 (equot); 2282 2283 /* test for nonzero remainder after roundoff bit */ 2284 p = &num[M]; 2285 j = 0; 2286 for (i = M; i < NI; i++) 2287 { 2288 j |= *p++; 2289 } 2290 if (j) 2291 j = 1; 2292 2293 2294 for (i = 0; i < NI; i++) 2295 num[i] = equot[i]; 2296 return ((int) j); 2297} 2298 2299 2300/* Multiply significands */ 2301 2302int 2303emulm (a, b) 2304 UEMUSHORT a[], b[]; 2305{ 2306 UEMUSHORT *p, *q; 2307 int i, j, k; 2308 2309 equot[0] = b[0]; 2310 equot[1] = b[1]; 2311 for (i = M; i < NI; i++) 2312 equot[i] = 0; 2313 2314 p = &a[NI - 2]; 2315 k = NBITS; 2316 while (*p == 0) /* significand is not supposed to be zero */ 2317 { 2318 eshdn6 (a); 2319 k -= 16; 2320 } 2321 if ((*p & 0xff) == 0) 2322 { 2323 eshdn8 (a); 2324 k -= 8; 2325 } 2326 2327 q = &equot[NI - 1]; 2328 j = 0; 2329 for (i = 0; i < k; i++) 2330 { 2331 if (*p & 1) 2332 eaddm (b, equot); 2333 /* remember if there were any nonzero bits shifted out */ 2334 if (*q & 1) 2335 j |= 1; 2336 eshdn1 (a); 2337 eshdn1 (equot); 2338 } 2339 2340 for (i = 0; i < NI; i++) 2341 b[i] = equot[i]; 2342 2343 /* return flag for lost nonzero bits */ 2344 return (j); 2345} 2346 2347#else 2348 2349/* Radix 65536 versions of multiply and divide. */ 2350 2351/* Multiply significand of e-type number B 2352 by 16-bit quantity A, return e-type result to C. */ 2353 2354static void 2355m16m (a, b, c) 2356 unsigned int a; 2357 const UEMUSHORT b[]; 2358 UEMUSHORT c[]; 2359{ 2360 UEMUSHORT *pp; 2361 unsigned EMULONG carry; 2362 const UEMUSHORT *ps; 2363 UEMUSHORT p[NI]; 2364 unsigned EMULONG aa, m; 2365 int i; 2366 2367 aa = a; 2368 pp = &p[NI-2]; 2369 *pp++ = 0; 2370 *pp = 0; 2371 ps = &b[NI-1]; 2372 2373 for (i=M+1; i<NI; i++) 2374 { 2375 if (*ps == 0) 2376 { 2377 --ps; 2378 --pp; 2379 *(pp-1) = 0; 2380 } 2381 else 2382 { 2383 m = (unsigned EMULONG) aa * *ps--; 2384 carry = (m & 0xffff) + *pp; 2385 *pp-- = (UEMUSHORT) carry; 2386 carry = (carry >> 16) + (m >> 16) + *pp; 2387 *pp = (UEMUSHORT) carry; 2388 *(pp-1) = carry >> 16; 2389 } 2390 } 2391 for (i=M; i<NI; i++) 2392 c[i] = p[i]; 2393} 2394 2395/* Divide significands of exploded e-types NUM / DEN. Neither the 2396 numerator NUM nor the denominator DEN is permitted to have its high guard 2397 word nonzero. */ 2398 2399static int 2400edivm (den, num) 2401 const UEMUSHORT den[]; 2402 UEMUSHORT num[]; 2403{ 2404 int i; 2405 UEMUSHORT *p; 2406 unsigned EMULONG tnum; 2407 UEMUSHORT j, tdenm, tquot; 2408 UEMUSHORT tprod[NI+1]; 2409 2410 p = &equot[0]; 2411 *p++ = num[0]; 2412 *p++ = num[1]; 2413 2414 for (i=M; i<NI; i++) 2415 { 2416 *p++ = 0; 2417 } 2418 eshdn1 (num); 2419 tdenm = den[M+1]; 2420 for (i=M; i<NI; i++) 2421 { 2422 /* Find trial quotient digit (the radix is 65536). */ 2423 tnum = (((unsigned EMULONG) num[M]) << 16) + num[M+1]; 2424 2425 /* Do not execute the divide instruction if it will overflow. */ 2426 if ((tdenm * (unsigned long) 0xffff) < tnum) 2427 tquot = 0xffff; 2428 else 2429 tquot = tnum / tdenm; 2430 /* Multiply denominator by trial quotient digit. */ 2431 m16m ((unsigned int) tquot, den, tprod); 2432 /* The quotient digit may have been overestimated. */ 2433 if (ecmpm (tprod, num) > 0) 2434 { 2435 tquot -= 1; 2436 esubm (den, tprod); 2437 if (ecmpm (tprod, num) > 0) 2438 { 2439 tquot -= 1; 2440 esubm (den, tprod); 2441 } 2442 } 2443 esubm (tprod, num); 2444 equot[i] = tquot; 2445 eshup6 (num); 2446 } 2447 /* test for nonzero remainder after roundoff bit */ 2448 p = &num[M]; 2449 j = 0; 2450 for (i=M; i<NI; i++) 2451 { 2452 j |= *p++; 2453 } 2454 if (j) 2455 j = 1; 2456 2457 for (i=0; i<NI; i++) 2458 num[i] = equot[i]; 2459 2460 return ((int) j); 2461} 2462 2463/* Multiply significands of exploded e-type A and B, result in B. */ 2464 2465static int 2466emulm (a, b) 2467 const UEMUSHORT a[]; 2468 UEMUSHORT b[]; 2469{ 2470 const UEMUSHORT *p; 2471 UEMUSHORT *q; 2472 UEMUSHORT pprod[NI]; 2473 UEMUSHORT j; 2474 int i; 2475 2476 equot[0] = b[0]; 2477 equot[1] = b[1]; 2478 for (i=M; i<NI; i++) 2479 equot[i] = 0; 2480 2481 j = 0; 2482 p = &a[NI-1]; 2483 q = &equot[NI-1]; 2484 for (i=M+1; i<NI; i++) 2485 { 2486 if (*p == 0) 2487 { 2488 --p; 2489 } 2490 else 2491 { 2492 m16m ((unsigned int) *p--, b, pprod); 2493 eaddm (pprod, equot); 2494 } 2495 j |= *q; 2496 eshdn6 (equot); 2497 } 2498 2499 for (i=0; i<NI; i++) 2500 b[i] = equot[i]; 2501 2502 /* return flag for lost nonzero bits */ 2503 return ((int) j); 2504} 2505#endif 2506 2507 2508/* Normalize and round off. 2509 2510 The internal format number to be rounded is S. 2511 Input LOST is 0 if the value is exact. This is the so-called sticky bit. 2512 2513 Input SUBFLG indicates whether the number was obtained 2514 by a subtraction operation. In that case if LOST is nonzero 2515 then the number is slightly smaller than indicated. 2516 2517 Input EXP is the biased exponent, which may be negative. 2518 the exponent field of S is ignored but is replaced by 2519 EXP as adjusted by normalization and rounding. 2520 2521 Input RCNTRL is the rounding control. If it is nonzero, the 2522 returned value will be rounded to RNDPRC bits. 2523 2524 For future reference: In order for emdnorm to round off denormal 2525 significands at the right point, the input exponent must be 2526 adjusted to be the actual value it would have after conversion to 2527 the final floating point type. This adjustment has been 2528 implemented for all type conversions (etoe53, etc.) and decimal 2529 conversions, but not for the arithmetic functions (eadd, etc.). 2530 Data types having standard 15-bit exponents are not affected by 2531 this, but SFmode and DFmode are affected. For example, ediv with 2532 rndprc = 24 will not round correctly to 24-bit precision if the 2533 result is denormal. */ 2534 2535static int rlast = -1; 2536static int rw = 0; 2537static UEMUSHORT rmsk = 0; 2538static UEMUSHORT rmbit = 0; 2539static UEMUSHORT rebit = 0; 2540static int re = 0; 2541static UEMUSHORT rbit[NI]; 2542 2543static void 2544emdnorm (s, lost, subflg, exp, rcntrl) 2545 UEMUSHORT s[]; 2546 int lost; 2547 int subflg; 2548 EMULONG exp; 2549 int rcntrl; 2550{ 2551 int i, j; 2552 UEMUSHORT r; 2553 2554 /* Normalize */ 2555 j = enormlz (s); 2556 2557 /* a blank significand could mean either zero or infinity. */ 2558#ifndef INFINITY 2559 if (j > NBITS) 2560 { 2561 ecleazs (s); 2562 return; 2563 } 2564#endif 2565 exp -= j; 2566#ifndef INFINITY 2567 if (exp >= 32767L) 2568 goto overf; 2569#else 2570 if ((j > NBITS) && (exp < 32767)) 2571 { 2572 ecleazs (s); 2573 return; 2574 } 2575#endif 2576 if (exp < 0L) 2577 { 2578 if (exp > (EMULONG) (-NBITS - 1)) 2579 { 2580 j = (int) exp; 2581 i = eshift (s, j); 2582 if (i) 2583 lost = 1; 2584 } 2585 else 2586 { 2587 ecleazs (s); 2588 return; 2589 } 2590 } 2591 /* Round off, unless told not to by rcntrl. */ 2592 if (rcntrl == 0) 2593 goto mdfin; 2594 /* Set up rounding parameters if the control register changed. */ 2595 if (rndprc != rlast) 2596 { 2597 ecleaz (rbit); 2598 switch (rndprc) 2599 { 2600 default: 2601 case NBITS: 2602 rw = NI - 1; /* low guard word */ 2603 rmsk = 0xffff; 2604 rmbit = 0x8000; 2605 re = rw - 1; 2606 rebit = 1; 2607 break; 2608 2609 case 113: 2610 rw = 10; 2611 rmsk = 0x7fff; 2612 rmbit = 0x4000; 2613 rebit = 0x8000; 2614 re = rw; 2615 break; 2616 2617 case 64: 2618 rw = 7; 2619 rmsk = 0xffff; 2620 rmbit = 0x8000; 2621 re = rw - 1; 2622 rebit = 1; 2623 break; 2624 2625 /* For DEC or IBM arithmetic */ 2626 case 56: 2627 rw = 6; 2628 rmsk = 0xff; 2629 rmbit = 0x80; 2630 rebit = 0x100; 2631 re = rw; 2632 break; 2633 2634 case 53: 2635 rw = 6; 2636 rmsk = 0x7ff; 2637 rmbit = 0x0400; 2638 rebit = 0x800; 2639 re = rw; 2640 break; 2641 2642 /* For C4x arithmetic */ 2643 case 32: 2644 rw = 5; 2645 rmsk = 0xffff; 2646 rmbit = 0x8000; 2647 rebit = 1; 2648 re = rw - 1; 2649 break; 2650 2651 case 24: 2652 rw = 4; 2653 rmsk = 0xff; 2654 rmbit = 0x80; 2655 rebit = 0x100; 2656 re = rw; 2657 break; 2658 } 2659 rbit[re] = rebit; 2660 rlast = rndprc; 2661 } 2662 2663 /* Shift down 1 temporarily if the data structure has an implied 2664 most significant bit and the number is denormal. 2665 Intel long double denormals also lose one bit of precision. */ 2666 if ((exp <= 0) && (rndprc != NBITS) 2667 && ((rndprc != 64) || ((rndprc == 64) && ! REAL_WORDS_BIG_ENDIAN))) 2668 { 2669 lost |= s[NI - 1] & 1; 2670 eshdn1 (s); 2671 } 2672 /* Clear out all bits below the rounding bit, 2673 remembering in r if any were nonzero. */ 2674 r = s[rw] & rmsk; 2675 if (rndprc < NBITS) 2676 { 2677 i = rw + 1; 2678 while (i < NI) 2679 { 2680 if (s[i]) 2681 r |= 1; 2682 s[i] = 0; 2683 ++i; 2684 } 2685 } 2686 s[rw] &= ~rmsk; 2687 if ((r & rmbit) != 0) 2688 { 2689#ifndef C4X 2690 if (r == rmbit) 2691 { 2692 if (lost == 0) 2693 { /* round to even */ 2694 if ((s[re] & rebit) == 0) 2695 goto mddone; 2696 } 2697 else 2698 { 2699 if (subflg != 0) 2700 goto mddone; 2701 } 2702 } 2703#endif 2704 eaddm (rbit, s); 2705 } 2706 mddone: 2707/* Undo the temporary shift for denormal values. */ 2708 if ((exp <= 0) && (rndprc != NBITS) 2709 && ((rndprc != 64) || ((rndprc == 64) && ! REAL_WORDS_BIG_ENDIAN))) 2710 { 2711 eshup1 (s); 2712 } 2713 if (s[2] != 0) 2714 { /* overflow on roundoff */ 2715 eshdn1 (s); 2716 exp += 1; 2717 } 2718 mdfin: 2719 s[NI - 1] = 0; 2720 if (exp >= 32767L) 2721 { 2722#ifndef INFINITY 2723 overf: 2724#endif 2725#ifdef INFINITY 2726 s[1] = 32767; 2727 for (i = 2; i < NI - 1; i++) 2728 s[i] = 0; 2729 if (extra_warnings) 2730 warning ("floating point overflow"); 2731#else 2732 s[1] = 32766; 2733 s[2] = 0; 2734 for (i = M + 1; i < NI - 1; i++) 2735 s[i] = 0xffff; 2736 s[NI - 1] = 0; 2737 if ((rndprc < 64) || (rndprc == 113)) 2738 { 2739 s[rw] &= ~rmsk; 2740 if (rndprc == 24) 2741 { 2742 s[5] = 0; 2743 s[6] = 0; 2744 } 2745 } 2746#endif 2747 return; 2748 } 2749 if (exp < 0) 2750 s[1] = 0; 2751 else 2752 s[1] = (UEMUSHORT) exp; 2753} 2754 2755/* Subtract. C = B - A, all e type numbers. */ 2756 2757static int subflg = 0; 2758 2759static void 2760esub (a, b, c) 2761 const UEMUSHORT *a, *b; 2762 UEMUSHORT *c; 2763{ 2764 2765#ifdef NANS 2766 if (eisnan (a)) 2767 { 2768 emov (a, c); 2769 return; 2770 } 2771 if (eisnan (b)) 2772 { 2773 emov (b, c); 2774 return; 2775 } 2776/* Infinity minus infinity is a NaN. 2777 Test for subtracting infinities of the same sign. */ 2778 if (eisinf (a) && eisinf (b) 2779 && ((eisneg (a) ^ eisneg (b)) == 0)) 2780 { 2781 mtherr ("esub", INVALID); 2782 enan (c, 0); 2783 return; 2784 } 2785#endif 2786 subflg = 1; 2787 eadd1 (a, b, c); 2788} 2789 2790/* Add. C = A + B, all e type. */ 2791 2792static void 2793eadd (a, b, c) 2794 const UEMUSHORT *a, *b; 2795 UEMUSHORT *c; 2796{ 2797 2798#ifdef NANS 2799/* NaN plus anything is a NaN. */ 2800 if (eisnan (a)) 2801 { 2802 emov (a, c); 2803 return; 2804 } 2805 if (eisnan (b)) 2806 { 2807 emov (b, c); 2808 return; 2809 } 2810/* Infinity minus infinity is a NaN. 2811 Test for adding infinities of opposite signs. */ 2812 if (eisinf (a) && eisinf (b) 2813 && ((eisneg (a) ^ eisneg (b)) != 0)) 2814 { 2815 mtherr ("esub", INVALID); 2816 enan (c, 0); 2817 return; 2818 } 2819#endif 2820 subflg = 0; 2821 eadd1 (a, b, c); 2822} 2823 2824/* Arithmetic common to both addition and subtraction. */ 2825 2826static void 2827eadd1 (a, b, c) 2828 const UEMUSHORT *a, *b; 2829 UEMUSHORT *c; 2830{ 2831 UEMUSHORT ai[NI], bi[NI], ci[NI]; 2832 int i, lost, j, k; 2833 EMULONG lt, lta, ltb; 2834 2835#ifdef INFINITY 2836 if (eisinf (a)) 2837 { 2838 emov (a, c); 2839 if (subflg) 2840 eneg (c); 2841 return; 2842 } 2843 if (eisinf (b)) 2844 { 2845 emov (b, c); 2846 return; 2847 } 2848#endif 2849 emovi (a, ai); 2850 emovi (b, bi); 2851 if (subflg) 2852 ai[0] = ~ai[0]; 2853 2854 /* compare exponents */ 2855 lta = ai[E]; 2856 ltb = bi[E]; 2857 lt = lta - ltb; 2858 if (lt > 0L) 2859 { /* put the larger number in bi */ 2860 emovz (bi, ci); 2861 emovz (ai, bi); 2862 emovz (ci, ai); 2863 ltb = bi[E]; 2864 lt = -lt; 2865 } 2866 lost = 0; 2867 if (lt != 0L) 2868 { 2869 if (lt < (EMULONG) (-NBITS - 1)) 2870 goto done; /* answer same as larger addend */ 2871 k = (int) lt; 2872 lost = eshift (ai, k); /* shift the smaller number down */ 2873 } 2874 else 2875 { 2876 /* exponents were the same, so must compare significands */ 2877 i = ecmpm (ai, bi); 2878 if (i == 0) 2879 { /* the numbers are identical in magnitude */ 2880 /* if different signs, result is zero */ 2881 if (ai[0] != bi[0]) 2882 { 2883 eclear (c); 2884 return; 2885 } 2886 /* if same sign, result is double */ 2887 /* double denormalized tiny number */ 2888 if ((bi[E] == 0) && ((bi[3] & 0x8000) == 0)) 2889 { 2890 eshup1 (bi); 2891 goto done; 2892 } 2893 /* add 1 to exponent unless both are zero! */ 2894 for (j = 1; j < NI - 1; j++) 2895 { 2896 if (bi[j] != 0) 2897 { 2898 ltb += 1; 2899 if (ltb >= 0x7fff) 2900 { 2901 eclear (c); 2902 if (ai[0] != 0) 2903 eneg (c); 2904 einfin (c); 2905 return; 2906 } 2907 break; 2908 } 2909 } 2910 bi[E] = (UEMUSHORT) ltb; 2911 goto done; 2912 } 2913 if (i > 0) 2914 { /* put the larger number in bi */ 2915 emovz (bi, ci); 2916 emovz (ai, bi); 2917 emovz (ci, ai); 2918 } 2919 } 2920 if (ai[0] == bi[0]) 2921 { 2922 eaddm (ai, bi); 2923 subflg = 0; 2924 } 2925 else 2926 { 2927 esubm (ai, bi); 2928 subflg = 1; 2929 } 2930 emdnorm (bi, lost, subflg, ltb, 64); 2931 2932 done: 2933 emovo (bi, c); 2934} 2935 2936/* Divide: C = B/A, all e type. */ 2937 2938static void 2939ediv (a, b, c) 2940 const UEMUSHORT *a, *b; 2941 UEMUSHORT *c; 2942{ 2943 UEMUSHORT ai[NI], bi[NI]; 2944 int i, sign; 2945 EMULONG lt, lta, ltb; 2946 2947/* IEEE says if result is not a NaN, the sign is "-" if and only if 2948 operands have opposite signs -- but flush -0 to 0 later if not IEEE. */ 2949 sign = eisneg (a) ^ eisneg (b); 2950 2951#ifdef NANS 2952/* Return any NaN input. */ 2953 if (eisnan (a)) 2954 { 2955 emov (a, c); 2956 return; 2957 } 2958 if (eisnan (b)) 2959 { 2960 emov (b, c); 2961 return; 2962 } 2963/* Zero over zero, or infinity over infinity, is a NaN. */ 2964 if (((ecmp (a, ezero) == 0) && (ecmp (b, ezero) == 0)) 2965 || (eisinf (a) && eisinf (b))) 2966 { 2967 mtherr ("ediv", INVALID); 2968 enan (c, sign); 2969 return; 2970 } 2971#endif 2972/* Infinity over anything else is infinity. */ 2973#ifdef INFINITY 2974 if (eisinf (b)) 2975 { 2976 einfin (c); 2977 goto divsign; 2978 } 2979/* Anything else over infinity is zero. */ 2980 if (eisinf (a)) 2981 { 2982 eclear (c); 2983 goto divsign; 2984 } 2985#endif 2986 emovi (a, ai); 2987 emovi (b, bi); 2988 lta = ai[E]; 2989 ltb = bi[E]; 2990 if (bi[E] == 0) 2991 { /* See if numerator is zero. */ 2992 for (i = 1; i < NI - 1; i++) 2993 { 2994 if (bi[i] != 0) 2995 { 2996 ltb -= enormlz (bi); 2997 goto dnzro1; 2998 } 2999 } 3000 eclear (c); 3001 goto divsign; 3002 } 3003 dnzro1: 3004 3005 if (ai[E] == 0) 3006 { /* possible divide by zero */ 3007 for (i = 1; i < NI - 1; i++) 3008 { 3009 if (ai[i] != 0) 3010 { 3011 lta -= enormlz (ai); 3012 goto dnzro2; 3013 } 3014 } 3015/* Divide by zero is not an invalid operation. 3016 It is a divide-by-zero operation! */ 3017 einfin (c); 3018 mtherr ("ediv", SING); 3019 goto divsign; 3020 } 3021 dnzro2: 3022 3023 i = edivm (ai, bi); 3024 /* calculate exponent */ 3025 lt = ltb - lta + EXONE; 3026 emdnorm (bi, i, 0, lt, 64); 3027 emovo (bi, c); 3028 3029 divsign: 3030 3031 if (sign 3032#ifndef IEEE 3033 && (ecmp (c, ezero) != 0) 3034#endif 3035 ) 3036 *(c+(NE-1)) |= 0x8000; 3037 else 3038 *(c+(NE-1)) &= ~0x8000; 3039} 3040 3041/* Multiply e-types A and B, return e-type product C. */ 3042 3043static void 3044emul (a, b, c) 3045 const UEMUSHORT *a, *b; 3046 UEMUSHORT *c; 3047{ 3048 UEMUSHORT ai[NI], bi[NI]; 3049 int i, j, sign; 3050 EMULONG lt, lta, ltb; 3051 3052/* IEEE says if result is not a NaN, the sign is "-" if and only if 3053 operands have opposite signs -- but flush -0 to 0 later if not IEEE. */ 3054 sign = eisneg (a) ^ eisneg (b); 3055 3056#ifdef NANS 3057/* NaN times anything is the same NaN. */ 3058 if (eisnan (a)) 3059 { 3060 emov (a, c); 3061 return; 3062 } 3063 if (eisnan (b)) 3064 { 3065 emov (b, c); 3066 return; 3067 } 3068/* Zero times infinity is a NaN. */ 3069 if ((eisinf (a) && (ecmp (b, ezero) == 0)) 3070 || (eisinf (b) && (ecmp (a, ezero) == 0))) 3071 { 3072 mtherr ("emul", INVALID); 3073 enan (c, sign); 3074 return; 3075 } 3076#endif 3077/* Infinity times anything else is infinity. */ 3078#ifdef INFINITY 3079 if (eisinf (a) || eisinf (b)) 3080 { 3081 einfin (c); 3082 goto mulsign; 3083 } 3084#endif 3085 emovi (a, ai); 3086 emovi (b, bi); 3087 lta = ai[E]; 3088 ltb = bi[E]; 3089 if (ai[E] == 0) 3090 { 3091 for (i = 1; i < NI - 1; i++) 3092 { 3093 if (ai[i] != 0) 3094 { 3095 lta -= enormlz (ai); 3096 goto mnzer1; 3097 } 3098 } 3099 eclear (c); 3100 goto mulsign; 3101 } 3102 mnzer1: 3103 3104 if (bi[E] == 0) 3105 { 3106 for (i = 1; i < NI - 1; i++) 3107 { 3108 if (bi[i] != 0) 3109 { 3110 ltb -= enormlz (bi); 3111 goto mnzer2; 3112 } 3113 } 3114 eclear (c); 3115 goto mulsign; 3116 } 3117 mnzer2: 3118 3119 /* Multiply significands */ 3120 j = emulm (ai, bi); 3121 /* calculate exponent */ 3122 lt = lta + ltb - (EXONE - 1); 3123 emdnorm (bi, j, 0, lt, 64); 3124 emovo (bi, c); 3125 3126 mulsign: 3127 3128 if (sign 3129#ifndef IEEE 3130 && (ecmp (c, ezero) != 0) 3131#endif 3132 ) 3133 *(c+(NE-1)) |= 0x8000; 3134 else 3135 *(c+(NE-1)) &= ~0x8000; 3136} 3137 3138/* Convert double precision PE to e-type Y. */ 3139 3140static void 3141e53toe (pe, y) 3142 const UEMUSHORT *pe; 3143 UEMUSHORT *y; 3144{ 3145#ifdef DEC 3146 3147 dectoe (pe, y); 3148 3149#else 3150#ifdef IBM 3151 3152 ibmtoe (pe, y, DFmode); 3153 3154#else 3155#ifdef C4X 3156 3157 c4xtoe (pe, y, HFmode); 3158 3159#else 3160 UEMUSHORT r; 3161 const UEMUSHORT *e; 3162 UEMUSHORT *p; 3163 UEMUSHORT yy[NI]; 3164 int denorm, k; 3165 3166 e = pe; 3167 denorm = 0; /* flag if denormalized number */ 3168 ecleaz (yy); 3169 if (! REAL_WORDS_BIG_ENDIAN) 3170 e += 3; 3171 r = *e; 3172 yy[0] = 0; 3173 if (r & 0x8000) 3174 yy[0] = 0xffff; 3175 yy[M] = (r & 0x0f) | 0x10; 3176 r &= ~0x800f; /* strip sign and 4 significand bits */ 3177#ifdef INFINITY 3178 if (r == 0x7ff0) 3179 { 3180#ifdef NANS 3181 if (! REAL_WORDS_BIG_ENDIAN) 3182 { 3183 if (((pe[3] & 0xf) != 0) || (pe[2] != 0) 3184 || (pe[1] != 0) || (pe[0] != 0)) 3185 { 3186 enan (y, yy[0] != 0); 3187 return; 3188 } 3189 } 3190 else 3191 { 3192 if (((pe[0] & 0xf) != 0) || (pe[1] != 0) 3193 || (pe[2] != 0) || (pe[3] != 0)) 3194 { 3195 enan (y, yy[0] != 0); 3196 return; 3197 } 3198 } 3199#endif /* NANS */ 3200 eclear (y); 3201 einfin (y); 3202 if (yy[0]) 3203 eneg (y); 3204 return; 3205 } 3206#endif /* INFINITY */ 3207 r >>= 4; 3208 /* If zero exponent, then the significand is denormalized. 3209 So take back the understood high significand bit. */ 3210 3211 if (r == 0) 3212 { 3213 denorm = 1; 3214 yy[M] &= ~0x10; 3215 } 3216 r += EXONE - 01777; 3217 yy[E] = r; 3218 p = &yy[M + 1]; 3219#ifdef IEEE 3220 if (! REAL_WORDS_BIG_ENDIAN) 3221 { 3222 *p++ = *(--e); 3223 *p++ = *(--e); 3224 *p++ = *(--e); 3225 } 3226 else 3227 { 3228 ++e; 3229 *p++ = *e++; 3230 *p++ = *e++; 3231 *p++ = *e++; 3232 } 3233#endif 3234 eshift (yy, -5); 3235 if (denorm) 3236 { 3237 /* If zero exponent, then normalize the significand. */ 3238 if ((k = enormlz (yy)) > NBITS) 3239 ecleazs (yy); 3240 else 3241 yy[E] -= (UEMUSHORT) (k - 1); 3242 } 3243 emovo (yy, y); 3244#endif /* not C4X */ 3245#endif /* not IBM */ 3246#endif /* not DEC */ 3247} 3248 3249/* Convert double extended precision float PE to e type Y. */ 3250 3251static void 3252e64toe (pe, y) 3253 const UEMUSHORT *pe; 3254 UEMUSHORT *y; 3255{ 3256 UEMUSHORT yy[NI]; 3257 const UEMUSHORT *e; 3258 UEMUSHORT *p, *q; 3259 int i; 3260 3261 e = pe; 3262 p = yy; 3263 for (i = 0; i < NE - 5; i++) 3264 *p++ = 0; 3265/* This precision is not ordinarily supported on DEC or IBM. */ 3266#ifdef DEC 3267 for (i = 0; i < 5; i++) 3268 *p++ = *e++; 3269#endif 3270#ifdef IBM 3271 p = &yy[0] + (NE - 1); 3272 *p-- = *e++; 3273 ++e; 3274 for (i = 0; i < 5; i++) 3275 *p-- = *e++; 3276#endif 3277#ifdef IEEE 3278 if (! REAL_WORDS_BIG_ENDIAN) 3279 { 3280 for (i = 0; i < 5; i++) 3281 *p++ = *e++; 3282 3283 /* For denormal long double Intel format, shift significand up one 3284 -- but only if the top significand bit is zero. A top bit of 1 3285 is "pseudodenormal" when the exponent is zero. */ 3286 if ((yy[NE-1] & 0x7fff) == 0 && (yy[NE-2] & 0x8000) == 0) 3287 { 3288 UEMUSHORT temp[NI]; 3289 3290 emovi (yy, temp); 3291 eshup1 (temp); 3292 emovo (temp,y); 3293 return; 3294 } 3295 } 3296 else 3297 { 3298 p = &yy[0] + (NE - 1); 3299#ifdef ARM_EXTENDED_IEEE_FORMAT 3300 /* For ARMs, the exponent is in the lowest 15 bits of the word. */ 3301 *p-- = (e[0] & 0x8000) | (e[1] & 0x7ffff); 3302 e += 2; 3303#else 3304 *p-- = *e++; 3305 ++e; 3306#endif 3307 for (i = 0; i < 4; i++) 3308 *p-- = *e++; 3309 } 3310#endif 3311#ifdef INFINITY 3312 /* Point to the exponent field and check max exponent cases. */ 3313 p = &yy[NE - 1]; 3314 if ((*p & 0x7fff) == 0x7fff) 3315 { 3316#ifdef NANS 3317 if (! REAL_WORDS_BIG_ENDIAN) 3318 { 3319 for (i = 0; i < 4; i++) 3320 { 3321 if ((i != 3 && pe[i] != 0) 3322 /* Anything but 0x8000 here, including 0, is a NaN. */ 3323 || (i == 3 && pe[i] != 0x8000)) 3324 { 3325 enan (y, (*p & 0x8000) != 0); 3326 return; 3327 } 3328 } 3329 } 3330 else 3331 { 3332#ifdef ARM_EXTENDED_IEEE_FORMAT 3333 for (i = 2; i <= 5; i++) 3334 { 3335 if (pe[i] != 0) 3336 { 3337 enan (y, (*p & 0x8000) != 0); 3338 return; 3339 } 3340 } 3341#else /* not ARM */ 3342 /* In Motorola extended precision format, the most significant 3343 bit of an infinity mantissa could be either 1 or 0. It is 3344 the lower order bits that tell whether the value is a NaN. */ 3345 if ((pe[2] & 0x7fff) != 0) 3346 goto bigend_nan; 3347 3348 for (i = 3; i <= 5; i++) 3349 { 3350 if (pe[i] != 0) 3351 { 3352bigend_nan: 3353 enan (y, (*p & 0x8000) != 0); 3354 return; 3355 } 3356 } 3357#endif /* not ARM */ 3358 } 3359#endif /* NANS */ 3360 eclear (y); 3361 einfin (y); 3362 if (*p & 0x8000) 3363 eneg (y); 3364 return; 3365 } 3366#endif /* INFINITY */ 3367 p = yy; 3368 q = y; 3369 for (i = 0; i < NE; i++) 3370 *q++ = *p++; 3371} 3372 3373#if (INTEL_EXTENDED_IEEE_FORMAT == 0) 3374/* Convert 128-bit long double precision float PE to e type Y. */ 3375 3376static void 3377e113toe (pe, y) 3378 const UEMUSHORT *pe; 3379 UEMUSHORT *y; 3380{ 3381 UEMUSHORT r; 3382 const UEMUSHORT *e; 3383 UEMUSHORT *p; 3384 UEMUSHORT yy[NI]; 3385 int denorm, i; 3386 3387 e = pe; 3388 denorm = 0; 3389 ecleaz (yy); 3390#ifdef IEEE 3391 if (! REAL_WORDS_BIG_ENDIAN) 3392 e += 7; 3393#endif 3394 r = *e; 3395 yy[0] = 0; 3396 if (r & 0x8000) 3397 yy[0] = 0xffff; 3398 r &= 0x7fff; 3399#ifdef INFINITY 3400 if (r == 0x7fff) 3401 { 3402#ifdef NANS 3403 if (! REAL_WORDS_BIG_ENDIAN) 3404 { 3405 for (i = 0; i < 7; i++) 3406 { 3407 if (pe[i] != 0) 3408 { 3409 enan (y, yy[0] != 0); 3410 return; 3411 } 3412 } 3413 } 3414 else 3415 { 3416 for (i = 1; i < 8; i++) 3417 { 3418 if (pe[i] != 0) 3419 { 3420 enan (y, yy[0] != 0); 3421 return; 3422 } 3423 } 3424 } 3425#endif /* NANS */ 3426 eclear (y); 3427 einfin (y); 3428 if (yy[0]) 3429 eneg (y); 3430 return; 3431 } 3432#endif /* INFINITY */ 3433 yy[E] = r; 3434 p = &yy[M + 1]; 3435#ifdef IEEE 3436 if (! REAL_WORDS_BIG_ENDIAN) 3437 { 3438 for (i = 0; i < 7; i++) 3439 *p++ = *(--e); 3440 } 3441 else 3442 { 3443 ++e; 3444 for (i = 0; i < 7; i++) 3445 *p++ = *e++; 3446 } 3447#endif 3448/* If denormal, remove the implied bit; else shift down 1. */ 3449 if (r == 0) 3450 { 3451 yy[M] = 0; 3452 } 3453 else 3454 { 3455 yy[M] = 1; 3456 eshift (yy, -1); 3457 } 3458 emovo (yy, y); 3459} 3460#endif 3461 3462/* Convert single precision float PE to e type Y. */ 3463 3464static void 3465e24toe (pe, y) 3466 const UEMUSHORT *pe; 3467 UEMUSHORT *y; 3468{ 3469#ifdef IBM 3470 3471 ibmtoe (pe, y, SFmode); 3472 3473#else 3474 3475#ifdef C4X 3476 3477 c4xtoe (pe, y, QFmode); 3478 3479#else 3480 3481 UEMUSHORT r; 3482 const UEMUSHORT *e; 3483 UEMUSHORT *p; 3484 UEMUSHORT yy[NI]; 3485 int denorm, k; 3486 3487 e = pe; 3488 denorm = 0; /* flag if denormalized number */ 3489 ecleaz (yy); 3490#ifdef IEEE 3491 if (! REAL_WORDS_BIG_ENDIAN) 3492 e += 1; 3493#endif 3494#ifdef DEC 3495 e += 1; 3496#endif 3497 r = *e; 3498 yy[0] = 0; 3499 if (r & 0x8000) 3500 yy[0] = 0xffff; 3501 yy[M] = (r & 0x7f) | 0200; 3502 r &= ~0x807f; /* strip sign and 7 significand bits */ 3503#ifdef INFINITY 3504 if (r == 0x7f80) 3505 { 3506#ifdef NANS 3507 if (REAL_WORDS_BIG_ENDIAN) 3508 { 3509 if (((pe[0] & 0x7f) != 0) || (pe[1] != 0)) 3510 { 3511 enan (y, yy[0] != 0); 3512 return; 3513 } 3514 } 3515 else 3516 { 3517 if (((pe[1] & 0x7f) != 0) || (pe[0] != 0)) 3518 { 3519 enan (y, yy[0] != 0); 3520 return; 3521 } 3522 } 3523#endif /* NANS */ 3524 eclear (y); 3525 einfin (y); 3526 if (yy[0]) 3527 eneg (y); 3528 return; 3529 } 3530#endif /* INFINITY */ 3531 r >>= 7; 3532 /* If zero exponent, then the significand is denormalized. 3533 So take back the understood high significand bit. */ 3534 if (r == 0) 3535 { 3536 denorm = 1; 3537 yy[M] &= ~0200; 3538 } 3539 r += EXONE - 0177; 3540 yy[E] = r; 3541 p = &yy[M + 1]; 3542#ifdef DEC 3543 *p++ = *(--e); 3544#endif 3545#ifdef IEEE 3546 if (! REAL_WORDS_BIG_ENDIAN) 3547 *p++ = *(--e); 3548 else 3549 { 3550 ++e; 3551 *p++ = *e++; 3552 } 3553#endif 3554 eshift (yy, -8); 3555 if (denorm) 3556 { /* if zero exponent, then normalize the significand */ 3557 if ((k = enormlz (yy)) > NBITS) 3558 ecleazs (yy); 3559 else 3560 yy[E] -= (UEMUSHORT) (k - 1); 3561 } 3562 emovo (yy, y); 3563#endif /* not C4X */ 3564#endif /* not IBM */ 3565} 3566 3567#if (INTEL_EXTENDED_IEEE_FORMAT == 0) 3568/* Convert e-type X to IEEE 128-bit long double format E. */ 3569 3570static void 3571etoe113 (x, e) 3572 const UEMUSHORT *x; 3573 UEMUSHORT *e; 3574{ 3575 UEMUSHORT xi[NI]; 3576 EMULONG exp; 3577 int rndsav; 3578 3579#ifdef NANS 3580 if (eisnan (x)) 3581 { 3582 make_nan (e, eisneg (x), TFmode); 3583 return; 3584 } 3585#endif 3586 emovi (x, xi); 3587 exp = (EMULONG) xi[E]; 3588#ifdef INFINITY 3589 if (eisinf (x)) 3590 goto nonorm; 3591#endif 3592 /* round off to nearest or even */ 3593 rndsav = rndprc; 3594 rndprc = 113; 3595 emdnorm (xi, 0, 0, exp, 64); 3596 rndprc = rndsav; 3597#ifdef INFINITY 3598 nonorm: 3599#endif 3600 toe113 (xi, e); 3601} 3602 3603/* Convert exploded e-type X, that has already been rounded to 3604 113-bit precision, to IEEE 128-bit long double format Y. */ 3605 3606static void 3607toe113 (a, b) 3608 UEMUSHORT *a, *b; 3609{ 3610 UEMUSHORT *p, *q; 3611 UEMUSHORT i; 3612 3613#ifdef NANS 3614 if (eiisnan (a)) 3615 { 3616 make_nan (b, eiisneg (a), TFmode); 3617 return; 3618 } 3619#endif 3620 p = a; 3621 if (REAL_WORDS_BIG_ENDIAN) 3622 q = b; 3623 else 3624 q = b + 7; /* point to output exponent */ 3625 3626 /* If not denormal, delete the implied bit. */ 3627 if (a[E] != 0) 3628 { 3629 eshup1 (a); 3630 } 3631 /* combine sign and exponent */ 3632 i = *p++; 3633 if (REAL_WORDS_BIG_ENDIAN) 3634 { 3635 if (i) 3636 *q++ = *p++ | 0x8000; 3637 else 3638 *q++ = *p++; 3639 } 3640 else 3641 { 3642 if (i) 3643 *q-- = *p++ | 0x8000; 3644 else 3645 *q-- = *p++; 3646 } 3647 /* skip over guard word */ 3648 ++p; 3649 /* move the significand */ 3650 if (REAL_WORDS_BIG_ENDIAN) 3651 { 3652 for (i = 0; i < 7; i++) 3653 *q++ = *p++; 3654 } 3655 else 3656 { 3657 for (i = 0; i < 7; i++) 3658 *q-- = *p++; 3659 } 3660} 3661#endif 3662 3663/* Convert e-type X to IEEE double extended format E. */ 3664 3665static void 3666etoe64 (x, e) 3667 const UEMUSHORT *x; 3668 UEMUSHORT *e; 3669{ 3670 UEMUSHORT xi[NI]; 3671 EMULONG exp; 3672 int rndsav; 3673 3674#ifdef NANS 3675 if (eisnan (x)) 3676 { 3677 make_nan (e, eisneg (x), XFmode); 3678 return; 3679 } 3680#endif 3681 emovi (x, xi); 3682 /* adjust exponent for offset */ 3683 exp = (EMULONG) xi[E]; 3684#ifdef INFINITY 3685 if (eisinf (x)) 3686 goto nonorm; 3687#endif 3688 /* round off to nearest or even */ 3689 rndsav = rndprc; 3690 rndprc = 64; 3691 emdnorm (xi, 0, 0, exp, 64); 3692 rndprc = rndsav; 3693#ifdef INFINITY 3694 nonorm: 3695#endif 3696 toe64 (xi, e); 3697} 3698 3699/* Convert exploded e-type X, that has already been rounded to 3700 64-bit precision, to IEEE double extended format Y. */ 3701 3702static void 3703toe64 (a, b) 3704 UEMUSHORT *a, *b; 3705{ 3706 UEMUSHORT *p, *q; 3707 UEMUSHORT i; 3708 3709#ifdef NANS 3710 if (eiisnan (a)) 3711 { 3712 make_nan (b, eiisneg (a), XFmode); 3713 return; 3714 } 3715#endif 3716 /* Shift denormal long double Intel format significand down one bit. */ 3717 if ((a[E] == 0) && ! REAL_WORDS_BIG_ENDIAN) 3718 eshdn1 (a); 3719 p = a; 3720#ifdef IBM 3721 q = b; 3722#endif 3723#ifdef DEC 3724 q = b + 4; 3725#endif 3726#ifdef IEEE 3727 if (REAL_WORDS_BIG_ENDIAN) 3728 q = b; 3729 else 3730 { 3731 q = b + 4; /* point to output exponent */ 3732 /* Clear the last two bytes of 12-byte Intel format. q is pointing 3733 into an array of size 6 (e.g. x[NE]), so the last two bytes are 3734 always there, and there are never more bytes, even when we are using 3735 INTEL_EXTENDED_IEEE_FORMAT. */ 3736 *(q+1) = 0; 3737 } 3738#endif 3739 3740 /* combine sign and exponent */ 3741 i = *p++; 3742#ifdef IBM 3743 if (i) 3744 *q++ = *p++ | 0x8000; 3745 else 3746 *q++ = *p++; 3747 *q++ = 0; 3748#endif 3749#ifdef DEC 3750 if (i) 3751 *q-- = *p++ | 0x8000; 3752 else 3753 *q-- = *p++; 3754#endif 3755#ifdef IEEE 3756 if (REAL_WORDS_BIG_ENDIAN) 3757 { 3758#ifdef ARM_EXTENDED_IEEE_FORMAT 3759 /* The exponent is in the lowest 15 bits of the first word. */ 3760 *q++ = i ? 0x8000 : 0; 3761 *q++ = *p++; 3762#else 3763 if (i) 3764 *q++ = *p++ | 0x8000; 3765 else 3766 *q++ = *p++; 3767 *q++ = 0; 3768#endif 3769 } 3770 else 3771 { 3772 if (i) 3773 *q-- = *p++ | 0x8000; 3774 else 3775 *q-- = *p++; 3776 } 3777#endif 3778 /* skip over guard word */ 3779 ++p; 3780 /* move the significand */ 3781#ifdef IBM 3782 for (i = 0; i < 4; i++) 3783 *q++ = *p++; 3784#endif 3785#ifdef DEC 3786 for (i = 0; i < 4; i++) 3787 *q-- = *p++; 3788#endif 3789#ifdef IEEE 3790 if (REAL_WORDS_BIG_ENDIAN) 3791 { 3792 for (i = 0; i < 4; i++) 3793 *q++ = *p++; 3794 } 3795 else 3796 { 3797#ifdef INFINITY 3798 if (eiisinf (a)) 3799 { 3800 /* Intel long double infinity significand. */ 3801 *q-- = 0x8000; 3802 *q-- = 0; 3803 *q-- = 0; 3804 *q = 0; 3805 return; 3806 } 3807#endif 3808 for (i = 0; i < 4; i++) 3809 *q-- = *p++; 3810 } 3811#endif 3812} 3813 3814/* e type to double precision. */ 3815 3816#ifdef DEC 3817/* Convert e-type X to DEC-format double E. */ 3818 3819static void 3820etoe53 (x, e) 3821 const UEMUSHORT *x; 3822 UEMUSHORT *e; 3823{ 3824 etodec (x, e); /* see etodec.c */ 3825} 3826 3827/* Convert exploded e-type X, that has already been rounded to 3828 56-bit double precision, to DEC double Y. */ 3829 3830static void 3831toe53 (x, y) 3832 UEMUSHORT *x, *y; 3833{ 3834 todec (x, y); 3835} 3836 3837#else 3838#ifdef IBM 3839/* Convert e-type X to IBM 370-format double E. */ 3840 3841static void 3842etoe53 (x, e) 3843 const UEMUSHORT *x; 3844 UEMUSHORT *e; 3845{ 3846 etoibm (x, e, DFmode); 3847} 3848 3849/* Convert exploded e-type X, that has already been rounded to 3850 56-bit precision, to IBM 370 double Y. */ 3851 3852static void 3853toe53 (x, y) 3854 UEMUSHORT *x, *y; 3855{ 3856 toibm (x, y, DFmode); 3857} 3858 3859#else /* it's neither DEC nor IBM */ 3860#ifdef C4X 3861/* Convert e-type X to C4X-format long double E. */ 3862 3863static void 3864etoe53 (x, e) 3865 const UEMUSHORT *x; 3866 UEMUSHORT *e; 3867{ 3868 etoc4x (x, e, HFmode); 3869} 3870 3871/* Convert exploded e-type X, that has already been rounded to 3872 56-bit precision, to IBM 370 double Y. */ 3873 3874static void 3875toe53 (x, y) 3876 UEMUSHORT *x, *y; 3877{ 3878 toc4x (x, y, HFmode); 3879} 3880 3881#else /* it's neither DEC nor IBM nor C4X */ 3882 3883/* Convert e-type X to IEEE double E. */ 3884 3885static void 3886etoe53 (x, e) 3887 const UEMUSHORT *x; 3888 UEMUSHORT *e; 3889{ 3890 UEMUSHORT xi[NI]; 3891 EMULONG exp; 3892 int rndsav; 3893 3894#ifdef NANS 3895 if (eisnan (x)) 3896 { 3897 make_nan (e, eisneg (x), DFmode); 3898 return; 3899 } 3900#endif 3901 emovi (x, xi); 3902 /* adjust exponent for offsets */ 3903 exp = (EMULONG) xi[E] - (EXONE - 0x3ff); 3904#ifdef INFINITY 3905 if (eisinf (x)) 3906 goto nonorm; 3907#endif 3908 /* round off to nearest or even */ 3909 rndsav = rndprc; 3910 rndprc = 53; 3911 emdnorm (xi, 0, 0, exp, 64); 3912 rndprc = rndsav; 3913#ifdef INFINITY 3914 nonorm: 3915#endif 3916 toe53 (xi, e); 3917} 3918 3919/* Convert exploded e-type X, that has already been rounded to 3920 53-bit precision, to IEEE double Y. */ 3921 3922static void 3923toe53 (x, y) 3924 UEMUSHORT *x, *y; 3925{ 3926 UEMUSHORT i; 3927 UEMUSHORT *p; 3928 3929#ifdef NANS 3930 if (eiisnan (x)) 3931 { 3932 make_nan (y, eiisneg (x), DFmode); 3933 return; 3934 } 3935#endif 3936 p = &x[0]; 3937#ifdef IEEE 3938 if (! REAL_WORDS_BIG_ENDIAN) 3939 y += 3; 3940#endif 3941 *y = 0; /* output high order */ 3942 if (*p++) 3943 *y = 0x8000; /* output sign bit */ 3944 3945 i = *p++; 3946 if (i >= (unsigned int) 2047) 3947 { 3948 /* Saturate at largest number less than infinity. */ 3949#ifdef INFINITY 3950 *y |= 0x7ff0; 3951 if (! REAL_WORDS_BIG_ENDIAN) 3952 { 3953 *(--y) = 0; 3954 *(--y) = 0; 3955 *(--y) = 0; 3956 } 3957 else 3958 { 3959 ++y; 3960 *y++ = 0; 3961 *y++ = 0; 3962 *y++ = 0; 3963 } 3964#else 3965 *y |= (UEMUSHORT) 0x7fef; 3966 if (! REAL_WORDS_BIG_ENDIAN) 3967 { 3968 *(--y) = 0xffff; 3969 *(--y) = 0xffff; 3970 *(--y) = 0xffff; 3971 } 3972 else 3973 { 3974 ++y; 3975 *y++ = 0xffff; 3976 *y++ = 0xffff; 3977 *y++ = 0xffff; 3978 } 3979#endif 3980 return; 3981 } 3982 if (i == 0) 3983 { 3984 eshift (x, 4); 3985 } 3986 else 3987 { 3988 i <<= 4; 3989 eshift (x, 5); 3990 } 3991 i |= *p++ & (UEMUSHORT) 0x0f; /* *p = xi[M] */ 3992 *y |= (UEMUSHORT) i; /* high order output already has sign bit set */ 3993 if (! REAL_WORDS_BIG_ENDIAN) 3994 { 3995 *(--y) = *p++; 3996 *(--y) = *p++; 3997 *(--y) = *p; 3998 } 3999 else 4000 { 4001 ++y; 4002 *y++ = *p++; 4003 *y++ = *p++; 4004 *y++ = *p++; 4005 } 4006} 4007 4008#endif /* not C4X */ 4009#endif /* not IBM */ 4010#endif /* not DEC */ 4011 4012 4013 4014/* e type to single precision. */ 4015 4016#ifdef IBM 4017/* Convert e-type X to IBM 370 float E. */ 4018 4019static void 4020etoe24 (x, e) 4021 const UEMUSHORT *x; 4022 UEMUSHORT *e; 4023{ 4024 etoibm (x, e, SFmode); 4025} 4026 4027/* Convert exploded e-type X, that has already been rounded to 4028 float precision, to IBM 370 float Y. */ 4029 4030static void 4031toe24 (x, y) 4032 UEMUSHORT *x, *y; 4033{ 4034 toibm (x, y, SFmode); 4035} 4036 4037#else 4038 4039#ifdef C4X 4040/* Convert e-type X to C4X float E. */ 4041 4042static void 4043etoe24 (x, e) 4044 const UEMUSHORT *x; 4045 UEMUSHORT *e; 4046{ 4047 etoc4x (x, e, QFmode); 4048} 4049 4050/* Convert exploded e-type X, that has already been rounded to 4051 float precision, to IBM 370 float Y. */ 4052 4053static void 4054toe24 (x, y) 4055 UEMUSHORT *x, *y; 4056{ 4057 toc4x (x, y, QFmode); 4058} 4059 4060#else 4061 4062/* Convert e-type X to IEEE float E. DEC float is the same as IEEE float. */ 4063 4064static void 4065etoe24 (x, e) 4066 const UEMUSHORT *x; 4067 UEMUSHORT *e; 4068{ 4069 EMULONG exp; 4070 UEMUSHORT xi[NI]; 4071 int rndsav; 4072 4073#ifdef NANS 4074 if (eisnan (x)) 4075 { 4076 make_nan (e, eisneg (x), SFmode); 4077 return; 4078 } 4079#endif 4080 emovi (x, xi); 4081 /* adjust exponent for offsets */ 4082 exp = (EMULONG) xi[E] - (EXONE - 0177); 4083#ifdef INFINITY 4084 if (eisinf (x)) 4085 goto nonorm; 4086#endif 4087 /* round off to nearest or even */ 4088 rndsav = rndprc; 4089 rndprc = 24; 4090 emdnorm (xi, 0, 0, exp, 64); 4091 rndprc = rndsav; 4092#ifdef INFINITY 4093 nonorm: 4094#endif 4095 toe24 (xi, e); 4096} 4097 4098/* Convert exploded e-type X, that has already been rounded to 4099 float precision, to IEEE float Y. */ 4100 4101static void 4102toe24 (x, y) 4103 UEMUSHORT *x, *y; 4104{ 4105 UEMUSHORT i; 4106 UEMUSHORT *p; 4107 4108#ifdef NANS 4109 if (eiisnan (x)) 4110 { 4111 make_nan (y, eiisneg (x), SFmode); 4112 return; 4113 } 4114#endif 4115 p = &x[0]; 4116#ifdef IEEE 4117 if (! REAL_WORDS_BIG_ENDIAN) 4118 y += 1; 4119#endif 4120#ifdef DEC 4121 y += 1; 4122#endif 4123 *y = 0; /* output high order */ 4124 if (*p++) 4125 *y = 0x8000; /* output sign bit */ 4126 4127 i = *p++; 4128/* Handle overflow cases. */ 4129 if (i >= 255) 4130 { 4131#ifdef INFINITY 4132 *y |= (UEMUSHORT) 0x7f80; 4133#ifdef DEC 4134 *(--y) = 0; 4135#endif 4136#ifdef IEEE 4137 if (! REAL_WORDS_BIG_ENDIAN) 4138 *(--y) = 0; 4139 else 4140 { 4141 ++y; 4142 *y = 0; 4143 } 4144#endif 4145#else /* no INFINITY */ 4146 *y |= (UEMUSHORT) 0x7f7f; 4147#ifdef DEC 4148 *(--y) = 0xffff; 4149#endif 4150#ifdef IEEE 4151 if (! REAL_WORDS_BIG_ENDIAN) 4152 *(--y) = 0xffff; 4153 else 4154 { 4155 ++y; 4156 *y = 0xffff; 4157 } 4158#endif 4159#ifdef ERANGE 4160 errno = ERANGE; 4161#endif 4162#endif /* no INFINITY */ 4163 return; 4164 } 4165 if (i == 0) 4166 { 4167 eshift (x, 7); 4168 } 4169 else 4170 { 4171 i <<= 7; 4172 eshift (x, 8); 4173 } 4174 i |= *p++ & (UEMUSHORT) 0x7f; /* *p = xi[M] */ 4175 /* High order output already has sign bit set. */ 4176 *y |= i; 4177#ifdef DEC 4178 *(--y) = *p; 4179#endif 4180#ifdef IEEE 4181 if (! REAL_WORDS_BIG_ENDIAN) 4182 *(--y) = *p; 4183 else 4184 { 4185 ++y; 4186 *y = *p; 4187 } 4188#endif 4189} 4190#endif /* not C4X */ 4191#endif /* not IBM */ 4192 4193/* Compare two e type numbers. 4194 Return +1 if a > b 4195 0 if a == b 4196 -1 if a < b 4197 -2 if either a or b is a NaN. */ 4198 4199static int 4200ecmp (a, b) 4201 const UEMUSHORT *a, *b; 4202{ 4203 UEMUSHORT ai[NI], bi[NI]; 4204 UEMUSHORT *p, *q; 4205 int i; 4206 int msign; 4207 4208#ifdef NANS 4209 if (eisnan (a) || eisnan (b)) 4210 return (-2); 4211#endif 4212 emovi (a, ai); 4213 p = ai; 4214 emovi (b, bi); 4215 q = bi; 4216 4217 if (*p != *q) 4218 { /* the signs are different */ 4219 /* -0 equals + 0 */ 4220 for (i = 1; i < NI - 1; i++) 4221 { 4222 if (ai[i] != 0) 4223 goto nzro; 4224 if (bi[i] != 0) 4225 goto nzro; 4226 } 4227 return (0); 4228 nzro: 4229 if (*p == 0) 4230 return (1); 4231 else 4232 return (-1); 4233 } 4234 /* both are the same sign */ 4235 if (*p == 0) 4236 msign = 1; 4237 else 4238 msign = -1; 4239 i = NI - 1; 4240 do 4241 { 4242 if (*p++ != *q++) 4243 { 4244 goto diff; 4245 } 4246 } 4247 while (--i > 0); 4248 4249 return (0); /* equality */ 4250 4251 diff: 4252 4253 if (*(--p) > *(--q)) 4254 return (msign); /* p is bigger */ 4255 else 4256 return (-msign); /* p is littler */ 4257} 4258 4259#if 0 4260/* Find e-type nearest integer to X, as floor (X + 0.5). */ 4261 4262static void 4263eround (x, y) 4264 const UEMUSHORT *x; 4265 UEMUSHORT *y; 4266{ 4267 eadd (ehalf, x, y); 4268 efloor (y, y); 4269} 4270#endif /* 0 */ 4271 4272/* Convert HOST_WIDE_INT LP to e type Y. */ 4273 4274static void 4275ltoe (lp, y) 4276 const HOST_WIDE_INT *lp; 4277 UEMUSHORT *y; 4278{ 4279 UEMUSHORT yi[NI]; 4280 unsigned HOST_WIDE_INT ll; 4281 int k; 4282 4283 ecleaz (yi); 4284 if (*lp < 0) 4285 { 4286 /* make it positive */ 4287 ll = (unsigned HOST_WIDE_INT) (-(*lp)); 4288 yi[0] = 0xffff; /* put correct sign in the e type number */ 4289 } 4290 else 4291 { 4292 ll = (unsigned HOST_WIDE_INT) (*lp); 4293 } 4294 /* move the long integer to yi significand area */ 4295#if HOST_BITS_PER_WIDE_INT == 64 4296 yi[M] = (UEMUSHORT) (ll >> 48); 4297 yi[M + 1] = (UEMUSHORT) (ll >> 32); 4298 yi[M + 2] = (UEMUSHORT) (ll >> 16); 4299 yi[M + 3] = (UEMUSHORT) ll; 4300 yi[E] = EXONE + 47; /* exponent if normalize shift count were 0 */ 4301#else 4302 yi[M] = (UEMUSHORT) (ll >> 16); 4303 yi[M + 1] = (UEMUSHORT) ll; 4304 yi[E] = EXONE + 15; /* exponent if normalize shift count were 0 */ 4305#endif 4306 4307 if ((k = enormlz (yi)) > NBITS)/* normalize the significand */ 4308 ecleaz (yi); /* it was zero */ 4309 else 4310 yi[E] -= (UEMUSHORT) k;/* subtract shift count from exponent */ 4311 emovo (yi, y); /* output the answer */ 4312} 4313 4314/* Convert unsigned HOST_WIDE_INT LP to e type Y. */ 4315 4316static void 4317ultoe (lp, y) 4318 const unsigned HOST_WIDE_INT *lp; 4319 UEMUSHORT *y; 4320{ 4321 UEMUSHORT yi[NI]; 4322 unsigned HOST_WIDE_INT ll; 4323 int k; 4324 4325 ecleaz (yi); 4326 ll = *lp; 4327 4328 /* move the long integer to ayi significand area */ 4329#if HOST_BITS_PER_WIDE_INT == 64 4330 yi[M] = (UEMUSHORT) (ll >> 48); 4331 yi[M + 1] = (UEMUSHORT) (ll >> 32); 4332 yi[M + 2] = (UEMUSHORT) (ll >> 16); 4333 yi[M + 3] = (UEMUSHORT) ll; 4334 yi[E] = EXONE + 47; /* exponent if normalize shift count were 0 */ 4335#else 4336 yi[M] = (UEMUSHORT) (ll >> 16); 4337 yi[M + 1] = (UEMUSHORT) ll; 4338 yi[E] = EXONE + 15; /* exponent if normalize shift count were 0 */ 4339#endif 4340 4341 if ((k = enormlz (yi)) > NBITS)/* normalize the significand */ 4342 ecleaz (yi); /* it was zero */ 4343 else 4344 yi[E] -= (UEMUSHORT) k; /* subtract shift count from exponent */ 4345 emovo (yi, y); /* output the answer */ 4346} 4347 4348 4349/* Find signed HOST_WIDE_INT integer I and floating point fractional 4350 part FRAC of e-type (packed internal format) floating point input X. 4351 The integer output I has the sign of the input, except that 4352 positive overflow is permitted if FIXUNS_TRUNC_LIKE_FIX_TRUNC. 4353 The output e-type fraction FRAC is the positive fractional 4354 part of abs (X). */ 4355 4356static void 4357eifrac (x, i, frac) 4358 const UEMUSHORT *x; 4359 HOST_WIDE_INT *i; 4360 UEMUSHORT *frac; 4361{ 4362 UEMUSHORT xi[NI]; 4363 int j, k; 4364 unsigned HOST_WIDE_INT ll; 4365 4366 emovi (x, xi); 4367 k = (int) xi[E] - (EXONE - 1); 4368 if (k <= 0) 4369 { 4370 /* if exponent <= 0, integer = 0 and real output is fraction */ 4371 *i = 0L; 4372 emovo (xi, frac); 4373 return; 4374 } 4375 if (k > (HOST_BITS_PER_WIDE_INT - 1)) 4376 { 4377 /* long integer overflow: output large integer 4378 and correct fraction */ 4379 if (xi[0]) 4380 *i = ((unsigned HOST_WIDE_INT) 1) << (HOST_BITS_PER_WIDE_INT - 1); 4381 else 4382 { 4383#ifdef FIXUNS_TRUNC_LIKE_FIX_TRUNC 4384 /* In this case, let it overflow and convert as if unsigned. */ 4385 euifrac (x, &ll, frac); 4386 *i = (HOST_WIDE_INT) ll; 4387 return; 4388#else 4389 /* In other cases, return the largest positive integer. */ 4390 *i = (((unsigned HOST_WIDE_INT) 1) << (HOST_BITS_PER_WIDE_INT - 1)) - 1; 4391#endif 4392 } 4393 eshift (xi, k); 4394 if (extra_warnings) 4395 warning ("overflow on truncation to integer"); 4396 } 4397 else if (k > 16) 4398 { 4399 /* Shift more than 16 bits: first shift up k-16 mod 16, 4400 then shift up by 16's. */ 4401 j = k - ((k >> 4) << 4); 4402 eshift (xi, j); 4403 ll = xi[M]; 4404 k -= j; 4405 do 4406 { 4407 eshup6 (xi); 4408 ll = (ll << 16) | xi[M]; 4409 } 4410 while ((k -= 16) > 0); 4411 *i = ll; 4412 if (xi[0]) 4413 *i = -(*i); 4414 } 4415 else 4416 { 4417 /* shift not more than 16 bits */ 4418 eshift (xi, k); 4419 *i = (HOST_WIDE_INT) xi[M] & 0xffff; 4420 if (xi[0]) 4421 *i = -(*i); 4422 } 4423 xi[0] = 0; 4424 xi[E] = EXONE - 1; 4425 xi[M] = 0; 4426 if ((k = enormlz (xi)) > NBITS) 4427 ecleaz (xi); 4428 else 4429 xi[E] -= (UEMUSHORT) k; 4430 4431 emovo (xi, frac); 4432} 4433 4434 4435/* Find unsigned HOST_WIDE_INT integer I and floating point fractional part 4436 FRAC of e-type X. A negative input yields integer output = 0 but 4437 correct fraction. */ 4438 4439static void 4440euifrac (x, i, frac) 4441 const UEMUSHORT *x; 4442 unsigned HOST_WIDE_INT *i; 4443 UEMUSHORT *frac; 4444{ 4445 unsigned HOST_WIDE_INT ll; 4446 UEMUSHORT xi[NI]; 4447 int j, k; 4448 4449 emovi (x, xi); 4450 k = (int) xi[E] - (EXONE - 1); 4451 if (k <= 0) 4452 { 4453 /* if exponent <= 0, integer = 0 and argument is fraction */ 4454 *i = 0L; 4455 emovo (xi, frac); 4456 return; 4457 } 4458 if (k > HOST_BITS_PER_WIDE_INT) 4459 { 4460 /* Long integer overflow: output large integer 4461 and correct fraction. 4462 Note, the BSD MicroVAX compiler says that ~(0UL) 4463 is a syntax error. */ 4464 *i = ~(0L); 4465 eshift (xi, k); 4466 if (extra_warnings) 4467 warning ("overflow on truncation to unsigned integer"); 4468 } 4469 else if (k > 16) 4470 { 4471 /* Shift more than 16 bits: first shift up k-16 mod 16, 4472 then shift up by 16's. */ 4473 j = k - ((k >> 4) << 4); 4474 eshift (xi, j); 4475 ll = xi[M]; 4476 k -= j; 4477 do 4478 { 4479 eshup6 (xi); 4480 ll = (ll << 16) | xi[M]; 4481 } 4482 while ((k -= 16) > 0); 4483 *i = ll; 4484 } 4485 else 4486 { 4487 /* shift not more than 16 bits */ 4488 eshift (xi, k); 4489 *i = (HOST_WIDE_INT) xi[M] & 0xffff; 4490 } 4491 4492 if (xi[0]) /* A negative value yields unsigned integer 0. */ 4493 *i = 0L; 4494 4495 xi[0] = 0; 4496 xi[E] = EXONE - 1; 4497 xi[M] = 0; 4498 if ((k = enormlz (xi)) > NBITS) 4499 ecleaz (xi); 4500 else 4501 xi[E] -= (UEMUSHORT) k; 4502 4503 emovo (xi, frac); 4504} 4505 4506/* Shift the significand of exploded e-type X up or down by SC bits. */ 4507 4508static int 4509eshift (x, sc) 4510 UEMUSHORT *x; 4511 int sc; 4512{ 4513 UEMUSHORT lost; 4514 UEMUSHORT *p; 4515 4516 if (sc == 0) 4517 return (0); 4518 4519 lost = 0; 4520 p = x + NI - 1; 4521 4522 if (sc < 0) 4523 { 4524 sc = -sc; 4525 while (sc >= 16) 4526 { 4527 lost |= *p; /* remember lost bits */ 4528 eshdn6 (x); 4529 sc -= 16; 4530 } 4531 4532 while (sc >= 8) 4533 { 4534 lost |= *p & 0xff; 4535 eshdn8 (x); 4536 sc -= 8; 4537 } 4538 4539 while (sc > 0) 4540 { 4541 lost |= *p & 1; 4542 eshdn1 (x); 4543 sc -= 1; 4544 } 4545 } 4546 else 4547 { 4548 while (sc >= 16) 4549 { 4550 eshup6 (x); 4551 sc -= 16; 4552 } 4553 4554 while (sc >= 8) 4555 { 4556 eshup8 (x); 4557 sc -= 8; 4558 } 4559 4560 while (sc > 0) 4561 { 4562 eshup1 (x); 4563 sc -= 1; 4564 } 4565 } 4566 if (lost) 4567 lost = 1; 4568 return ((int) lost); 4569} 4570 4571/* Shift normalize the significand area of exploded e-type X. 4572 Return the shift count (up = positive). */ 4573 4574static int 4575enormlz (x) 4576 UEMUSHORT x[]; 4577{ 4578 UEMUSHORT *p; 4579 int sc; 4580 4581 sc = 0; 4582 p = &x[M]; 4583 if (*p != 0) 4584 goto normdn; 4585 ++p; 4586 if (*p & 0x8000) 4587 return (0); /* already normalized */ 4588 while (*p == 0) 4589 { 4590 eshup6 (x); 4591 sc += 16; 4592 4593 /* With guard word, there are NBITS+16 bits available. 4594 Return true if all are zero. */ 4595 if (sc > NBITS) 4596 return (sc); 4597 } 4598 /* see if high byte is zero */ 4599 while ((*p & 0xff00) == 0) 4600 { 4601 eshup8 (x); 4602 sc += 8; 4603 } 4604 /* now shift 1 bit at a time */ 4605 while ((*p & 0x8000) == 0) 4606 { 4607 eshup1 (x); 4608 sc += 1; 4609 if (sc > NBITS) 4610 { 4611 mtherr ("enormlz", UNDERFLOW); 4612 return (sc); 4613 } 4614 } 4615 return (sc); 4616 4617 /* Normalize by shifting down out of the high guard word 4618 of the significand */ 4619 normdn: 4620 4621 if (*p & 0xff00) 4622 { 4623 eshdn8 (x); 4624 sc -= 8; 4625 } 4626 while (*p != 0) 4627 { 4628 eshdn1 (x); 4629 sc -= 1; 4630 4631 if (sc < -NBITS) 4632 { 4633 mtherr ("enormlz", OVERFLOW); 4634 return (sc); 4635 } 4636 } 4637 return (sc); 4638} 4639 4640/* Powers of ten used in decimal <-> binary conversions. */ 4641 4642#define NTEN 12 4643#define MAXP 4096 4644 4645#if MAX_LONG_DOUBLE_TYPE_SIZE == 128 && (INTEL_EXTENDED_IEEE_FORMAT == 0) 4646static const UEMUSHORT etens[NTEN + 1][NE] = 4647{ 4648 {0x6576, 0x4a92, 0x804a, 0x153f, 4649 0xc94c, 0x979a, 0x8a20, 0x5202, 0xc460, 0x7525,}, /* 10**4096 */ 4650 {0x6a32, 0xce52, 0x329a, 0x28ce, 4651 0xa74d, 0x5de4, 0xc53d, 0x3b5d, 0x9e8b, 0x5a92,}, /* 10**2048 */ 4652 {0x526c, 0x50ce, 0xf18b, 0x3d28, 4653 0x650d, 0x0c17, 0x8175, 0x7586, 0xc976, 0x4d48,}, 4654 {0x9c66, 0x58f8, 0xbc50, 0x5c54, 4655 0xcc65, 0x91c6, 0xa60e, 0xa0ae, 0xe319, 0x46a3,}, 4656 {0x851e, 0xeab7, 0x98fe, 0x901b, 4657 0xddbb, 0xde8d, 0x9df9, 0xebfb, 0xaa7e, 0x4351,}, 4658 {0x0235, 0x0137, 0x36b1, 0x336c, 4659 0xc66f, 0x8cdf, 0x80e9, 0x47c9, 0x93ba, 0x41a8,}, 4660 {0x50f8, 0x25fb, 0xc76b, 0x6b71, 4661 0x3cbf, 0xa6d5, 0xffcf, 0x1f49, 0xc278, 0x40d3,}, 4662 {0x0000, 0x0000, 0x0000, 0x0000, 4663 0xf020, 0xb59d, 0x2b70, 0xada8, 0x9dc5, 0x4069,}, 4664 {0x0000, 0x0000, 0x0000, 0x0000, 4665 0x0000, 0x0000, 0x0400, 0xc9bf, 0x8e1b, 0x4034,}, 4666 {0x0000, 0x0000, 0x0000, 0x0000, 4667 0x0000, 0x0000, 0x0000, 0x2000, 0xbebc, 0x4019,}, 4668 {0x0000, 0x0000, 0x0000, 0x0000, 4669 0x0000, 0x0000, 0x0000, 0x0000, 0x9c40, 0x400c,}, 4670 {0x0000, 0x0000, 0x0000, 0x0000, 4671 0x0000, 0x0000, 0x0000, 0x0000, 0xc800, 0x4005,}, 4672 {0x0000, 0x0000, 0x0000, 0x0000, 4673 0x0000, 0x0000, 0x0000, 0x0000, 0xa000, 0x4002,}, /* 10**1 */ 4674}; 4675 4676static const UEMUSHORT emtens[NTEN + 1][NE] = 4677{ 4678 {0x2030, 0xcffc, 0xa1c3, 0x8123, 4679 0x2de3, 0x9fde, 0xd2ce, 0x04c8, 0xa6dd, 0x0ad8,}, /* 10**-4096 */ 4680 {0x8264, 0xd2cb, 0xf2ea, 0x12d4, 4681 0x4925, 0x2de4, 0x3436, 0x534f, 0xceae, 0x256b,}, /* 10**-2048 */ 4682 {0xf53f, 0xf698, 0x6bd3, 0x0158, 4683 0x87a6, 0xc0bd, 0xda57, 0x82a5, 0xa2a6, 0x32b5,}, 4684 {0xe731, 0x04d4, 0xe3f2, 0xd332, 4685 0x7132, 0xd21c, 0xdb23, 0xee32, 0x9049, 0x395a,}, 4686 {0xa23e, 0x5308, 0xfefb, 0x1155, 4687 0xfa91, 0x1939, 0x637a, 0x4325, 0xc031, 0x3cac,}, 4688 {0xe26d, 0xdbde, 0xd05d, 0xb3f6, 4689 0xac7c, 0xe4a0, 0x64bc, 0x467c, 0xddd0, 0x3e55,}, 4690 {0x2a20, 0x6224, 0x47b3, 0x98d7, 4691 0x3f23, 0xe9a5, 0xa539, 0xea27, 0xa87f, 0x3f2a,}, 4692 {0x0b5b, 0x4af2, 0xa581, 0x18ed, 4693 0x67de, 0x94ba, 0x4539, 0x1ead, 0xcfb1, 0x3f94,}, 4694 {0xbf71, 0xa9b3, 0x7989, 0xbe68, 4695 0x4c2e, 0xe15b, 0xc44d, 0x94be, 0xe695, 0x3fc9,}, 4696 {0x3d4d, 0x7c3d, 0x36ba, 0x0d2b, 4697 0xfdc2, 0xcefc, 0x8461, 0x7711, 0xabcc, 0x3fe4,}, 4698 {0xc155, 0xa4a8, 0x404e, 0x6113, 4699 0xd3c3, 0x652b, 0xe219, 0x1758, 0xd1b7, 0x3ff1,}, 4700 {0xd70a, 0x70a3, 0x0a3d, 0xa3d7, 4701 0x3d70, 0xd70a, 0x70a3, 0x0a3d, 0xa3d7, 0x3ff8,}, 4702 {0xcccd, 0xcccc, 0xcccc, 0xcccc, 4703 0xcccc, 0xcccc, 0xcccc, 0xcccc, 0xcccc, 0x3ffb,}, /* 10**-1 */ 4704}; 4705#else 4706/* LONG_DOUBLE_TYPE_SIZE is other than 128 */ 4707static const UEMUSHORT etens[NTEN + 1][NE] = 4708{ 4709 {0xc94c, 0x979a, 0x8a20, 0x5202, 0xc460, 0x7525,}, /* 10**4096 */ 4710 {0xa74d, 0x5de4, 0xc53d, 0x3b5d, 0x9e8b, 0x5a92,}, /* 10**2048 */ 4711 {0x650d, 0x0c17, 0x8175, 0x7586, 0xc976, 0x4d48,}, 4712 {0xcc65, 0x91c6, 0xa60e, 0xa0ae, 0xe319, 0x46a3,}, 4713 {0xddbc, 0xde8d, 0x9df9, 0xebfb, 0xaa7e, 0x4351,}, 4714 {0xc66f, 0x8cdf, 0x80e9, 0x47c9, 0x93ba, 0x41a8,}, 4715 {0x3cbf, 0xa6d5, 0xffcf, 0x1f49, 0xc278, 0x40d3,}, 4716 {0xf020, 0xb59d, 0x2b70, 0xada8, 0x9dc5, 0x4069,}, 4717 {0x0000, 0x0000, 0x0400, 0xc9bf, 0x8e1b, 0x4034,}, 4718 {0x0000, 0x0000, 0x0000, 0x2000, 0xbebc, 0x4019,}, 4719 {0x0000, 0x0000, 0x0000, 0x0000, 0x9c40, 0x400c,}, 4720 {0x0000, 0x0000, 0x0000, 0x0000, 0xc800, 0x4005,}, 4721 {0x0000, 0x0000, 0x0000, 0x0000, 0xa000, 0x4002,}, /* 10**1 */ 4722}; 4723 4724static const UEMUSHORT emtens[NTEN + 1][NE] = 4725{ 4726 {0x2de4, 0x9fde, 0xd2ce, 0x04c8, 0xa6dd, 0x0ad8,}, /* 10**-4096 */ 4727 {0x4925, 0x2de4, 0x3436, 0x534f, 0xceae, 0x256b,}, /* 10**-2048 */ 4728 {0x87a6, 0xc0bd, 0xda57, 0x82a5, 0xa2a6, 0x32b5,}, 4729 {0x7133, 0xd21c, 0xdb23, 0xee32, 0x9049, 0x395a,}, 4730 {0xfa91, 0x1939, 0x637a, 0x4325, 0xc031, 0x3cac,}, 4731 {0xac7d, 0xe4a0, 0x64bc, 0x467c, 0xddd0, 0x3e55,}, 4732 {0x3f24, 0xe9a5, 0xa539, 0xea27, 0xa87f, 0x3f2a,}, 4733 {0x67de, 0x94ba, 0x4539, 0x1ead, 0xcfb1, 0x3f94,}, 4734 {0x4c2f, 0xe15b, 0xc44d, 0x94be, 0xe695, 0x3fc9,}, 4735 {0xfdc2, 0xcefc, 0x8461, 0x7711, 0xabcc, 0x3fe4,}, 4736 {0xd3c3, 0x652b, 0xe219, 0x1758, 0xd1b7, 0x3ff1,}, 4737 {0x3d71, 0xd70a, 0x70a3, 0x0a3d, 0xa3d7, 0x3ff8,}, 4738 {0xcccd, 0xcccc, 0xcccc, 0xcccc, 0xcccc, 0x3ffb,}, /* 10**-1 */ 4739}; 4740#endif 4741 4742#if 0 4743/* Convert float value X to ASCII string STRING with NDIG digits after 4744 the decimal point. */ 4745 4746static void 4747e24toasc (x, string, ndigs) 4748 const UEMUSHORT x[]; 4749 char *string; 4750 int ndigs; 4751{ 4752 UEMUSHORT w[NI]; 4753 4754 e24toe (x, w); 4755 etoasc (w, string, ndigs); 4756} 4757 4758/* Convert double value X to ASCII string STRING with NDIG digits after 4759 the decimal point. */ 4760 4761static void 4762e53toasc (x, string, ndigs) 4763 const UEMUSHORT x[]; 4764 char *string; 4765 int ndigs; 4766{ 4767 UEMUSHORT w[NI]; 4768 4769 e53toe (x, w); 4770 etoasc (w, string, ndigs); 4771} 4772 4773/* Convert double extended value X to ASCII string STRING with NDIG digits 4774 after the decimal point. */ 4775 4776static void 4777e64toasc (x, string, ndigs) 4778 const UEMUSHORT x[]; 4779 char *string; 4780 int ndigs; 4781{ 4782 UEMUSHORT w[NI]; 4783 4784 e64toe (x, w); 4785 etoasc (w, string, ndigs); 4786} 4787 4788/* Convert 128-bit long double value X to ASCII string STRING with NDIG digits 4789 after the decimal point. */ 4790 4791static void 4792e113toasc (x, string, ndigs) 4793 const UEMUSHORT x[]; 4794 char *string; 4795 int ndigs; 4796{ 4797 UEMUSHORT w[NI]; 4798 4799 e113toe (x, w); 4800 etoasc (w, string, ndigs); 4801} 4802#endif /* 0 */ 4803 4804/* Convert e-type X to ASCII string STRING with NDIGS digits after 4805 the decimal point. */ 4806 4807static char wstring[80]; /* working storage for ASCII output */ 4808 4809static void 4810etoasc (x, string, ndigs) 4811 const UEMUSHORT x[]; 4812 char *string; 4813 int ndigs; 4814{ 4815 EMUSHORT digit; 4816 UEMUSHORT y[NI], t[NI], u[NI], w[NI]; 4817 const UEMUSHORT *p, *r, *ten; 4818 UEMUSHORT sign; 4819 int i, j, k, expon, rndsav; 4820 char *s, *ss; 4821 UEMUSHORT m; 4822 4823 4824 rndsav = rndprc; 4825 ss = string; 4826 s = wstring; 4827 *ss = '\0'; 4828 *s = '\0'; 4829#ifdef NANS 4830 if (eisnan (x)) 4831 { 4832 sprintf (wstring, " NaN "); 4833 goto bxit; 4834 } 4835#endif 4836 rndprc = NBITS; /* set to full precision */ 4837 emov (x, y); /* retain external format */ 4838 if (y[NE - 1] & 0x8000) 4839 { 4840 sign = 0xffff; 4841 y[NE - 1] &= 0x7fff; 4842 } 4843 else 4844 { 4845 sign = 0; 4846 } 4847 expon = 0; 4848 ten = &etens[NTEN][0]; 4849 emov (eone, t); 4850 /* Test for zero exponent */ 4851 if (y[NE - 1] == 0) 4852 { 4853 for (k = 0; k < NE - 1; k++) 4854 { 4855 if (y[k] != 0) 4856 goto tnzro; /* denormalized number */ 4857 } 4858 goto isone; /* valid all zeros */ 4859 } 4860 tnzro: 4861 4862 /* Test for infinity. */ 4863 if (y[NE - 1] == 0x7fff) 4864 { 4865 if (sign) 4866 sprintf (wstring, " -Infinity "); 4867 else 4868 sprintf (wstring, " Infinity "); 4869 goto bxit; 4870 } 4871 4872 /* Test for exponent nonzero but significand denormalized. 4873 * This is an error condition. 4874 */ 4875 if ((y[NE - 1] != 0) && ((y[NE - 2] & 0x8000) == 0)) 4876 { 4877 mtherr ("etoasc", DOMAIN); 4878 sprintf (wstring, "NaN"); 4879 goto bxit; 4880 } 4881 4882 /* Compare to 1.0 */ 4883 i = ecmp (eone, y); 4884 if (i == 0) 4885 goto isone; 4886 4887 if (i == -2) 4888 abort (); 4889 4890 if (i < 0) 4891 { /* Number is greater than 1 */ 4892 /* Convert significand to an integer and strip trailing decimal zeros. */ 4893 emov (y, u); 4894 u[NE - 1] = EXONE + NBITS - 1; 4895 4896 p = &etens[NTEN - 4][0]; 4897 m = 16; 4898 do 4899 { 4900 ediv (p, u, t); 4901 efloor (t, w); 4902 for (j = 0; j < NE - 1; j++) 4903 { 4904 if (t[j] != w[j]) 4905 goto noint; 4906 } 4907 emov (t, u); 4908 expon += (int) m; 4909 noint: 4910 p += NE; 4911 m >>= 1; 4912 } 4913 while (m != 0); 4914 4915 /* Rescale from integer significand */ 4916 u[NE - 1] += y[NE - 1] - (unsigned int) (EXONE + NBITS - 1); 4917 emov (u, y); 4918 /* Find power of 10 */ 4919 emov (eone, t); 4920 m = MAXP; 4921 p = &etens[0][0]; 4922 /* An unordered compare result shouldn't happen here. */ 4923 while (ecmp (ten, u) <= 0) 4924 { 4925 if (ecmp (p, u) <= 0) 4926 { 4927 ediv (p, u, u); 4928 emul (p, t, t); 4929 expon += (int) m; 4930 } 4931 m >>= 1; 4932 if (m == 0) 4933 break; 4934 p += NE; 4935 } 4936 } 4937 else 4938 { /* Number is less than 1.0 */ 4939 /* Pad significand with trailing decimal zeros. */ 4940 if (y[NE - 1] == 0) 4941 { 4942 while ((y[NE - 2] & 0x8000) == 0) 4943 { 4944 emul (ten, y, y); 4945 expon -= 1; 4946 } 4947 } 4948 else 4949 { 4950 emovi (y, w); 4951 for (i = 0; i < NDEC + 1; i++) 4952 { 4953 if ((w[NI - 1] & 0x7) != 0) 4954 break; 4955 /* multiply by 10 */ 4956 emovz (w, u); 4957 eshdn1 (u); 4958 eshdn1 (u); 4959 eaddm (w, u); 4960 u[1] += 3; 4961 while (u[2] != 0) 4962 { 4963 eshdn1 (u); 4964 u[1] += 1; 4965 } 4966 if (u[NI - 1] != 0) 4967 break; 4968 if (eone[NE - 1] <= u[1]) 4969 break; 4970 emovz (u, w); 4971 expon -= 1; 4972 } 4973 emovo (w, y); 4974 } 4975 k = -MAXP; 4976 p = &emtens[0][0]; 4977 r = &etens[0][0]; 4978 emov (y, w); 4979 emov (eone, t); 4980 while (ecmp (eone, w) > 0) 4981 { 4982 if (ecmp (p, w) >= 0) 4983 { 4984 emul (r, w, w); 4985 emul (r, t, t); 4986 expon += k; 4987 } 4988 k /= 2; 4989 if (k == 0) 4990 break; 4991 p += NE; 4992 r += NE; 4993 } 4994 ediv (t, eone, t); 4995 } 4996 isone: 4997 /* Find the first (leading) digit. */ 4998 emovi (t, w); 4999 emovz (w, t); 5000 emovi (y, w); 5001 emovz (w, y); 5002 eiremain (t, y); 5003 digit = equot[NI - 1]; 5004 while ((digit == 0) && (ecmp (y, ezero) != 0)) 5005 { 5006 eshup1 (y); 5007 emovz (y, u); 5008 eshup1 (u); 5009 eshup1 (u); 5010 eaddm (u, y); 5011 eiremain (t, y); 5012 digit = equot[NI - 1]; 5013 expon -= 1; 5014 } 5015 s = wstring; 5016 if (sign) 5017 *s++ = '-'; 5018 else 5019 *s++ = ' '; 5020 /* Examine number of digits requested by caller. */ 5021 if (ndigs < 0) 5022 ndigs = 0; 5023 if (ndigs > NDEC) 5024 ndigs = NDEC; 5025 if (digit == 10) 5026 { 5027 *s++ = '1'; 5028 *s++ = '.'; 5029 if (ndigs > 0) 5030 { 5031 *s++ = '0'; 5032 ndigs -= 1; 5033 } 5034 expon += 1; 5035 } 5036 else 5037 { 5038 *s++ = (char) digit + '0'; 5039 *s++ = '.'; 5040 } 5041 /* Generate digits after the decimal point. */ 5042 for (k = 0; k <= ndigs; k++) 5043 { 5044 /* multiply current number by 10, without normalizing */ 5045 eshup1 (y); 5046 emovz (y, u); 5047 eshup1 (u); 5048 eshup1 (u); 5049 eaddm (u, y); 5050 eiremain (t, y); 5051 *s++ = (char) equot[NI - 1] + '0'; 5052 } 5053 digit = equot[NI - 1]; 5054 --s; 5055 ss = s; 5056 /* round off the ASCII string */ 5057 if (digit > 4) 5058 { 5059 /* Test for critical rounding case in ASCII output. */ 5060 if (digit == 5) 5061 { 5062 emovo (y, t); 5063 if (ecmp (t, ezero) != 0) 5064 goto roun; /* round to nearest */ 5065#ifndef C4X 5066 if ((*(s - 1) & 1) == 0) 5067 goto doexp; /* round to even */ 5068#endif 5069 } 5070 /* Round up and propagate carry-outs */ 5071 roun: 5072 --s; 5073 k = *s & CHARMASK; 5074 /* Carry out to most significant digit? */ 5075 if (k == '.') 5076 { 5077 --s; 5078 k = *s; 5079 k += 1; 5080 *s = (char) k; 5081 /* Most significant digit carries to 10? */ 5082 if (k > '9') 5083 { 5084 expon += 1; 5085 *s = '1'; 5086 } 5087 goto doexp; 5088 } 5089 /* Round up and carry out from less significant digits */ 5090 k += 1; 5091 *s = (char) k; 5092 if (k > '9') 5093 { 5094 *s = '0'; 5095 goto roun; 5096 } 5097 } 5098 doexp: 5099 /* 5100 if (expon >= 0) 5101 sprintf (ss, "e+%d", expon); 5102 else 5103 sprintf (ss, "e%d", expon); 5104 */ 5105 sprintf (ss, "e%d", expon); 5106 bxit: 5107 rndprc = rndsav; 5108 /* copy out the working string */ 5109 s = string; 5110 ss = wstring; 5111 while (*ss == ' ') /* strip possible leading space */ 5112 ++ss; 5113 while ((*s++ = *ss++) != '\0') 5114 ; 5115} 5116 5117 5118/* Convert ASCII string to floating point. 5119 5120 Numeric input is a free format decimal number of any length, with 5121 or without decimal point. Entering E after the number followed by an 5122 integer number causes the second number to be interpreted as a power of 5123 10 to be multiplied by the first number (i.e., "scientific" notation). */ 5124 5125/* Convert ASCII string S to single precision float value Y. */ 5126 5127static void 5128asctoe24 (s, y) 5129 const char *s; 5130 UEMUSHORT *y; 5131{ 5132 asctoeg (s, y, 24); 5133} 5134 5135 5136/* Convert ASCII string S to double precision value Y. */ 5137 5138static void 5139asctoe53 (s, y) 5140 const char *s; 5141 UEMUSHORT *y; 5142{ 5143#if defined(DEC) || defined(IBM) 5144 asctoeg (s, y, 56); 5145#else 5146#if defined(C4X) 5147 asctoeg (s, y, 32); 5148#else 5149 asctoeg (s, y, 53); 5150#endif 5151#endif 5152} 5153 5154 5155/* Convert ASCII string S to double extended value Y. */ 5156 5157static void 5158asctoe64 (s, y) 5159 const char *s; 5160 UEMUSHORT *y; 5161{ 5162 asctoeg (s, y, 64); 5163} 5164 5165#if (INTEL_EXTENDED_IEEE_FORMAT == 0) 5166/* Convert ASCII string S to 128-bit long double Y. */ 5167 5168static void 5169asctoe113 (s, y) 5170 const char *s; 5171 UEMUSHORT *y; 5172{ 5173 asctoeg (s, y, 113); 5174} 5175#endif 5176 5177/* Convert ASCII string S to e type Y. */ 5178 5179static void 5180asctoe (s, y) 5181 const char *s; 5182 UEMUSHORT *y; 5183{ 5184 asctoeg (s, y, NBITS); 5185} 5186 5187/* Convert ASCII string SS to e type Y, with a specified rounding precision 5188 of OPREC bits. BASE is 16 for C99 hexadecimal floating constants. */ 5189 5190static void 5191asctoeg (ss, y, oprec) 5192 const char *ss; 5193 UEMUSHORT *y; 5194 int oprec; 5195{ 5196 UEMUSHORT yy[NI], xt[NI], tt[NI]; 5197 int esign, decflg, sgnflg, nexp, exp, prec, lost; 5198 int i, k, trail, c, rndsav; 5199 EMULONG lexp; 5200 UEMUSHORT nsign; 5201 char *sp, *s, *lstr; 5202 int base = 10; 5203 5204 /* Copy the input string. */ 5205 lstr = (char *) alloca (strlen (ss) + 1); 5206 5207 while (*ss == ' ') /* skip leading spaces */ 5208 ++ss; 5209 5210 sp = lstr; 5211 while ((*sp++ = *ss++) != '\0') 5212 ; 5213 s = lstr; 5214 5215 if (s[0] == '0' && (s[1] == 'x' || s[1] == 'X')) 5216 { 5217 base = 16; 5218 s += 2; 5219 } 5220 5221 rndsav = rndprc; 5222 rndprc = NBITS; /* Set to full precision */ 5223 lost = 0; 5224 nsign = 0; 5225 decflg = 0; 5226 sgnflg = 0; 5227 nexp = 0; 5228 exp = 0; 5229 prec = 0; 5230 ecleaz (yy); 5231 trail = 0; 5232 5233 nxtcom: 5234 k = hex_value (*s); 5235 if ((k >= 0) && (k < base)) 5236 { 5237 /* Ignore leading zeros */ 5238 if ((prec == 0) && (decflg == 0) && (k == 0)) 5239 goto donchr; 5240 /* Identify and strip trailing zeros after the decimal point. */ 5241 if ((trail == 0) && (decflg != 0)) 5242 { 5243 sp = s; 5244 while (ISDIGIT (*sp) || (base == 16 && ISXDIGIT (*sp))) 5245 ++sp; 5246 /* Check for syntax error */ 5247 c = *sp & CHARMASK; 5248 if ((base != 10 || ((c != 'e') && (c != 'E'))) 5249 && (base != 16 || ((c != 'p') && (c != 'P'))) 5250 && (c != '\0') 5251 && (c != '\n') && (c != '\r') && (c != ' ') 5252 && (c != ',')) 5253 goto unexpected_char_error; 5254 --sp; 5255 while (*sp == '0') 5256 *sp-- = 'z'; 5257 trail = 1; 5258 if (*s == 'z') 5259 goto donchr; 5260 } 5261 5262 /* If enough digits were given to more than fill up the yy register, 5263 continuing until overflow into the high guard word yy[2] 5264 guarantees that there will be a roundoff bit at the top 5265 of the low guard word after normalization. */ 5266 5267 if (yy[2] == 0) 5268 { 5269 if (base == 16) 5270 { 5271 if (decflg) 5272 nexp += 4; /* count digits after decimal point */ 5273 5274 eshup1 (yy); /* multiply current number by 16 */ 5275 eshup1 (yy); 5276 eshup1 (yy); 5277 eshup1 (yy); 5278 } 5279 else 5280 { 5281 if (decflg) 5282 nexp += 1; /* count digits after decimal point */ 5283 5284 eshup1 (yy); /* multiply current number by 10 */ 5285 emovz (yy, xt); 5286 eshup1 (xt); 5287 eshup1 (xt); 5288 eaddm (xt, yy); 5289 } 5290 /* Insert the current digit. */ 5291 ecleaz (xt); 5292 xt[NI - 2] = (UEMUSHORT) k; 5293 eaddm (xt, yy); 5294 } 5295 else 5296 { 5297 /* Mark any lost non-zero digit. */ 5298 lost |= k; 5299 /* Count lost digits before the decimal point. */ 5300 if (decflg == 0) 5301 { 5302 if (base == 10) 5303 nexp -= 1; 5304 else 5305 nexp -= 4; 5306 } 5307 } 5308 prec += 1; 5309 goto donchr; 5310 } 5311 5312 switch (*s) 5313 { 5314 case 'z': 5315 break; 5316 case 'E': 5317 case 'e': 5318 case 'P': 5319 case 'p': 5320 goto expnt; 5321 case '.': /* decimal point */ 5322 if (decflg) 5323 goto unexpected_char_error; 5324 ++decflg; 5325 break; 5326 case '-': 5327 nsign = 0xffff; 5328 if (sgnflg) 5329 goto unexpected_char_error; 5330 ++sgnflg; 5331 break; 5332 case '+': 5333 if (sgnflg) 5334 goto unexpected_char_error; 5335 ++sgnflg; 5336 break; 5337 case ',': 5338 case ' ': 5339 case '\0': 5340 case '\n': 5341 case '\r': 5342 goto daldone; 5343 case 'i': 5344 case 'I': 5345 goto infinite; 5346 default: 5347 unexpected_char_error: 5348#ifdef NANS 5349 einan (yy); 5350#else 5351 mtherr ("asctoe", DOMAIN); 5352 eclear (yy); 5353#endif 5354 goto aexit; 5355 } 5356 donchr: 5357 ++s; 5358 goto nxtcom; 5359 5360 /* Exponent interpretation */ 5361 expnt: 5362 /* 0.0eXXX is zero, regardless of XXX. Check for the 0.0. */ 5363 for (k = 0; k < NI; k++) 5364 { 5365 if (yy[k] != 0) 5366 goto read_expnt; 5367 } 5368 goto aexit; 5369 5370read_expnt: 5371 esign = 1; 5372 exp = 0; 5373 ++s; 5374 /* check for + or - */ 5375 if (*s == '-') 5376 { 5377 esign = -1; 5378 ++s; 5379 } 5380 if (*s == '+') 5381 ++s; 5382 while (ISDIGIT (*s)) 5383 { 5384 exp *= 10; 5385 exp += *s++ - '0'; 5386 if (exp > 999999) 5387 break; 5388 } 5389 if (esign < 0) 5390 exp = -exp; 5391 if ((exp > MAXDECEXP) && (base == 10)) 5392 { 5393 infinite: 5394 ecleaz (yy); 5395 yy[E] = 0x7fff; /* infinity */ 5396 goto aexit; 5397 } 5398 if ((exp < MINDECEXP) && (base == 10)) 5399 { 5400 zero: 5401 ecleaz (yy); 5402 goto aexit; 5403 } 5404 5405 daldone: 5406 if (base == 16) 5407 { 5408 /* Base 16 hexadecimal floating constant. */ 5409 if ((k = enormlz (yy)) > NBITS) 5410 { 5411 ecleaz (yy); 5412 goto aexit; 5413 } 5414 /* Adjust the exponent. NEXP is the number of hex digits, 5415 EXP is a power of 2. */ 5416 lexp = (EXONE - 1 + NBITS) - k + yy[E] + exp - nexp; 5417 if (lexp > 0x7fff) 5418 goto infinite; 5419 if (lexp < 0) 5420 goto zero; 5421 yy[E] = lexp; 5422 goto expdon; 5423 } 5424 5425 nexp = exp - nexp; 5426 /* Pad trailing zeros to minimize power of 10, per IEEE spec. */ 5427 while ((nexp > 0) && (yy[2] == 0)) 5428 { 5429 emovz (yy, xt); 5430 eshup1 (xt); 5431 eshup1 (xt); 5432 eaddm (yy, xt); 5433 eshup1 (xt); 5434 if (xt[2] != 0) 5435 break; 5436 nexp -= 1; 5437 emovz (xt, yy); 5438 } 5439 if ((k = enormlz (yy)) > NBITS) 5440 { 5441 ecleaz (yy); 5442 goto aexit; 5443 } 5444 lexp = (EXONE - 1 + NBITS) - k; 5445 emdnorm (yy, lost, 0, lexp, 64); 5446 lost = 0; 5447 5448 /* Convert to external format: 5449 5450 Multiply by 10**nexp. If precision is 64 bits, 5451 the maximum relative error incurred in forming 10**n 5452 for 0 <= n <= 324 is 8.2e-20, at 10**180. 5453 For 0 <= n <= 999, the peak relative error is 1.4e-19 at 10**947. 5454 For 0 >= n >= -999, it is -1.55e-19 at 10**-435. */ 5455 5456 lexp = yy[E]; 5457 if (nexp == 0) 5458 { 5459 k = 0; 5460 goto expdon; 5461 } 5462 esign = 1; 5463 if (nexp < 0) 5464 { 5465 nexp = -nexp; 5466 esign = -1; 5467 if (nexp > 4096) 5468 { 5469 /* Punt. Can't handle this without 2 divides. */ 5470 emovi (etens[0], tt); 5471 lexp -= tt[E]; 5472 k = edivm (tt, yy); 5473 lexp += EXONE; 5474 nexp -= 4096; 5475 } 5476 } 5477 emov (eone, xt); 5478 exp = 1; 5479 i = NTEN; 5480 do 5481 { 5482 if (exp & nexp) 5483 emul (etens[i], xt, xt); 5484 i--; 5485 exp = exp + exp; 5486 } 5487 while (exp <= MAXP); 5488 5489 emovi (xt, tt); 5490 if (esign < 0) 5491 { 5492 lexp -= tt[E]; 5493 k = edivm (tt, yy); 5494 lexp += EXONE; 5495 } 5496 else 5497 { 5498 lexp += tt[E]; 5499 k = emulm (tt, yy); 5500 lexp -= EXONE - 1; 5501 } 5502 lost = k; 5503 5504 expdon: 5505 5506 /* Round and convert directly to the destination type */ 5507 if (oprec == 53) 5508 lexp -= EXONE - 0x3ff; 5509#ifdef C4X 5510 else if (oprec == 24 || oprec == 32) 5511 lexp -= (EXONE - 0x7f); 5512#else 5513#ifdef IBM 5514 else if (oprec == 24 || oprec == 56) 5515 lexp -= EXONE - (0x41 << 2); 5516#else 5517 else if (oprec == 24) 5518 lexp -= EXONE - 0177; 5519#endif /* IBM */ 5520#endif /* C4X */ 5521#ifdef DEC 5522 else if (oprec == 56) 5523 lexp -= EXONE - 0201; 5524#endif 5525 rndprc = oprec; 5526 emdnorm (yy, lost, 0, lexp, 64); 5527 5528 aexit: 5529 5530 rndprc = rndsav; 5531 yy[0] = nsign; 5532 switch (oprec) 5533 { 5534#ifdef DEC 5535 case 56: 5536 todec (yy, y); /* see etodec.c */ 5537 break; 5538#endif 5539#ifdef IBM 5540 case 56: 5541 toibm (yy, y, DFmode); 5542 break; 5543#endif 5544#ifdef C4X 5545 case 32: 5546 toc4x (yy, y, HFmode); 5547 break; 5548#endif 5549 5550 case 53: 5551 toe53 (yy, y); 5552 break; 5553 case 24: 5554 toe24 (yy, y); 5555 break; 5556 case 64: 5557 toe64 (yy, y); 5558 break; 5559#if (INTEL_EXTENDED_IEEE_FORMAT == 0) 5560 case 113: 5561 toe113 (yy, y); 5562 break; 5563#endif 5564 case NBITS: 5565 emovo (yy, y); 5566 break; 5567 } 5568} 5569 5570 5571 5572/* Return Y = largest integer not greater than X (truncated toward minus 5573 infinity). */ 5574 5575static const UEMUSHORT bmask[] = 5576{ 5577 0xffff, 5578 0xfffe, 5579 0xfffc, 5580 0xfff8, 5581 0xfff0, 5582 0xffe0, 5583 0xffc0, 5584 0xff80, 5585 0xff00, 5586 0xfe00, 5587 0xfc00, 5588 0xf800, 5589 0xf000, 5590 0xe000, 5591 0xc000, 5592 0x8000, 5593 0x0000, 5594}; 5595 5596static void 5597efloor (x, y) 5598 const UEMUSHORT x[]; 5599 UEMUSHORT y[]; 5600{ 5601 UEMUSHORT *p; 5602 int e, expon, i; 5603 UEMUSHORT f[NE]; 5604 5605 emov (x, f); /* leave in external format */ 5606 expon = (int) f[NE - 1]; 5607 e = (expon & 0x7fff) - (EXONE - 1); 5608 if (e <= 0) 5609 { 5610 eclear (y); 5611 goto isitneg; 5612 } 5613 /* number of bits to clear out */ 5614 e = NBITS - e; 5615 emov (f, y); 5616 if (e <= 0) 5617 return; 5618 5619 p = &y[0]; 5620 while (e >= 16) 5621 { 5622 *p++ = 0; 5623 e -= 16; 5624 } 5625 /* clear the remaining bits */ 5626 *p &= bmask[e]; 5627 /* truncate negatives toward minus infinity */ 5628 isitneg: 5629 5630 if ((UEMUSHORT) expon & (UEMUSHORT) 0x8000) 5631 { 5632 for (i = 0; i < NE - 1; i++) 5633 { 5634 if (f[i] != y[i]) 5635 { 5636 esub (eone, y, y); 5637 break; 5638 } 5639 } 5640 } 5641} 5642 5643 5644#if 0 5645/* Return S and EXP such that S * 2^EXP = X and .5 <= S < 1. 5646 For example, 1.1 = 0.55 * 2^1. */ 5647 5648static void 5649efrexp (x, exp, s) 5650 const UEMUSHORT x[]; 5651 int *exp; 5652 UEMUSHORT s[]; 5653{ 5654 UEMUSHORT xi[NI]; 5655 EMULONG li; 5656 5657 emovi (x, xi); 5658 /* Handle denormalized numbers properly using long integer exponent. */ 5659 li = (EMULONG) ((EMUSHORT) xi[1]); 5660 5661 if (li == 0) 5662 { 5663 li -= enormlz (xi); 5664 } 5665 xi[1] = 0x3ffe; 5666 emovo (xi, s); 5667 *exp = (int) (li - 0x3ffe); 5668} 5669#endif 5670 5671/* Return e type Y = X * 2^PWR2. */ 5672 5673static void 5674eldexp (x, pwr2, y) 5675 const UEMUSHORT x[]; 5676 int pwr2; 5677 UEMUSHORT y[]; 5678{ 5679 UEMUSHORT xi[NI]; 5680 EMULONG li; 5681 int i; 5682 5683 emovi (x, xi); 5684 li = xi[1]; 5685 li += pwr2; 5686 i = 0; 5687 emdnorm (xi, i, i, li, 64); 5688 emovo (xi, y); 5689} 5690 5691 5692#if 0 5693/* C = remainder after dividing B by A, all e type values. 5694 Least significant integer quotient bits left in EQUOT. */ 5695 5696static void 5697eremain (a, b, c) 5698 const UEMUSHORT a[], b[]; 5699 UEMUSHORT c[]; 5700{ 5701 UEMUSHORT den[NI], num[NI]; 5702 5703#ifdef NANS 5704 if (eisinf (b) 5705 || (ecmp (a, ezero) == 0) 5706 || eisnan (a) 5707 || eisnan (b)) 5708 { 5709 enan (c, 0); 5710 return; 5711 } 5712#endif 5713 if (ecmp (a, ezero) == 0) 5714 { 5715 mtherr ("eremain", SING); 5716 eclear (c); 5717 return; 5718 } 5719 emovi (a, den); 5720 emovi (b, num); 5721 eiremain (den, num); 5722 /* Sign of remainder = sign of quotient */ 5723 if (a[0] == b[0]) 5724 num[0] = 0; 5725 else 5726 num[0] = 0xffff; 5727 emovo (num, c); 5728} 5729#endif 5730 5731/* Return quotient of exploded e-types NUM / DEN in EQUOT, 5732 remainder in NUM. */ 5733 5734static void 5735eiremain (den, num) 5736 UEMUSHORT den[], num[]; 5737{ 5738 EMULONG ld, ln; 5739 UEMUSHORT j; 5740 5741 ld = den[E]; 5742 ld -= enormlz (den); 5743 ln = num[E]; 5744 ln -= enormlz (num); 5745 ecleaz (equot); 5746 while (ln >= ld) 5747 { 5748 if (ecmpm (den, num) <= 0) 5749 { 5750 esubm (den, num); 5751 j = 1; 5752 } 5753 else 5754 j = 0; 5755 eshup1 (equot); 5756 equot[NI - 1] |= j; 5757 eshup1 (num); 5758 ln -= 1; 5759 } 5760 emdnorm (num, 0, 0, ln, 0); 5761} 5762 5763/* Report an error condition CODE encountered in function NAME. 5764 5765 Mnemonic Value Significance 5766 5767 DOMAIN 1 argument domain error 5768 SING 2 function singularity 5769 OVERFLOW 3 overflow range error 5770 UNDERFLOW 4 underflow range error 5771 TLOSS 5 total loss of precision 5772 PLOSS 6 partial loss of precision 5773 INVALID 7 NaN - producing operation 5774 EDOM 33 Unix domain error code 5775 ERANGE 34 Unix range error code 5776 5777 The order of appearance of the following messages is bound to the 5778 error codes defined above. */ 5779 5780int merror = 0; 5781extern int merror; 5782 5783static void 5784mtherr (name, code) 5785 const char *name; 5786 int code; 5787{ 5788 /* The string passed by the calling program is supposed to be the 5789 name of the function in which the error occurred. 5790 The code argument selects which error message string will be printed. */ 5791 5792 if (strcmp (name, "esub") == 0) 5793 name = "subtraction"; 5794 else if (strcmp (name, "ediv") == 0) 5795 name = "division"; 5796 else if (strcmp (name, "emul") == 0) 5797 name = "multiplication"; 5798 else if (strcmp (name, "enormlz") == 0) 5799 name = "normalization"; 5800 else if (strcmp (name, "etoasc") == 0) 5801 name = "conversion to text"; 5802 else if (strcmp (name, "asctoe") == 0) 5803 name = "parsing"; 5804 else if (strcmp (name, "eremain") == 0) 5805 name = "modulus"; 5806 else if (strcmp (name, "esqrt") == 0) 5807 name = "square root"; 5808 if (extra_warnings) 5809 { 5810 switch (code) 5811 { 5812 case DOMAIN: warning ("%s: argument domain error" , name); break; 5813 case SING: warning ("%s: function singularity" , name); break; 5814 case OVERFLOW: warning ("%s: overflow range error" , name); break; 5815 case UNDERFLOW: warning ("%s: underflow range error" , name); break; 5816 case TLOSS: warning ("%s: total loss of precision" , name); break; 5817 case PLOSS: warning ("%s: partial loss of precision", name); break; 5818 case INVALID: warning ("%s: NaN - producing operation", name); break; 5819 default: abort (); 5820 } 5821 } 5822 5823 /* Set global error message word */ 5824 merror = code + 1; 5825} 5826 5827#ifdef DEC 5828/* Convert DEC double precision D to e type E. */ 5829 5830static void 5831dectoe (d, e) 5832 const UEMUSHORT *d; 5833 UEMUSHORT *e; 5834{ 5835 UEMUSHORT y[NI]; 5836 UEMUSHORT r, *p; 5837 5838 ecleaz (y); /* start with a zero */ 5839 p = y; /* point to our number */ 5840 r = *d; /* get DEC exponent word */ 5841 if (*d & (unsigned int) 0x8000) 5842 *p = 0xffff; /* fill in our sign */ 5843 ++p; /* bump pointer to our exponent word */ 5844 r &= 0x7fff; /* strip the sign bit */ 5845 if (r == 0) /* answer = 0 if high order DEC word = 0 */ 5846 goto done; 5847 5848 5849 r >>= 7; /* shift exponent word down 7 bits */ 5850 r += EXONE - 0201; /* subtract DEC exponent offset */ 5851 /* add our e type exponent offset */ 5852 *p++ = r; /* to form our exponent */ 5853 5854 r = *d++; /* now do the high order mantissa */ 5855 r &= 0177; /* strip off the DEC exponent and sign bits */ 5856 r |= 0200; /* the DEC understood high order mantissa bit */ 5857 *p++ = r; /* put result in our high guard word */ 5858 5859 *p++ = *d++; /* fill in the rest of our mantissa */ 5860 *p++ = *d++; 5861 *p = *d; 5862 5863 eshdn8 (y); /* shift our mantissa down 8 bits */ 5864 done: 5865 emovo (y, e); 5866} 5867 5868/* Convert e type X to DEC double precision D. */ 5869 5870static void 5871etodec (x, d) 5872 const UEMUSHORT *x; 5873 UEMUSHORT *d; 5874{ 5875 UEMUSHORT xi[NI]; 5876 EMULONG exp; 5877 int rndsav; 5878 5879 emovi (x, xi); 5880 /* Adjust exponent for offsets. */ 5881 exp = (EMULONG) xi[E] - (EXONE - 0201); 5882 /* Round off to nearest or even. */ 5883 rndsav = rndprc; 5884 rndprc = 56; 5885 emdnorm (xi, 0, 0, exp, 64); 5886 rndprc = rndsav; 5887 todec (xi, d); 5888} 5889 5890/* Convert exploded e-type X, that has already been rounded to 5891 56-bit precision, to DEC format double Y. */ 5892 5893static void 5894todec (x, y) 5895 UEMUSHORT *x, *y; 5896{ 5897 UEMUSHORT i; 5898 UEMUSHORT *p; 5899 5900 p = x; 5901 *y = 0; 5902 if (*p++) 5903 *y = 0100000; 5904 i = *p++; 5905 if (i == 0) 5906 { 5907 *y++ = 0; 5908 *y++ = 0; 5909 *y++ = 0; 5910 *y++ = 0; 5911 return; 5912 } 5913 if (i > 0377) 5914 { 5915 *y++ |= 077777; 5916 *y++ = 0xffff; 5917 *y++ = 0xffff; 5918 *y++ = 0xffff; 5919#ifdef ERANGE 5920 errno = ERANGE; 5921#endif 5922 return; 5923 } 5924 i &= 0377; 5925 i <<= 7; 5926 eshup8 (x); 5927 x[M] &= 0177; 5928 i |= x[M]; 5929 *y++ |= i; 5930 *y++ = x[M + 1]; 5931 *y++ = x[M + 2]; 5932 *y++ = x[M + 3]; 5933} 5934#endif /* DEC */ 5935 5936#ifdef IBM 5937/* Convert IBM single/double precision to e type. */ 5938 5939static void 5940ibmtoe (d, e, mode) 5941 const UEMUSHORT *d; 5942 UEMUSHORT *e; 5943 enum machine_mode mode; 5944{ 5945 UEMUSHORT y[NI]; 5946 UEMUSHORT r, *p; 5947 5948 ecleaz (y); /* start with a zero */ 5949 p = y; /* point to our number */ 5950 r = *d; /* get IBM exponent word */ 5951 if (*d & (unsigned int) 0x8000) 5952 *p = 0xffff; /* fill in our sign */ 5953 ++p; /* bump pointer to our exponent word */ 5954 r &= 0x7f00; /* strip the sign bit */ 5955 r >>= 6; /* shift exponent word down 6 bits */ 5956 /* in fact shift by 8 right and 2 left */ 5957 r += EXONE - (0x41 << 2); /* subtract IBM exponent offset */ 5958 /* add our e type exponent offset */ 5959 *p++ = r; /* to form our exponent */ 5960 5961 *p++ = *d++ & 0xff; /* now do the high order mantissa */ 5962 /* strip off the IBM exponent and sign bits */ 5963 if (mode != SFmode) /* there are only 2 words in SFmode */ 5964 { 5965 *p++ = *d++; /* fill in the rest of our mantissa */ 5966 *p++ = *d++; 5967 } 5968 *p = *d; 5969 5970 if (y[M] == 0 && y[M+1] == 0 && y[M+2] == 0 && y[M+3] == 0) 5971 y[0] = y[E] = 0; 5972 else 5973 y[E] -= 5 + enormlz (y); /* now normalise the mantissa */ 5974 /* handle change in RADIX */ 5975 emovo (y, e); 5976} 5977 5978 5979 5980/* Convert e type to IBM single/double precision. */ 5981 5982static void 5983etoibm (x, d, mode) 5984 const UEMUSHORT *x; 5985 UEMUSHORT *d; 5986 enum machine_mode mode; 5987{ 5988 UEMUSHORT xi[NI]; 5989 EMULONG exp; 5990 int rndsav; 5991 5992 emovi (x, xi); 5993 exp = (EMULONG) xi[E] - (EXONE - (0x41 << 2)); /* adjust exponent for offsets */ 5994 /* round off to nearest or even */ 5995 rndsav = rndprc; 5996 rndprc = 56; 5997 emdnorm (xi, 0, 0, exp, 64); 5998 rndprc = rndsav; 5999 toibm (xi, d, mode); 6000} 6001 6002static void 6003toibm (x, y, mode) 6004 UEMUSHORT *x, *y; 6005 enum machine_mode mode; 6006{ 6007 UEMUSHORT i; 6008 UEMUSHORT *p; 6009 int r; 6010 6011 p = x; 6012 *y = 0; 6013 if (*p++) 6014 *y = 0x8000; 6015 i = *p++; 6016 if (i == 0) 6017 { 6018 *y++ = 0; 6019 *y++ = 0; 6020 if (mode != SFmode) 6021 { 6022 *y++ = 0; 6023 *y++ = 0; 6024 } 6025 return; 6026 } 6027 r = i & 0x3; 6028 i >>= 2; 6029 if (i > 0x7f) 6030 { 6031 *y++ |= 0x7fff; 6032 *y++ = 0xffff; 6033 if (mode != SFmode) 6034 { 6035 *y++ = 0xffff; 6036 *y++ = 0xffff; 6037 } 6038#ifdef ERANGE 6039 errno = ERANGE; 6040#endif 6041 return; 6042 } 6043 i &= 0x7f; 6044 *y |= (i << 8); 6045 eshift (x, r + 5); 6046 *y++ |= x[M]; 6047 *y++ = x[M + 1]; 6048 if (mode != SFmode) 6049 { 6050 *y++ = x[M + 2]; 6051 *y++ = x[M + 3]; 6052 } 6053} 6054#endif /* IBM */ 6055 6056 6057#ifdef C4X 6058/* Convert C4X single/double precision to e type. */ 6059 6060static void 6061c4xtoe (d, e, mode) 6062 const UEMUSHORT *d; 6063 UEMUSHORT *e; 6064 enum machine_mode mode; 6065{ 6066 UEMUSHORT y[NI]; 6067 UEMUSHORT dn[4]; 6068 int r; 6069 int isnegative; 6070 int size; 6071 int i; 6072 int carry; 6073 6074 dn[0] = d[0]; 6075 dn[1] = d[1]; 6076 if (mode != QFmode) 6077 { 6078 dn[2] = d[3] << 8; 6079 dn[3] = 0; 6080 } 6081 6082 /* Short-circuit the zero case. */ 6083 if ((dn[0] == 0x8000) 6084 && (dn[1] == 0x0000) 6085 && ((mode == QFmode) || ((dn[2] == 0x0000) && (dn[3] == 0x0000)))) 6086 { 6087 e[0] = 0; 6088 e[1] = 0; 6089 e[2] = 0; 6090 e[3] = 0; 6091 e[4] = 0; 6092 e[5] = 0; 6093 return; 6094 } 6095 6096 ecleaz (y); /* start with a zero */ 6097 r = dn[0]; /* get sign/exponent part */ 6098 if (r & (unsigned int) 0x0080) 6099 { 6100 y[0] = 0xffff; /* fill in our sign */ 6101 isnegative = TRUE; 6102 } 6103 else 6104 isnegative = FALSE; 6105 6106 r >>= 8; /* Shift exponent word down 8 bits. */ 6107 if (r & 0x80) /* Make the exponent negative if it is. */ 6108 r = r | (~0 & ~0xff); 6109 6110 if (isnegative) 6111 { 6112 /* Now do the high order mantissa. We don't "or" on the high bit 6113 because it is 2 (not 1) and is handled a little differently 6114 below. */ 6115 y[M] = dn[0] & 0x7f; 6116 6117 y[M+1] = dn[1]; 6118 if (mode != QFmode) /* There are only 2 words in QFmode. */ 6119 { 6120 y[M+2] = dn[2]; /* Fill in the rest of our mantissa. */ 6121 y[M+3] = dn[3]; 6122 size = 4; 6123 } 6124 else 6125 size = 2; 6126 eshift (y, -8); 6127 6128 /* Now do the two's complement on the data. */ 6129 6130 carry = 1; /* Initially add 1 for the two's complement. */ 6131 for (i=size + M; i > M; i--) 6132 { 6133 if (carry && (y[i] == 0x0000)) 6134 /* We overflowed into the next word, carry is the same. */ 6135 y[i] = carry ? 0x0000 : 0xffff; 6136 else 6137 { 6138 /* No overflow, just invert and add carry. */ 6139 y[i] = ((~y[i]) + carry) & 0xffff; 6140 carry = 0; 6141 } 6142 } 6143 6144 if (carry) 6145 { 6146 eshift (y, -1); 6147 y[M+1] |= 0x8000; 6148 r++; 6149 } 6150 y[1] = r + EXONE; 6151 } 6152 else 6153 { 6154 /* Add our e type exponent offset to form our exponent. */ 6155 r += EXONE; 6156 y[1] = r; 6157 6158 /* Now do the high order mantissa strip off the exponent and sign 6159 bits and add the high 1 bit. */ 6160 y[M] = (dn[0] & 0x7f) | 0x80; 6161 6162 y[M+1] = dn[1]; 6163 if (mode != QFmode) /* There are only 2 words in QFmode. */ 6164 { 6165 y[M+2] = dn[2]; /* Fill in the rest of our mantissa. */ 6166 y[M+3] = dn[3]; 6167 } 6168 eshift (y, -8); 6169 } 6170 6171 emovo (y, e); 6172} 6173 6174 6175/* Convert e type to C4X single/double precision. */ 6176 6177static void 6178etoc4x (x, d, mode) 6179 const UEMUSHORT *x; 6180 UEMUSHORT *d; 6181 enum machine_mode mode; 6182{ 6183 UEMUSHORT xi[NI]; 6184 EMULONG exp; 6185 int rndsav; 6186 6187 emovi (x, xi); 6188 6189 /* Adjust exponent for offsets. */ 6190 exp = (EMULONG) xi[E] - (EXONE - 0x7f); 6191 6192 /* Round off to nearest or even. */ 6193 rndsav = rndprc; 6194 rndprc = mode == QFmode ? 24 : 32; 6195 emdnorm (xi, 0, 0, exp, 64); 6196 rndprc = rndsav; 6197 toc4x (xi, d, mode); 6198} 6199 6200static void 6201toc4x (x, y, mode) 6202 UEMUSHORT *x, *y; 6203 enum machine_mode mode; 6204{ 6205 int i; 6206 int v; 6207 int carry; 6208 6209 /* Short-circuit the zero case */ 6210 if ((x[0] == 0) /* Zero exponent and sign */ 6211 && (x[1] == 0) 6212 && (x[M] == 0) /* The rest is for zero mantissa */ 6213 && (x[M+1] == 0) 6214 /* Only check for double if necessary */ 6215 && ((mode == QFmode) || ((x[M+2] == 0) && (x[M+3] == 0)))) 6216 { 6217 /* We have a zero. Put it into the output and return. */ 6218 *y++ = 0x8000; 6219 *y++ = 0x0000; 6220 if (mode != QFmode) 6221 { 6222 *y++ = 0x0000; 6223 *y++ = 0x0000; 6224 } 6225 return; 6226 } 6227 6228 *y = 0; 6229 6230 /* Negative number require a two's complement conversion of the 6231 mantissa. */ 6232 if (x[0]) 6233 { 6234 *y = 0x0080; 6235 6236 i = ((int) x[1]) - 0x7f; 6237 6238 /* Now add 1 to the inverted data to do the two's complement. */ 6239 if (mode != QFmode) 6240 v = 4 + M; 6241 else 6242 v = 2 + M; 6243 carry = 1; 6244 while (v > M) 6245 { 6246 if (x[v] == 0x0000) 6247 x[v] = carry ? 0x0000 : 0xffff; 6248 else 6249 { 6250 x[v] = ((~x[v]) + carry) & 0xffff; 6251 carry = 0; 6252 } 6253 v--; 6254 } 6255 6256 /* The following is a special case. The C4X negative float requires 6257 a zero in the high bit (because the format is (2 - x) x 2^m), so 6258 if a one is in that bit, we have to shift left one to get rid 6259 of it. This only occurs if the number is -1 x 2^m. */ 6260 if (x[M+1] & 0x8000) 6261 { 6262 /* This is the case of -1 x 2^m, we have to rid ourselves of the 6263 high sign bit and shift the exponent. */ 6264 eshift (x, 1); 6265 i--; 6266 } 6267 } 6268 else 6269 i = ((int) x[1]) - 0x7f; 6270 6271 if ((i < -128) || (i > 127)) 6272 { 6273 y[0] |= 0xff7f; 6274 y[1] = 0xffff; 6275 if (mode != QFmode) 6276 { 6277 y[2] = 0xffff; 6278 y[3] = 0xffff; 6279 y[3] = (y[1] << 8) | ((y[2] >> 8) & 0xff); 6280 y[2] = (y[0] << 8) | ((y[1] >> 8) & 0xff); 6281 } 6282#ifdef ERANGE 6283 errno = ERANGE; 6284#endif 6285 return; 6286 } 6287 6288 y[0] |= ((i & 0xff) << 8); 6289 6290 eshift (x, 8); 6291 6292 y[0] |= x[M] & 0x7f; 6293 y[1] = x[M + 1]; 6294 if (mode != QFmode) 6295 { 6296 y[2] = x[M + 2]; 6297 y[3] = x[M + 3]; 6298 y[3] = (y[1] << 8) | ((y[2] >> 8) & 0xff); 6299 y[2] = (y[0] << 8) | ((y[1] >> 8) & 0xff); 6300 } 6301} 6302#endif /* C4X */ 6303 6304/* Output a binary NaN bit pattern in the target machine's format. */ 6305 6306/* If special NaN bit patterns are required, define them in tm.h 6307 as arrays of unsigned 16-bit shorts. Otherwise, use the default 6308 patterns here. */ 6309#ifdef TFMODE_NAN 6310TFMODE_NAN; 6311#else 6312#ifdef IEEE 6313static const UEMUSHORT TFbignan[8] = 6314 {0x7fff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff}; 6315static const UEMUSHORT TFlittlenan[8] = {0, 0, 0, 0, 0, 0, 0x8000, 0xffff}; 6316#endif 6317#endif 6318 6319#ifdef XFMODE_NAN 6320XFMODE_NAN; 6321#else 6322#ifdef IEEE 6323static const UEMUSHORT XFbignan[6] = 6324 {0x7fff, 0xffff, 0xffff, 0xffff, 0xffff, 0xffff}; 6325static const UEMUSHORT XFlittlenan[6] = {0, 0, 0, 0xc000, 0xffff, 0}; 6326#endif 6327#endif 6328 6329#ifdef DFMODE_NAN 6330DFMODE_NAN; 6331#else 6332#ifdef IEEE 6333static const UEMUSHORT DFbignan[4] = {0x7fff, 0xffff, 0xffff, 0xffff}; 6334static const UEMUSHORT DFlittlenan[4] = {0, 0, 0, 0xfff8}; 6335#endif 6336#endif 6337 6338#ifdef SFMODE_NAN 6339SFMODE_NAN; 6340#else 6341#ifdef IEEE 6342static const UEMUSHORT SFbignan[2] = {0x7fff, 0xffff}; 6343static const UEMUSHORT SFlittlenan[2] = {0, 0xffc0}; 6344#endif 6345#endif 6346 6347 6348#ifdef NANS 6349static void 6350make_nan (nan, sign, mode) 6351 UEMUSHORT *nan; 6352 int sign; 6353 enum machine_mode mode; 6354{ 6355 int n; 6356 const UEMUSHORT *p; 6357 6358 switch (mode) 6359 { 6360/* Possibly the `reserved operand' patterns on a VAX can be 6361 used like NaN's, but probably not in the same way as IEEE. */ 6362#if !defined(DEC) && !defined(IBM) && !defined(C4X) 6363 case TFmode: 6364#if (INTEL_EXTENDED_IEEE_FORMAT == 0) 6365 n = 8; 6366 if (REAL_WORDS_BIG_ENDIAN) 6367 p = TFbignan; 6368 else 6369 p = TFlittlenan; 6370 break; 6371#endif 6372 /* FALLTHRU */ 6373 6374 case XFmode: 6375 n = 6; 6376 if (REAL_WORDS_BIG_ENDIAN) 6377 p = XFbignan; 6378 else 6379 p = XFlittlenan; 6380 break; 6381 6382 case DFmode: 6383 n = 4; 6384 if (REAL_WORDS_BIG_ENDIAN) 6385 p = DFbignan; 6386 else 6387 p = DFlittlenan; 6388 break; 6389 6390 case SFmode: 6391 case HFmode: 6392 n = 2; 6393 if (REAL_WORDS_BIG_ENDIAN) 6394 p = SFbignan; 6395 else 6396 p = SFlittlenan; 6397 break; 6398#endif 6399 6400 default: 6401 abort (); 6402 } 6403 if (REAL_WORDS_BIG_ENDIAN) 6404 *nan++ = (sign << 15) | (*p++ & 0x7fff); 6405 while (--n != 0) 6406 *nan++ = *p++; 6407 if (! REAL_WORDS_BIG_ENDIAN) 6408 *nan = (sign << 15) | (*p & 0x7fff); 6409} 6410#endif /* NANS */ 6411 6412/* This is the inverse of the function `etarsingle' invoked by 6413 REAL_VALUE_TO_TARGET_SINGLE. */ 6414 6415REAL_VALUE_TYPE 6416ereal_unto_float (f) 6417 long f; 6418{ 6419 REAL_VALUE_TYPE r; 6420 UEMUSHORT s[2]; 6421 UEMUSHORT e[NE]; 6422 6423 /* Convert 32 bit integer to array of 16 bit pieces in target machine order. 6424 This is the inverse operation to what the function `endian' does. */ 6425 if (REAL_WORDS_BIG_ENDIAN) 6426 { 6427 s[0] = (UEMUSHORT) (f >> 16); 6428 s[1] = (UEMUSHORT) f; 6429 } 6430 else 6431 { 6432 s[0] = (UEMUSHORT) f; 6433 s[1] = (UEMUSHORT) (f >> 16); 6434 } 6435 /* Convert and promote the target float to E-type. */ 6436 e24toe (s, e); 6437 /* Output E-type to REAL_VALUE_TYPE. */ 6438 PUT_REAL (e, &r); 6439 return r; 6440} 6441 6442 6443/* This is the inverse of the function `etardouble' invoked by 6444 REAL_VALUE_TO_TARGET_DOUBLE. */ 6445 6446REAL_VALUE_TYPE 6447ereal_unto_double (d) 6448 long d[]; 6449{ 6450 REAL_VALUE_TYPE r; 6451 UEMUSHORT s[4]; 6452 UEMUSHORT e[NE]; 6453 6454 /* Convert array of HOST_WIDE_INT to equivalent array of 16-bit pieces. */ 6455 if (REAL_WORDS_BIG_ENDIAN) 6456 { 6457 s[0] = (UEMUSHORT) (d[0] >> 16); 6458 s[1] = (UEMUSHORT) d[0]; 6459 s[2] = (UEMUSHORT) (d[1] >> 16); 6460 s[3] = (UEMUSHORT) d[1]; 6461 } 6462 else 6463 { 6464 /* Target float words are little-endian. */ 6465 s[0] = (UEMUSHORT) d[0]; 6466 s[1] = (UEMUSHORT) (d[0] >> 16); 6467 s[2] = (UEMUSHORT) d[1]; 6468 s[3] = (UEMUSHORT) (d[1] >> 16); 6469 } 6470 /* Convert target double to E-type. */ 6471 e53toe (s, e); 6472 /* Output E-type to REAL_VALUE_TYPE. */ 6473 PUT_REAL (e, &r); 6474 return r; 6475} 6476 6477 6478/* Convert an SFmode target `float' value to a REAL_VALUE_TYPE. 6479 This is somewhat like ereal_unto_float, but the input types 6480 for these are different. */ 6481 6482REAL_VALUE_TYPE 6483ereal_from_float (f) 6484 HOST_WIDE_INT f; 6485{ 6486 REAL_VALUE_TYPE r; 6487 UEMUSHORT s[2]; 6488 UEMUSHORT e[NE]; 6489 6490 /* Convert 32 bit integer to array of 16 bit pieces in target machine order. 6491 This is the inverse operation to what the function `endian' does. */ 6492 if (REAL_WORDS_BIG_ENDIAN) 6493 { 6494 s[0] = (UEMUSHORT) (f >> 16); 6495 s[1] = (UEMUSHORT) f; 6496 } 6497 else 6498 { 6499 s[0] = (UEMUSHORT) f; 6500 s[1] = (UEMUSHORT) (f >> 16); 6501 } 6502 /* Convert and promote the target float to E-type. */ 6503 e24toe (s, e); 6504 /* Output E-type to REAL_VALUE_TYPE. */ 6505 PUT_REAL (e, &r); 6506 return r; 6507} 6508 6509 6510/* Convert a DFmode target `double' value to a REAL_VALUE_TYPE. 6511 This is somewhat like ereal_unto_double, but the input types 6512 for these are different. 6513 6514 The DFmode is stored as an array of HOST_WIDE_INT in the target's 6515 data format, with no holes in the bit packing. The first element 6516 of the input array holds the bits that would come first in the 6517 target computer's memory. */ 6518 6519REAL_VALUE_TYPE 6520ereal_from_double (d) 6521 HOST_WIDE_INT d[]; 6522{ 6523 REAL_VALUE_TYPE r; 6524 UEMUSHORT s[4]; 6525 UEMUSHORT e[NE]; 6526 6527 /* Convert array of HOST_WIDE_INT to equivalent array of 16-bit pieces. */ 6528 if (REAL_WORDS_BIG_ENDIAN) 6529 { 6530#if HOST_BITS_PER_WIDE_INT == 32 6531 s[0] = (UEMUSHORT) (d[0] >> 16); 6532 s[1] = (UEMUSHORT) d[0]; 6533 s[2] = (UEMUSHORT) (d[1] >> 16); 6534 s[3] = (UEMUSHORT) d[1]; 6535#else 6536 /* In this case the entire target double is contained in the 6537 first array element. The second element of the input is 6538 ignored. */ 6539 s[0] = (UEMUSHORT) (d[0] >> 48); 6540 s[1] = (UEMUSHORT) (d[0] >> 32); 6541 s[2] = (UEMUSHORT) (d[0] >> 16); 6542 s[3] = (UEMUSHORT) d[0]; 6543#endif 6544 } 6545 else 6546 { 6547 /* Target float words are little-endian. */ 6548 s[0] = (UEMUSHORT) d[0]; 6549 s[1] = (UEMUSHORT) (d[0] >> 16); 6550#if HOST_BITS_PER_WIDE_INT == 32 6551 s[2] = (UEMUSHORT) d[1]; 6552 s[3] = (UEMUSHORT) (d[1] >> 16); 6553#else 6554 s[2] = (UEMUSHORT) (d[0] >> 32); 6555 s[3] = (UEMUSHORT) (d[0] >> 48); 6556#endif 6557 } 6558 /* Convert target double to E-type. */ 6559 e53toe (s, e); 6560 /* Output E-type to REAL_VALUE_TYPE. */ 6561 PUT_REAL (e, &r); 6562 return r; 6563} 6564 6565 6566#if 0 6567/* Convert target computer unsigned 64-bit integer to e-type. 6568 The endian-ness of DImode follows the convention for integers, 6569 so we use WORDS_BIG_ENDIAN here, not REAL_WORDS_BIG_ENDIAN. */ 6570 6571static void 6572uditoe (di, e) 6573 const UEMUSHORT *di; /* Address of the 64-bit int. */ 6574 UEMUSHORT *e; 6575{ 6576 UEMUSHORT yi[NI]; 6577 int k; 6578 6579 ecleaz (yi); 6580 if (WORDS_BIG_ENDIAN) 6581 { 6582 for (k = M; k < M + 4; k++) 6583 yi[k] = *di++; 6584 } 6585 else 6586 { 6587 for (k = M + 3; k >= M; k--) 6588 yi[k] = *di++; 6589 } 6590 yi[E] = EXONE + 47; /* exponent if normalize shift count were 0 */ 6591 if ((k = enormlz (yi)) > NBITS)/* normalize the significand */ 6592 ecleaz (yi); /* it was zero */ 6593 else 6594 yi[E] -= (UEMUSHORT) k;/* subtract shift count from exponent */ 6595 emovo (yi, e); 6596} 6597 6598/* Convert target computer signed 64-bit integer to e-type. */ 6599 6600static void 6601ditoe (di, e) 6602 const UEMUSHORT *di; /* Address of the 64-bit int. */ 6603 UEMUSHORT *e; 6604{ 6605 unsigned EMULONG acc; 6606 UEMUSHORT yi[NI]; 6607 UEMUSHORT carry; 6608 int k, sign; 6609 6610 ecleaz (yi); 6611 if (WORDS_BIG_ENDIAN) 6612 { 6613 for (k = M; k < M + 4; k++) 6614 yi[k] = *di++; 6615 } 6616 else 6617 { 6618 for (k = M + 3; k >= M; k--) 6619 yi[k] = *di++; 6620 } 6621 /* Take absolute value */ 6622 sign = 0; 6623 if (yi[M] & 0x8000) 6624 { 6625 sign = 1; 6626 carry = 0; 6627 for (k = M + 3; k >= M; k--) 6628 { 6629 acc = (unsigned EMULONG) (~yi[k] & 0xffff) + carry; 6630 yi[k] = acc; 6631 carry = 0; 6632 if (acc & 0x10000) 6633 carry = 1; 6634 } 6635 } 6636 yi[E] = EXONE + 47; /* exponent if normalize shift count were 0 */ 6637 if ((k = enormlz (yi)) > NBITS)/* normalize the significand */ 6638 ecleaz (yi); /* it was zero */ 6639 else 6640 yi[E] -= (UEMUSHORT) k;/* subtract shift count from exponent */ 6641 emovo (yi, e); 6642 if (sign) 6643 eneg (e); 6644} 6645 6646 6647/* Convert e-type to unsigned 64-bit int. */ 6648 6649static void 6650etoudi (x, i) 6651 const UEMUSHORT *x; 6652 UEMUSHORT *i; 6653{ 6654 UEMUSHORT xi[NI]; 6655 int j, k; 6656 6657 emovi (x, xi); 6658 if (xi[0]) 6659 { 6660 xi[M] = 0; 6661 goto noshift; 6662 } 6663 k = (int) xi[E] - (EXONE - 1); 6664 if (k <= 0) 6665 { 6666 for (j = 0; j < 4; j++) 6667 *i++ = 0; 6668 return; 6669 } 6670 if (k > 64) 6671 { 6672 for (j = 0; j < 4; j++) 6673 *i++ = 0xffff; 6674 if (extra_warnings) 6675 warning ("overflow on truncation to integer"); 6676 return; 6677 } 6678 if (k > 16) 6679 { 6680 /* Shift more than 16 bits: first shift up k-16 mod 16, 6681 then shift up by 16's. */ 6682 j = k - ((k >> 4) << 4); 6683 if (j == 0) 6684 j = 16; 6685 eshift (xi, j); 6686 if (WORDS_BIG_ENDIAN) 6687 *i++ = xi[M]; 6688 else 6689 { 6690 i += 3; 6691 *i-- = xi[M]; 6692 } 6693 k -= j; 6694 do 6695 { 6696 eshup6 (xi); 6697 if (WORDS_BIG_ENDIAN) 6698 *i++ = xi[M]; 6699 else 6700 *i-- = xi[M]; 6701 } 6702 while ((k -= 16) > 0); 6703 } 6704 else 6705 { 6706 /* shift not more than 16 bits */ 6707 eshift (xi, k); 6708 6709noshift: 6710 6711 if (WORDS_BIG_ENDIAN) 6712 { 6713 i += 3; 6714 *i-- = xi[M]; 6715 *i-- = 0; 6716 *i-- = 0; 6717 *i = 0; 6718 } 6719 else 6720 { 6721 *i++ = xi[M]; 6722 *i++ = 0; 6723 *i++ = 0; 6724 *i = 0; 6725 } 6726 } 6727} 6728 6729 6730/* Convert e-type to signed 64-bit int. */ 6731 6732static void 6733etodi (x, i) 6734 const UEMUSHORT *x; 6735 UEMUSHORT *i; 6736{ 6737 unsigned EMULONG acc; 6738 UEMUSHORT xi[NI]; 6739 UEMUSHORT carry; 6740 UEMUSHORT *isave; 6741 int j, k; 6742 6743 emovi (x, xi); 6744 k = (int) xi[E] - (EXONE - 1); 6745 if (k <= 0) 6746 { 6747 for (j = 0; j < 4; j++) 6748 *i++ = 0; 6749 return; 6750 } 6751 if (k > 64) 6752 { 6753 for (j = 0; j < 4; j++) 6754 *i++ = 0xffff; 6755 if (extra_warnings) 6756 warning ("overflow on truncation to integer"); 6757 return; 6758 } 6759 isave = i; 6760 if (k > 16) 6761 { 6762 /* Shift more than 16 bits: first shift up k-16 mod 16, 6763 then shift up by 16's. */ 6764 j = k - ((k >> 4) << 4); 6765 if (j == 0) 6766 j = 16; 6767 eshift (xi, j); 6768 if (WORDS_BIG_ENDIAN) 6769 *i++ = xi[M]; 6770 else 6771 { 6772 i += 3; 6773 *i-- = xi[M]; 6774 } 6775 k -= j; 6776 do 6777 { 6778 eshup6 (xi); 6779 if (WORDS_BIG_ENDIAN) 6780 *i++ = xi[M]; 6781 else 6782 *i-- = xi[M]; 6783 } 6784 while ((k -= 16) > 0); 6785 } 6786 else 6787 { 6788 /* shift not more than 16 bits */ 6789 eshift (xi, k); 6790 6791 if (WORDS_BIG_ENDIAN) 6792 { 6793 i += 3; 6794 *i = xi[M]; 6795 *i-- = 0; 6796 *i-- = 0; 6797 *i = 0; 6798 } 6799 else 6800 { 6801 *i++ = xi[M]; 6802 *i++ = 0; 6803 *i++ = 0; 6804 *i = 0; 6805 } 6806 } 6807 /* Negate if negative */ 6808 if (xi[0]) 6809 { 6810 carry = 0; 6811 if (WORDS_BIG_ENDIAN) 6812 isave += 3; 6813 for (k = 0; k < 4; k++) 6814 { 6815 acc = (unsigned EMULONG) (~(*isave) & 0xffff) + carry; 6816 if (WORDS_BIG_ENDIAN) 6817 *isave-- = acc; 6818 else 6819 *isave++ = acc; 6820 carry = 0; 6821 if (acc & 0x10000) 6822 carry = 1; 6823 } 6824 } 6825} 6826 6827 6828/* Longhand square root routine. */ 6829 6830 6831static int esqinited = 0; 6832static unsigned short sqrndbit[NI]; 6833 6834static void 6835esqrt (x, y) 6836 const UEMUSHORT *x; 6837 UEMUSHORT *y; 6838{ 6839 UEMUSHORT temp[NI], num[NI], sq[NI], xx[NI]; 6840 EMULONG m, exp; 6841 int i, j, k, n, nlups; 6842 6843 if (esqinited == 0) 6844 { 6845 ecleaz (sqrndbit); 6846 sqrndbit[NI - 2] = 1; 6847 esqinited = 1; 6848 } 6849 /* Check for arg <= 0 */ 6850 i = ecmp (x, ezero); 6851 if (i <= 0) 6852 { 6853 if (i == -1) 6854 { 6855 mtherr ("esqrt", DOMAIN); 6856 eclear (y); 6857 } 6858 else 6859 emov (x, y); 6860 return; 6861 } 6862 6863#ifdef INFINITY 6864 if (eisinf (x)) 6865 { 6866 eclear (y); 6867 einfin (y); 6868 return; 6869 } 6870#endif 6871 /* Bring in the arg and renormalize if it is denormal. */ 6872 emovi (x, xx); 6873 m = (EMULONG) xx[1]; /* local long word exponent */ 6874 if (m == 0) 6875 m -= enormlz (xx); 6876 6877 /* Divide exponent by 2 */ 6878 m -= 0x3ffe; 6879 exp = (unsigned short) ((m / 2) + 0x3ffe); 6880 6881 /* Adjust if exponent odd */ 6882 if ((m & 1) != 0) 6883 { 6884 if (m > 0) 6885 exp += 1; 6886 eshdn1 (xx); 6887 } 6888 6889 ecleaz (sq); 6890 ecleaz (num); 6891 n = 8; /* get 8 bits of result per inner loop */ 6892 nlups = rndprc; 6893 j = 0; 6894 6895 while (nlups > 0) 6896 { 6897 /* bring in next word of arg */ 6898 if (j < NE) 6899 num[NI - 1] = xx[j + 3]; 6900 /* Do additional bit on last outer loop, for roundoff. */ 6901 if (nlups <= 8) 6902 n = nlups + 1; 6903 for (i = 0; i < n; i++) 6904 { 6905 /* Next 2 bits of arg */ 6906 eshup1 (num); 6907 eshup1 (num); 6908 /* Shift up answer */ 6909 eshup1 (sq); 6910 /* Make trial divisor */ 6911 for (k = 0; k < NI; k++) 6912 temp[k] = sq[k]; 6913 eshup1 (temp); 6914 eaddm (sqrndbit, temp); 6915 /* Subtract and insert answer bit if it goes in */ 6916 if (ecmpm (temp, num) <= 0) 6917 { 6918 esubm (temp, num); 6919 sq[NI - 2] |= 1; 6920 } 6921 } 6922 nlups -= n; 6923 j += 1; 6924 } 6925 6926 /* Adjust for extra, roundoff loop done. */ 6927 exp += (NBITS - 1) - rndprc; 6928 6929 /* Sticky bit = 1 if the remainder is nonzero. */ 6930 k = 0; 6931 for (i = 3; i < NI; i++) 6932 k |= (int) num[i]; 6933 6934 /* Renormalize and round off. */ 6935 emdnorm (sq, k, 0, exp, 64); 6936 emovo (sq, y); 6937} 6938#endif 6939#endif /* EMU_NON_COMPILE not defined */ 6940 6941/* Return the binary precision of the significand for a given 6942 floating point mode. The mode can hold an integer value 6943 that many bits wide, without losing any bits. */ 6944 6945unsigned int 6946significand_size (mode) 6947 enum machine_mode mode; 6948{ 6949 6950/* Don't test the modes, but their sizes, lest this 6951 code won't work for BITS_PER_UNIT != 8 . */ 6952 6953switch (GET_MODE_BITSIZE (mode)) 6954 { 6955 case 32: 6956 6957#if TARGET_FLOAT_FORMAT == C4X_FLOAT_FORMAT 6958 return 56; 6959#endif 6960 6961 return 24; 6962 6963 case 64: 6964#if TARGET_FLOAT_FORMAT == IEEE_FLOAT_FORMAT 6965 return 53; 6966#else 6967#if TARGET_FLOAT_FORMAT == IBM_FLOAT_FORMAT 6968 return 56; 6969#else 6970#if TARGET_FLOAT_FORMAT == VAX_FLOAT_FORMAT 6971 return 56; 6972#else 6973#if TARGET_FLOAT_FORMAT == C4X_FLOAT_FORMAT 6974 return 56; 6975#else 6976 abort (); 6977#endif 6978#endif 6979#endif 6980#endif 6981 6982 case 96: 6983 return 64; 6984 6985 case 128: 6986#if (INTEL_EXTENDED_IEEE_FORMAT == 0) 6987 return 113; 6988#else 6989 return 64; 6990#endif 6991 6992 default: 6993 abort (); 6994 } 6995} 6996