double_decim.c revision 7c478bd95313f5f23a4c958a745db2134aa03244
/*
* CDDL HEADER START
*
* The contents of this file are subject to the terms of the
* Common Development and Distribution License, Version 1.0 only
* (the "License"). You may not use this file except in compliance
* with the License.
*
* You can obtain a copy of the license at usr/src/OPENSOLARIS.LICENSE
* See the License for the specific language governing permissions
* and limitations under the License.
*
* When distributing Covered Code, include this CDDL HEADER in each
* file and include the License file at usr/src/OPENSOLARIS.LICENSE.
* If applicable, add the following below this CDDL HEADER, with the
* fields enclosed by brackets "[]" replaced with your own identifying
* information: Portions Copyright [yyyy] [name of copyright owner]
*
* CDDL HEADER END
*/
/*
* Copyright 2005 Sun Microsystems, Inc. All rights reserved.
* Use is subject to license terms.
*/
#pragma ident "%Z%%M% %I% %E% SMI"
/*
* Conversion from binary to decimal floating point
*/
#include "synonyms.h"
#include <stdlib.h>
#include "base_conversion.h"
/*
* Any sensible programmer would inline the following routine where
* it is used below. Unfortunately, the Sun SPARC compilers are not
* consistent in generating efficient code for this, so inlining it
* as written can cause the *_to_decimal functions to take twice as
* long in some cases.
*
* We might be tempted, then, to rewrite the source to match the most
* efficient code the compilers generate and inline that. Alas, the
* most efficient code on SPARC uses 32x32->64 bit multiply, which
* can't be expressed directly in source code. We could use long long,
* which would imply 64x64->64 bit multiply; this would work perfectly
* well on SPARC in v8plus mode. But as of Solaris 10, libc for SPARC
* is still built in v8 mode, and of course, x86 is another story.
*
* We could also choose to use an inline template to get the most
* efficient code without incurring the full cost of a function call.
* Since I expect that would not buy much performance gain, and I
* prefer to avoid using inline templates for things that can be
* written in a perfectly straightforward way in C, I've settled
* for this implementation. I hope that someday the compilers will
* do this.
*/
static unsigned int
__quorem10000(unsigned int x, unsigned short *pr)
{
*pr = x % 10000;
return (x / 10000);
}
/*
* Convert the integer part of a nonzero base-2^16 _big_float *pb
* to base 10^4 in **ppd. The converted value is accurate to nsig
* significant digits. On exit, *sticky is nonzero if *pb had a
* nonzero fractional part. If pb->exponent > 0 and **ppd is not
* large enough to hold the final converted value (i.e., the con-
* verted significand scaled by 2^pb->exponent), then on exit,
* *ppd will point to a newly allocated _big_float, which must be
* freed by the caller. (The number of significant digits we need
* should fit in pd, but __big_float_times_power may allocate new
* storage anyway because we could be multiplying by as much as
* 2^16271, which would require more than 4000 digits.)
*
* This routine does not check that **ppd is large enough to hold
* the result of converting the significand of *pb.
*/
static void
int *sticky)
{
_big_float *pd;
int i, j, len, s;
unsigned int carry;
/* convert pb a digit at a time, most significant first */
&pd->bsignificand[0]);
/* multiply pd by 2^16 and add next digit */
for (j = 0; j < len; j++) {
<< 16;
&pd->bsignificand[j]);
}
while (carry != 0) {
&pd->bsignificand[j]);
j++;
}
len = j;
}
} else {
len = 0;
}
/* convert any partial digit */
/* multiply pd by 2^s and add partial digit */
for (j = 0; j < len; j++) {
}
while (carry != 0) {
j++;
}
len = j;
i--;
} else {
s = 0;
}
/* continue accumulating sticky flag */
while (i >= 0)
s |= pb->bsignificand[i--];
*sticky = s;
/* scale pd by 2^pb->bexponent */
}
}
/*
* Convert a base-10^4 _big_float *pf to a decimal string in *pd,
* rounding according to the modes in *pm and recording any exceptions
* in *ps. If sticky is nonzero, then additional nonzero digits are
* assumed to follow those in *pf. pd->sign must have already been
* filled in, and pd->fpclass is not modified. The resulting string
* is stored in pd->ds, terminated by a null byte. The length of this
* string is stored in pd->ndigits, and the corresponding exponent
* is stored in pd->exponent. If the converted value is not exact,
* the inexact flag is set in *ps.
*
* When pm->df == fixed_form, we may discover that the result would
* have more than DECIMAL_STRING_LENGTH - 1 digits. In this case,
* we put DECIMAL_STRING_LENGTH - 1 digits into *pd, adjusting both
* the exponent and the decimal place at which the value is rounded
* as need be, and we set the overflow flag in *ps. (Raising overflow
* is a bug, but we have to do it to maintain backward compatibility.)
*
* *pf may be modified.
*/
static void
{
unsigned short d;
char s[4], round;
/* set e = floor(log10(*pf)) */
if (i < 0) {
} else {
d = pf->bsignificand[i];
if (d >= 1000)
e += 3;
else if (d >= 100)
e += 2;
else if (d >= 10)
e++;
}
/*
* Determine the power of ten after which to round and the
* powers corresponding to the first and last digits desired
* in the result.
*/
/* F format */
if (er < 0) {
efirst = (e >= 0)? e : -1;
} else {
elast = 0;
}
/* check for possible overflow of pd->ds */
efirst = e;
}
} else {
/* E format */
efirst = e;
}
/* retrieve digits down to the (er - 1) place */
is = 0;
e >= er - 1; e--)
if (j > 0 && e >= er - 1) {
e--;
}
while (j <= 3)
sticky |= (s[j++] - '0');
i--;
}
is += 4;
e -= 4;
i--;
}
if (i >= 0) {
if (e >= er - 1) {
for (j = 0; e >= er - 1; j++) {
e--;
}
while (j <= 3)
sticky |= (s[j++] - '0');
i--;
}
} else {
while (e-- >= er - 1)
}
/* collect rounding information */
while (i >= 0)
/* add more trailing zeroes if need be */
/* round */
return;
case fp_nearest:
return;
break;
case fp_positive:
return;
break;
case fp_negative:
return;
break;
default:
return;
}
/* round up */
if (i >= 0) {
} else {
/* rounding carry out has occurred */
} else {
if (is > 0)
is++;
}
}
}
/*
* Convert a binary floating point value represented by *pf to a
* decimal record *pd according to the modes in *pm. Any exceptions
* incurred are passed back via *ps.
*/
static void
{
/*
* If pm->ndigits is too large or too small, set the overflow
* flag in *ps and do nothing. (Raising overflow is a bug,
* but we have to do it to maintain backward compatibility.)
*/
return;
}
powten = 0;
/* pre-scale to get the digits we want into the integer part */
/* F format */
/*
* Scale by 10^min(-bf->bexponent, pm->ndigits + 1).
*/
/*
* Take sigbits large enough to get all integral
* digits correct.
*/
if (sigbits < 1)
sigbits = 1;
}
} else {
/* E format */
/* i is a lower bound on log2(x) */
/*
* Scale by 10^min(-bf->bexponent,
* pm->ndigits + 1 + u) where u is
* an upper bound on -log10(x).
*/
if (i < 0)
else if (i > 0)
/*
* Take sigbits large enough to get
* all integral digits correct.
*/
sigbits = i + 16 +
}
}
}
/* convert to base 10^4 */
d.bsize = _BIG_FLOAT_SIZE;
pbd = &d;
/* adjust pbd->bexponent based on the scale factor above */
/* convert to ASCII */
if (pbd != &d)
}
/* remove trailing zeroes from the significand of p */
static void
__shorten(_big_float *p)
{
int zeros, i;
/* count trailing zeros */
;
if (zeros) {
for (i = 0; i < length; i++)
}
}
/*
* Unpack a normal or subnormal double into a _big_float.
*/
static void
{
x = (double_equivalence *)px;
} else {
}
}
/*
* Unpack a normal or subnormal extended into a _big_float.
*/
static void
{
x = (extended_equivalence *)px;
}
}
/*
* Unpack a normal or subnormal quad into a _big_float.
*/
static void
{
x = (quadruple_equivalence *)px;
} else {
}
}
/* PUBLIC ROUTINES */
void
{
double x;
/* decide what to do based on the class of x */
*ps = 0;
return;
} else {
#if defined(__sparc)
int i;
/*
* On SPARC, simply converting *px to double
* can flush a subnormal value to zero when
* nonstandard mode is enabled, so we have
* to go through all this nonsense instead.
*/
i = *(int *)px;
x = (double)(i & ~0x80000000);
if (i < 0)
x = -x;
x *= 1.401298464324817070923730e-45; /* 2^-149 */
ef = 0;
__double_to_bigfloat(&x, &bf);
}
if (ef != 0)
return;
#else
#endif
}
else
*ps = 0;
return;
} else {
}
ef = 0;
x = *px;
__double_to_bigfloat(&x, &bf);
}
if (ef != 0)
}
void
{
/* decide what to do based on the class of x */
kluge->f.significand2 == 0) {
*ps = 0;
return;
} else {
}
kluge->f.significand2 == 0)
else
*ps = 0;
return;
} else {
}
ef = 0;
}
if (ef != 0)
}
void
{
/* decide what to do based on the class of x */
*ps = 0;
return;
} else {
/*
* x could be a pseudo-denormal, but the distinction
* doesn't matter
*/
}
/*
* In Intel's extended format, if the exponent is
* nonzero but the explicit integer bit is zero, this
* is an "unsupported format" bit pattern; treat it
* like a signaling NaN.
*/
*ps = 0;
return;
kluge->f.significand2) == 0)
else
*ps = 0;
return;
} else {
}
ef = 0;
if (ef != 0)
}
void
{
/* decide what to do based on the class of x */
kluge->f.significand4) == 0) {
*ps = 0;
return;
} else {
}
kluge->f.significand4) == 0)
else
*ps = 0;
return;
} else {
}
ef = 0;
if (ef != 0)
}