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