1/*
2 * CDDL HEADER START
3 *
4 * The contents of this file are subject to the terms of the
5 * Common Development and Distribution License (the "License").
6 * You may not use this file except in compliance with the License.
7 *
8 * You can obtain a copy of the license at usr/src/OPENSOLARIS.LICENSE
9 * or http://www.opensolaris.org/os/licensing.
10 * See the License for the specific language governing permissions
11 * and limitations under the License.
12 *
13 * When distributing Covered Code, include this CDDL HEADER in each
14 * file and include the License file at usr/src/OPENSOLARIS.LICENSE.
15 * If applicable, add the following below this CDDL HEADER, with the
16 * fields enclosed by brackets "[]" replaced with your own identifying
17 * information: Portions Copyright [yyyy] [name of copyright owner]
18 *
19 * CDDL HEADER END
20 */
21
22/*
23 * Copyright 2008 Sun Microsystems, Inc.  All rights reserved.
24 * Use is subject to license terms.
25 */
26
27#pragma ident	"%Z%%M%	%I%	%E% SMI"
28
29/*
30 * Conversion from decimal to binary floating point
31 */
32
33#include "lint.h"
34#include <stdlib.h>
35#include "base_conversion.h"
36
37/*
38 * Convert the integer part of a nonzero base-10^4 _big_float *pd
39 * to base 2^16 in **ppb.  The converted value is accurate to nsig
40 * significant bits.  On exit, *sticky is nonzero if *pd had a
41 * nonzero fractional part.  If pd->exponent > 0 and **ppb is not
42 * large enough to hold the final converted value (i.e., the con-
43 * verted significand scaled by 10^pd->exponent), then on exit,
44 * *ppb will point to a newly allocated _big_float, which must be
45 * freed by the caller.  (The number of significant bits we need
46 * should fit in pb, but __big_float_times_power may allocate new
47 * storage anyway because the exact product could require more than
48 * 16000 bits.)
49 *
50 * This routine does not check that **ppb is large enough to hold
51 * the result of converting the significand of *pd.
52 */
53static void
54__big_decimal_to_big_binary(_big_float *pd, int nsig, _big_float **ppb,
55    int *sticky)
56{
57	_big_float	*pb;
58	int		i, j, len, s;
59	unsigned int	carry;
60
61	pb = *ppb;
62
63	/* convert pd a digit at a time, most significant first */
64	if (pd->bexponent + ((pd->blength - 1) << 2) >= 0) {
65		pb->bsignificand[0] = pd->bsignificand[pd->blength - 1];
66		len = 1;
67		for (i = pd->blength - 2; i >= 0 &&
68		    pd->bexponent + (i << 2) >= 0; i--) {
69			/* multiply pb by 10^4 and add next digit */
70			carry = pd->bsignificand[i];
71			for (j = 0; j < len; j++) {
72				carry += (unsigned int)pb->bsignificand[j]
73				    * 10000;
74				pb->bsignificand[j] = carry & 0xffff;
75				carry >>= 16;
76			}
77			if (carry)
78				pb->bsignificand[j++] = carry;
79			len = j;
80		}
81	} else {
82		i = pd->blength - 1;
83		len = 0;
84	}
85
86	/* convert any partial digit */
87	if (i >= 0 && pd->bexponent + (i << 2) > -4) {
88		s = pd->bexponent + (i << 2) + 4;
89		/* multiply pb by 10^s and add partial digit */
90		carry = pd->bsignificand[i];
91		if (s == 1) {
92			s = carry % 1000;
93			carry = carry / 1000;
94			for (j = 0; j < len; j++) {
95				carry += (unsigned int)pb->bsignificand[j]
96				    * 10;
97				pb->bsignificand[j] = carry & 0xffff;
98				carry >>= 16;
99			}
100		} else if (s == 2) {
101			s = carry % 100;
102			carry = carry / 100;
103			for (j = 0; j < len; j++) {
104				carry += (unsigned int)pb->bsignificand[j]
105				    * 100;
106				pb->bsignificand[j] = carry & 0xffff;
107				carry >>= 16;
108			}
109		} else {
110			s = carry % 10;
111			carry = carry / 10;
112			for (j = 0; j < len; j++) {
113				carry += (unsigned int)pb->bsignificand[j]
114				    * 1000;
115				pb->bsignificand[j] = carry & 0xffff;
116				carry >>= 16;
117			}
118		}
119		if (carry)
120			pb->bsignificand[j++] = carry;
121		len = j;
122		i--;
123	} else {
124		s = 0;
125	}
126
127	pb->blength = len;
128	pb->bexponent = 0;
129
130	/* continue accumulating sticky flag */
131	while (i >= 0)
132		s |= pd->bsignificand[i--];
133	*sticky = s;
134
135	if (pd->bexponent > 0) {
136		/* scale pb by 10^pd->exponent */
137		__big_float_times_power(pb, 10, pd->bexponent, nsig, ppb);
138	}
139}
140
141/*
142 * Convert the decimal_record *pd to an unpacked datum *px accurately
143 * enough that *px can be rounded correctly to sigbits significant bits.
144 * (We may assume sigbits <= 113.)
145 */
146static void
147__decimal_to_unpacked(unpacked *px, decimal_record *pd, int sigbits)
148{
149	_big_float	d, b, *pbd, *pbb;
150	char		*ds;
151	int		ids, i, ix, exp, ndigs;
152	int		sticky, powtwo, sigdigits;
153
154	px->sign = pd->sign;
155	px->fpclass = pd->fpclass;
156	ds = pd->ds;
157	ndigs = pd->ndigits;
158	exp = pd->exponent;
159
160	/* remove trailing zeroes */
161	while (ndigs > 0 && ds[ndigs - 1] == '0') {
162		exp++;
163		ndigs--;
164	}
165	if (ndigs < 1) {
166		/* nothing left */
167		px->fpclass = fp_zero;
168		return;
169	}
170
171	/* convert remaining digits to a base-10^4 _big_float */
172	d.bsize = _BIG_FLOAT_SIZE;
173	d.bexponent = exp;
174	d.blength = (ndigs + 3) >> 2;
175	i = d.blength - 1;
176	ids = ndigs - (d.blength << 2);
177	switch (ids) {
178	case -1:
179		d.bsignificand[i] = 100 * ds[ids + 1] +
180		    10 * ds[ids + 2] + ds[ids + 3] - 111 * '0';
181		i--;
182		ids += 4;
183		break;
184
185	case -2:
186		d.bsignificand[i] = 10 * ds[ids + 2] + ds[ids + 3] - 11 * '0';
187		i--;
188		ids += 4;
189		break;
190
191	case -3:
192		d.bsignificand[i] = ds[ids + 3] - '0';
193		i--;
194		ids += 4;
195		break;
196	}
197	while (i >= 0) {
198		d.bsignificand[i] = 1000 * ds[ids] + 100 * ds[ids + 1] +
199		    10 * ds[ids + 2] + ds[ids + 3] - 1111 * '0';
200		i--;
201		ids += 4;
202	}
203
204	pbd = &d;
205	powtwo = 0;
206
207	/* pre-scale to get the bits we want into the integer part */
208	if (exp < 0) {
209		/* i is a lower bound on log10(x) */
210		i = exp + ndigs - 1;
211		if (i <= 0 || ((i * 217705) >> 16) < sigbits + 2) {
212			/*
213			 * Scale by 2^(sigbits + 2 + u) where
214			 * u is an upper bound on -log2(x).
215			 */
216			powtwo = sigbits + 2;
217			if (i < 0)
218				powtwo += ((-i * 217706) + 65535) >> 16;
219			else if (i > 0)
220				powtwo -= (i * 217705) >> 16;
221			/*
222			 * Take sigdigits large enough to get
223			 * all integral digits correct.
224			 */
225			sigdigits = i + 1 + (((powtwo * 19729) + 65535) >> 16);
226			__big_float_times_power(&d, 2, powtwo, sigdigits, &pbd);
227		}
228	}
229
230	/* convert to base 2^16 */
231	b.bsize = _BIG_FLOAT_SIZE;
232	pbb = &b;
233	__big_decimal_to_big_binary(pbd, sigbits + 2, &pbb, &sticky);
234
235	/* adjust pbb->bexponent based on the scale factor above */
236	pbb->bexponent -= powtwo;
237
238	/* convert to unpacked */
239	ix = 0;
240	for (i = pbb->blength - 1; i > 0 && ix < 5; i -= 2) {
241		px->significand[ix++] = (pbb->bsignificand[i] << 16) |
242		    pbb->bsignificand[i - 1];
243	}
244	if (ix < 5) {
245		/* pad with zeroes */
246		if (i == 0)
247			px->significand[ix++] = pbb->bsignificand[i] << 16;
248		while (ix < 5)
249			px->significand[ix++] = 0;
250	} else {
251		/* truncate and set a sticky bit if necessary */
252		while (i >= 0 && pbb->bsignificand[i] == 0)
253			i--;
254		if (i >= 0)
255			px->significand[4] |= 1;
256	}
257	if (sticky | pd->more)
258		px->significand[4] |= 1;
259	px->exponent = pbb->bexponent + (pbb->blength << 4) - 1;
260
261	/* normalize so the most significant bit is set */
262	while (px->significand[0] < 0x80000000u) {
263		px->significand[0] = (px->significand[0] << 1) |
264		    (px->significand[1] >> 31);
265		px->significand[1] = (px->significand[1] << 1) |
266		    (px->significand[2] >> 31);
267		px->significand[2] = (px->significand[2] << 1) |
268		    (px->significand[3] >> 31);
269		px->significand[3] = (px->significand[3] << 1) |
270		    (px->significand[4] >> 31);
271		px->significand[4] <<= 1;
272		px->exponent--;
273	}
274
275	if (pbd != &d)
276		(void) free((void *)pbd);
277	if (pbb != &b)
278		(void) free((void *)pbb);
279}
280
281/*
282 * Convert a string s consisting of n <= 18 ASCII decimal digits
283 * to an integer value in double precision format, and set *pe
284 * to the number of rounding errors incurred (0 or 1).
285 */
286static double
287__digits_to_double(char *s, int n, int *pe)
288{
289	int	i, acc;
290	double	t, th, tl;
291
292	if (n <= 9) {
293		acc = s[0] - '0';
294		for (i = 1; i < n; i++) {
295			/* acc <- 10 * acc + next digit */
296			acc = (acc << 1) + (acc << 3) + s[i] - '0';
297		}
298		t = (double)acc;
299		*pe = 0;
300	} else {
301		acc = s[0] - '0';
302		for (i = 1; i < (n - 9); i++) {
303			/* acc <- 10 * acc + next digit */
304			acc = (acc << 1) + (acc << 3) + s[i] - '0';
305		}
306		th = 1.0e9 * (double)acc;	/* this will be exact */
307		acc = s[n - 9] - '0';
308		for (i = n - 8; i < n; i++) {
309			/* acc <- 10 * acc + next digit */
310			acc = (acc << 1) + (acc << 3) + s[i] - '0';
311		}
312		tl = (double)acc;
313		/* add and indicate whether or not the sum is exact */
314		t = th + tl;
315		*pe = ((t - th) == tl)? 0 : 1;
316	}
317	return (t);
318}
319
320static union {
321	int	i[2];
322	double	d;
323} C[] = {
324#ifdef _LITTLE_ENDIAN
325	{ 0x00000000, 0x3cc40000 },
326#else
327	{ 0x3cc40000, 0x00000000 },	/* 5 * 2^-53 */
328#endif
329};
330
331#define	five2m53	C[0].d
332
333static int
334__fast_decimal_to_single(single *px, decimal_mode *pm, decimal_record *pd,
335    fp_exception_field_type *ps)
336{
337	double			dds, delta, ddsplus, ddsminus, df1;
338	int			n, exp, rounded, e;
339	float			f1, f2;
340	__ieee_flags_type	fb;
341
342	if (pm->rd != fp_nearest)
343		return (0);
344
345	exp = pd->exponent;
346	if (pd->ndigits <= 18) {
347		rounded = 0;
348		n = pd->ndigits;
349	} else {
350		rounded = 1;
351		n = 18;
352		exp += pd->ndigits - 18;
353	}
354	/*
355	 * exp must be in the range of the table, and the result
356	 * must not underflow or overflow.
357	 */
358	if (exp < -__TBL_TENS_MAX || exp + n < -36 || exp + n > 38)
359		return (0);
360
361	__get_ieee_flags(&fb);
362	dds = __digits_to_double(pd->ds, n, &e);
363	if (e != 0)
364		rounded = 1;
365	if (exp > 0) {
366		/* small positive exponent */
367		if (exp > __TBL_TENS_EXACT)
368			rounded = 1;
369		if (rounded) {
370			dds *= __tbl_tens[exp];
371		} else {
372			dds = __mul_set(dds, __tbl_tens[exp], &e);
373			if (e)
374				rounded = 1;
375		}
376	} else if (exp < 0) {
377		/* small negative exponent */
378		if (-exp > __TBL_TENS_EXACT)
379			rounded = 1;
380		if (rounded) {
381			dds /= __tbl_tens[-exp];
382		} else {
383			dds = __div_set(dds, __tbl_tens[-exp], &e);
384			if (e)
385				rounded = 1;
386		}
387	}
388
389	/*
390	 * At this point dds may have four rounding errors due to
391	 * (i) truncation of pd->ds to 18 digits, (ii) inexact con-
392	 * version of pd->ds to binary, (iii) scaling by a power of
393	 * ten that is not exactly representable, and (iv) roundoff
394	 * error in the multiplication.  Below we will incur one more
395	 * rounding error when we add or subtract delta and dds.  We
396	 * construct delta so that even after this last rounding error,
397	 * ddsplus is an upper bound on the exact value and ddsminus
398	 * is a lower bound.  Then as long as both of these quantities
399	 * round to the same single precision number, that number
400	 * will be the correctly rounded single precision result.
401	 * (If any rounding errors have been committed, then we must
402	 * also be certain that the result can't be exact.)
403	 */
404	delta = five2m53 * dds;
405	ddsplus = dds + delta;
406	ddsminus = dds - delta;
407	f1 = (float)(ddsplus);
408	f2 = (float)(ddsminus);
409	df1 = f1;
410	__set_ieee_flags(&fb);
411	if (f1 != f2)
412		return (0);
413	if (rounded) {
414		/*
415		 * If ddsminus <= f1 <= ddsplus, the result might be
416		 * exact; we have to convert the long way to be sure.
417		 */
418		if (ddsminus <= df1 && df1 <= ddsplus)
419			return (0);
420		*ps = (1 << fp_inexact);
421	} else {
422		*ps = (df1 == dds)? 0 : (1 << fp_inexact);
423	}
424	*px = (pd->sign)? -f1 : f1;
425	return (1);
426}
427
428/*
429 * Attempt conversion to double using floating-point arithmetic.
430 * Return 1 if it works (at most one rounding error), 0 if it doesn't.
431 */
432static int
433__fast_decimal_to_double(double *px, decimal_mode *pm, decimal_record *pd,
434    fp_exception_field_type *ps)
435{
436	double			dds;
437	int			e;
438	__ieee_flags_type	fb;
439
440	if (pm->rd != fp_nearest || pd->ndigits > 18 || pd->exponent
441	    > __TBL_TENS_EXACT || pd->exponent < -__TBL_TENS_EXACT)
442		return (0);
443
444	__get_ieee_flags(&fb);
445	dds = __digits_to_double(pd->ds, pd->ndigits, &e);
446	if (e != 0) {
447		__set_ieee_flags(&fb);
448		return (0);
449	}
450	if (pd->exponent > 0)
451		dds = __mul_set(dds, __tbl_tens[pd->exponent], &e);
452	else if (pd->exponent < 0)
453		dds = __div_set(dds, __tbl_tens[-pd->exponent], &e);
454	*px = (pd->sign)? -dds : dds;
455	*ps = (e)? (1 << fp_inexact) : 0;
456	__set_ieee_flags(&fb);
457	return (1);
458}
459
460/* PUBLIC FUNCTIONS */
461
462/*
463 * The following routines convert the decimal record *pd to a floating
464 * point value *px observing the rounding mode specified in pm->rd and
465 * passing back any floating point exceptions that occur in *ps.
466 *
467 * pd->sign and pd->fpclass are always taken into account.  pd->exponent
468 * and pd->ds are used when pd->fpclass is fp_normal or fp_subnormal.
469 * In these cases, pd->ds must contain a null-terminated string of one
470 * or more ASCII digits, the first of which is not zero, and pd->ndigits
471 * must equal the length of this string.  If m is the integer represented
472 * by the string pd->ds, then *px will be set to a correctly rounded
473 * approximation to
474 *
475 *   -1**(pd->sign) * m * 10**(pd->exponent)
476 *
477 * (If pd->more != 0 then additional nonzero digits are assumed to follow
478 * those in pd->ds, so m is effectively replaced by m + epsilon in the
479 * expression above.)
480 *
481 * For example, if pd->exponent == -2 and pd->ds holds "1234", then *px
482 * will be a correctly rounded approximation to 12.34.
483 *
484 * Note that the only mode that matters is the rounding direction pm->rd;
485 * pm->df and pm->ndigits are never used.
486 */
487
488/* maximum decimal exponent we need to consider */
489#define	SINGLE_MAXE	  47
490#define	DOUBLE_MAXE	 326
491#define	EXTENDED_MAXE	4968
492#define	QUAD_MAXE	4968
493
494void
495decimal_to_single(single *px, decimal_mode *pm, decimal_record *pd,
496    fp_exception_field_type *ps)
497{
498	single_equivalence	*kluge;
499	unpacked		u;
500	fp_exception_field_type	ef;
501	int			i;
502
503	/* special values */
504	kluge = (single_equivalence *)px;
505	switch (pd->fpclass) {
506	case fp_zero:
507		kluge->f.msw.sign = (pd->sign)? 1 : 0;
508		kluge->f.msw.exponent = 0;
509		kluge->f.msw.significand = 0;
510		*ps = 0;
511		return;
512
513	case fp_infinity:
514		kluge->f.msw.sign = (pd->sign)? 1 : 0;
515		kluge->f.msw.exponent = 0xff;
516		kluge->f.msw.significand = 0;
517		*ps = 0;
518		return;
519
520	case fp_quiet:
521		kluge->f.msw.sign = (pd->sign)? 1 : 0;
522		kluge->f.msw.exponent = 0xff;
523		kluge->f.msw.significand = 0x7fffff;
524		*ps = 0;
525		return;
526
527	case fp_signaling:
528		kluge->f.msw.sign = (pd->sign)? 1 : 0;
529		kluge->f.msw.exponent = 0xff;
530		kluge->f.msw.significand = 0x3fffff;
531		*ps = 0;
532		return;
533	}
534
535	/* numeric values */
536	ef = 0;
537	if (pd->exponent + pd->ndigits > SINGLE_MAXE) {
538		/* result must overflow */
539		u.sign = (pd->sign != 0);
540		u.fpclass = fp_normal;
541		u.exponent = 0x000fffff;
542		u.significand[0] = 0x80000000;
543		for (i = 1; i < UNPACKED_SIZE; i++)
544			u.significand[i] = 0;
545	} else if (pd->exponent + pd->ndigits < -SINGLE_MAXE) {
546		/* result must underflow completely */
547		u.sign = (pd->sign != 0);
548		u.fpclass = fp_normal;
549		u.exponent = -0x000fffff;
550		u.significand[0] = 0x80000000;
551		for (i = 1; i < UNPACKED_SIZE; i++)
552			u.significand[i] = 0;
553	} else {
554		/* result may be in range */
555		if (__fast_decimal_to_single(px, pm, pd, &ef) == 1) {
556			*ps = ef;
557			if (ef != 0)
558				__base_conversion_set_exception(ef);
559			return;
560		}
561		__decimal_to_unpacked(&u, pd, 24);
562	}
563	__pack_single(&u, px, pm->rd, &ef);
564	*ps = ef;
565	if (ef != 0)
566		__base_conversion_set_exception(ef);
567}
568
569void
570decimal_to_double(double *px, decimal_mode *pm, decimal_record *pd,
571    fp_exception_field_type *ps)
572{
573	double_equivalence	*kluge;
574	unpacked		u;
575	fp_exception_field_type	ef;
576	int			i;
577
578	/* special values */
579	kluge = (double_equivalence *)px;
580	switch (pd->fpclass) {
581	case fp_zero:
582		kluge->f.msw.sign = (pd->sign)? 1 : 0;
583		kluge->f.msw.exponent = 0;
584		kluge->f.msw.significand = 0;
585		kluge->f.significand2 = 0;
586		*ps = 0;
587		return;
588
589	case fp_infinity:
590		kluge->f.msw.sign = (pd->sign)? 1 : 0;
591		kluge->f.msw.exponent = 0x7ff;
592		kluge->f.msw.significand = 0;
593		kluge->f.significand2 = 0;
594		*ps = 0;
595		return;
596
597	case fp_quiet:
598		kluge->f.msw.sign = (pd->sign)? 1 : 0;
599		kluge->f.msw.exponent = 0x7ff;
600		kluge->f.msw.significand = 0xfffff;
601		kluge->f.significand2 = 0xffffffff;
602		*ps = 0;
603		return;
604
605	case fp_signaling:
606		kluge->f.msw.sign = (pd->sign)? 1 : 0;
607		kluge->f.msw.exponent = 0x7ff;
608		kluge->f.msw.significand = 0x7ffff;
609		kluge->f.significand2 = 0xffffffff;
610		*ps = 0;
611		return;
612	}
613
614	/* numeric values */
615	ef = 0;
616	if (pd->exponent + pd->ndigits > DOUBLE_MAXE) {
617		/* result must overflow */
618		u.sign = (pd->sign != 0);
619		u.fpclass = fp_normal;
620		u.exponent = 0x000fffff;
621		u.significand[0] = 0x80000000;
622		for (i = 1; i < UNPACKED_SIZE; i++)
623			u.significand[i] = 0;
624	} else if (pd->exponent + pd->ndigits < -DOUBLE_MAXE) {
625		/* result must underflow completely */
626		u.sign = (pd->sign != 0);
627		u.fpclass = fp_normal;
628		u.exponent = -0x000fffff;
629		u.significand[0] = 0x80000000;
630		for (i = 1; i < UNPACKED_SIZE; i++)
631			u.significand[i] = 0;
632	} else {
633		/* result may be in range */
634		if (__fast_decimal_to_double(px, pm, pd, &ef) == 1) {
635			*ps = ef;
636			if (ef != 0)
637				__base_conversion_set_exception(ef);
638			return;
639		}
640		__decimal_to_unpacked(&u, pd, 53);
641	}
642	__pack_double(&u, px, pm->rd, &ef);
643	*ps = ef;
644	if (ef != 0)
645		__base_conversion_set_exception(ef);
646}
647
648void
649decimal_to_extended(extended *px, decimal_mode *pm, decimal_record *pd,
650    fp_exception_field_type *ps)
651{
652	extended_equivalence	*kluge;
653	unpacked		u;
654	double_equivalence	dd;
655	fp_exception_field_type ef;
656	int			i;
657
658	/* special values */
659	kluge = (extended_equivalence *)px;
660	switch (pd->fpclass) {
661	case fp_zero:
662		kluge->f.msw.sign = (pd->sign)? 1 : 0;
663		kluge->f.msw.exponent = 0;
664		kluge->f.significand = 0;
665		kluge->f.significand2 = 0;
666		*ps = 0;
667		return;
668
669	case fp_infinity:
670		kluge->f.msw.sign = (pd->sign)? 1 : 0;
671		kluge->f.msw.exponent = 0x7fff;
672		kluge->f.significand = 0x80000000;
673		kluge->f.significand2 = 0;
674		*ps = 0;
675		return;
676
677	case fp_quiet:
678		kluge->f.msw.sign = (pd->sign)? 1 : 0;
679		kluge->f.msw.exponent = 0x7fff;
680		kluge->f.significand = 0xffffffff;
681		kluge->f.significand2 = 0xffffffff;
682		*ps = 0;
683		return;
684
685	case fp_signaling:
686		kluge->f.msw.sign = (pd->sign)? 1 : 0;
687		kluge->f.msw.exponent = 0x7fff;
688		kluge->f.significand = 0xbfffffff;
689		kluge->f.significand2 = 0xffffffff;
690		*ps = 0;
691		return;
692	}
693
694	/* numeric values */
695	ef = 0;
696	if (pd->exponent + pd->ndigits > EXTENDED_MAXE) {
697		/* result must overflow */
698		u.sign = (pd->sign != 0);
699		u.fpclass = fp_normal;
700		u.exponent = 0x000fffff;
701		u.significand[0] = 0x80000000;
702		for (i = 1; i < UNPACKED_SIZE; i++)
703			u.significand[i] = 0;
704	} else if (pd->exponent + pd->ndigits < -EXTENDED_MAXE) {
705		/* result must underflow completely */
706		u.sign = (pd->sign != 0);
707		u.fpclass = fp_normal;
708		u.exponent = -0x000fffff;
709		u.significand[0] = 0x80000000;
710		for (i = 1; i < UNPACKED_SIZE; i++)
711			u.significand[i] = 0;
712	} else {
713		/* result may be in range */
714		if (__fast_decimal_to_double(&dd.x, pm, pd, &ef) == 1 &&
715		    ef == 0) {
716			u.sign = dd.f.msw.sign;
717			u.fpclass = fp_normal;
718			u.exponent = dd.f.msw.exponent - DOUBLE_BIAS;
719			u.significand[0] = ((0x100000 |
720			    dd.f.msw.significand) << 11) |
721			    (dd.f.significand2 >> 21);
722			u.significand[1] = dd.f.significand2 << 11;
723			for (i = 2; i < UNPACKED_SIZE; i++)
724				u.significand[i] = 0;
725		} else {
726			__decimal_to_unpacked(&u, pd, 64);
727		}
728	}
729	__pack_extended(&u, px, pm->rd, &ef);
730	*ps = ef;
731	if (ef != 0)
732		__base_conversion_set_exception(ef);
733}
734
735void
736decimal_to_quadruple(quadruple *px, decimal_mode *pm, decimal_record *pd,
737    fp_exception_field_type *ps)
738{
739	quadruple_equivalence	*kluge;
740	unpacked		u;
741	double_equivalence	dd;
742	fp_exception_field_type ef;
743	int			i;
744
745	/* special values */
746	kluge = (quadruple_equivalence *)px;
747	switch (pd->fpclass) {
748	case fp_zero:
749		kluge->f.msw.sign = (pd->sign)? 1 : 0;
750		kluge->f.msw.exponent = 0;
751		kluge->f.msw.significand = 0;
752		kluge->f.significand2 = 0;
753		kluge->f.significand3 = 0;
754		kluge->f.significand4 = 0;
755		*ps = 0;
756		return;
757
758	case fp_infinity:
759		kluge->f.msw.sign = (pd->sign)? 1 : 0;
760		kluge->f.msw.exponent = 0x7fff;
761		kluge->f.msw.significand = 0;
762		kluge->f.significand2 = 0;
763		kluge->f.significand3 = 0;
764		kluge->f.significand4 = 0;
765		*ps = 0;
766		return;
767
768	case fp_quiet:
769		kluge->f.msw.sign = (pd->sign)? 1 : 0;
770		kluge->f.msw.exponent = 0x7fff;
771		kluge->f.msw.significand = 0xffff;
772		kluge->f.significand2 = 0xffffffff;
773		kluge->f.significand3 = 0xffffffff;
774		kluge->f.significand4 = 0xffffffff;
775		*ps = 0;
776		return;
777
778	case fp_signaling:
779		kluge->f.msw.sign = (pd->sign)? 1 : 0;
780		kluge->f.msw.exponent = 0x7fff;
781		kluge->f.msw.significand = 0x7fff;
782		kluge->f.significand2 = 0xffffffff;
783		kluge->f.significand3 = 0xffffffff;
784		kluge->f.significand4 = 0xffffffff;
785		*ps = 0;
786		return;
787	}
788
789	/* numeric values */
790	ef = 0;
791	if (pd->exponent + pd->ndigits > QUAD_MAXE) {
792		/* result must overflow */
793		u.sign = (pd->sign != 0);
794		u.fpclass = fp_normal;
795		u.exponent = 0x000fffff;
796		u.significand[0] = 0x80000000;
797		for (i = 1; i < UNPACKED_SIZE; i++)
798			u.significand[i] = 0;
799	} else if (pd->exponent + pd->ndigits < -QUAD_MAXE) {
800		/* result must underflow completely */
801		u.sign = (pd->sign != 0);
802		u.fpclass = fp_normal;
803		u.exponent = -0x000fffff;
804		u.significand[0] = 0x80000000;
805		for (i = 1; i < UNPACKED_SIZE; i++)
806			u.significand[i] = 0;
807	} else {
808		/* result may be in range */
809		if (__fast_decimal_to_double(&dd.x, pm, pd, &ef) == 1 &&
810		    ef == 0) {
811			u.sign = dd.f.msw.sign;
812			u.fpclass = fp_normal;
813			u.exponent = dd.f.msw.exponent - DOUBLE_BIAS;
814			u.significand[0] = ((0x100000 |
815			    dd.f.msw.significand) << 11) |
816			    (dd.f.significand2 >> 21);
817			u.significand[1] = dd.f.significand2 << 11;
818			for (i = 2; i < UNPACKED_SIZE; i++)
819				u.significand[i] = 0;
820		} else {
821			__decimal_to_unpacked(&u, pd, 113);
822		}
823	}
824	__pack_quadruple(&u, px, pm->rd, &ef);
825	*ps = ef;
826	if (ef != 0)
827		__base_conversion_set_exception(ef);
828}
829