500 lines
10 KiB
C
500 lines
10 KiB
C
/* Floating-point printing for `printf'.
|
|
This is an implementation of a restricted form of the `Dragon4'
|
|
algorithm described in "How to Print Floating-Point Numbers Accurately",
|
|
by Guy L. Steele, Jr. and Jon L. White, presented at the ACM SIGPLAN '90
|
|
Conference on Programming Language Design and Implementation.
|
|
|
|
Copyright (C) 1992 Free Software Foundation, Inc.
|
|
This file is part of the GNU C Library.
|
|
|
|
The GNU C Library is free software; you can redistribute it and/or
|
|
modify it under the terms of the GNU Library General Public License as
|
|
published by the Free Software Foundation; either version 2 of the
|
|
License, or (at your option) any later version.
|
|
|
|
The GNU C Library is distributed in the hope that it will be useful,
|
|
but WITHOUT ANY WARRANTY; without even the implied warranty of
|
|
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
|
|
Library General Public License for more details.
|
|
|
|
You should have received a copy of the GNU Library General Public
|
|
License along with the GNU C Library; see the file COPYING.LIB. If
|
|
not, write to the Free Software Foundation, Inc., 675 Mass Ave,
|
|
Cambridge, MA 02139, USA. */
|
|
|
|
#include <ansidecl.h>
|
|
#include <ctype.h>
|
|
#include <stdio.h>
|
|
#include <float.h>
|
|
#include <math.h>
|
|
#include <stdarg.h>
|
|
#include <stdlib.h>
|
|
#include <localeinfo.h>
|
|
|
|
#include <printf.h>
|
|
|
|
#define outchar(x) \
|
|
do \
|
|
{ \
|
|
register CONST int outc = (x); \
|
|
if (putc(outc, s) == EOF) \
|
|
return -1; \
|
|
else \
|
|
++done; \
|
|
} while (0)
|
|
|
|
#if FLT_RADIX != 2
|
|
|
|
double
|
|
frexp (double f, int *e)
|
|
{
|
|
#error "Don't know how to extract fraction and exponent from `double'."
|
|
}
|
|
|
|
#undef ldexp
|
|
#ifdef __GNUC__
|
|
inline
|
|
#endif
|
|
static double
|
|
ldexp (double f, int e)
|
|
{
|
|
while (e > 0)
|
|
{
|
|
f *= FLT_RADIX;
|
|
--e;
|
|
}
|
|
while (e < 0)
|
|
{
|
|
f /= FLT_RADIX;
|
|
++e;
|
|
}
|
|
}
|
|
|
|
#endif
|
|
|
|
|
|
int
|
|
DEFUN(__printf_fp, (s, info, args),
|
|
FILE *s AND CONST struct printf_info *info AND va_list *args)
|
|
{
|
|
int done = 0;
|
|
|
|
/* Decimal point character. */
|
|
CONST char *CONST decimal = _numeric_info->decimal_point;
|
|
|
|
LONG_DOUBLE fpnum; /* Input. */
|
|
int is_neg;
|
|
|
|
LONG_DOUBLE f; /* Fraction. */
|
|
|
|
int e; /* Base-2 exponent of the input. */
|
|
CONST int p = DBL_MANT_DIG; /* Internal precision. */
|
|
LONG_DOUBLE scale, scale10; /* Scale factor. */
|
|
LONG_DOUBLE loerr, hierr; /* Potential error in the fraction. */
|
|
int k; /* Digits to the left of the decimal point. */
|
|
int cutoff; /* Where to stop generating digits. */
|
|
LONG_DOUBLE r, r2, r10; /* Remainder. */
|
|
int roundup;
|
|
int low, high;
|
|
char digit;
|
|
|
|
int j;
|
|
|
|
char type = tolower (info->spec);
|
|
int prec = info->prec;
|
|
int width = info->width;
|
|
|
|
/* This algorithm has the nice property of not needing a buffer.
|
|
However, to get the padding right for %g format, we need to know
|
|
the length of the number before printing it. */
|
|
char *buf = __alloca ((prec > LDBL_DIG ? prec : LDBL_DIG) +
|
|
LDBL_MAX_10_EXP + 3); /* Dot, e, exp. sign. */
|
|
register char *bp = buf;
|
|
#define put(c) *bp++ = (c)
|
|
|
|
/* Fetch the argument value. */
|
|
if (info->is_long_double)
|
|
fpnum = va_arg (*args, LONG_DOUBLE);
|
|
else
|
|
fpnum = (LONG_DOUBLE) va_arg (*args, double);
|
|
|
|
#ifdef HANDLE_SPECIAL
|
|
/* Allow for machine-dependent (or floating point format-dependent) code. */
|
|
HANDLE_SPECIAL (done, s, info, fpnum);
|
|
#endif
|
|
|
|
#ifndef IS_NEGATIVE
|
|
#define IS_NEGATIVE(num) ((num) < 0)
|
|
#endif
|
|
|
|
is_neg = IS_NEGATIVE (fpnum);
|
|
if (is_neg)
|
|
fpnum = - fpnum;
|
|
|
|
if (prec == -1)
|
|
prec = 6;
|
|
|
|
if (type == 'g')
|
|
{
|
|
if (prec == 0)
|
|
prec = 1;
|
|
|
|
if (fpnum != 0)
|
|
{
|
|
if (fpnum < 1e-4)
|
|
type = 'e';
|
|
else
|
|
{
|
|
f = 10;
|
|
j = prec;
|
|
if (j > p)
|
|
j = p;
|
|
while (--j > 0)
|
|
{
|
|
f *= 10;
|
|
if (f > fpnum)
|
|
{
|
|
type = 'e';
|
|
break;
|
|
}
|
|
}
|
|
}
|
|
}
|
|
|
|
/* For 'g'/'G' format, the precision specifies "significant digits",
|
|
not digits to come after the decimal point. */
|
|
--prec;
|
|
}
|
|
|
|
if (fpnum == 0)
|
|
/* Special case for zero.
|
|
The general algorithm does not work for zero. */
|
|
{
|
|
put ('0');
|
|
if (tolower (info->spec) != 'g' || info->alt)
|
|
{
|
|
if (prec > 0 || info->alt)
|
|
put (*decimal);
|
|
while (--prec > 0)
|
|
put ('0');
|
|
}
|
|
if (type == 'e')
|
|
{
|
|
put (info->spec);
|
|
put ('+');
|
|
put ('0');
|
|
put ('0');
|
|
}
|
|
}
|
|
else
|
|
{
|
|
/* Split the number into a fraction and base-2 exponent. */
|
|
f = frexp (fpnum, &e);
|
|
|
|
/* Scale the fractional part by the highest possible number of
|
|
significant bits of fraction. We want to represent the
|
|
fractional part as a (very) large integer. */
|
|
f = ldexp (f, p);
|
|
|
|
cutoff = -prec;
|
|
|
|
roundup = 0;
|
|
|
|
if (e > p)
|
|
{
|
|
/* The exponent is bigger than the number of fractional digits. */
|
|
r = ldexp (f, e - p);
|
|
scale = 1;
|
|
/* The number is (E - P) factors of two larger than
|
|
the fraction can represent; this is the potential error. */
|
|
loerr = ldexp (1.0, e - p);
|
|
}
|
|
else
|
|
{
|
|
/* The number of fractional digits is greater than the exponent.
|
|
Scale by the difference factors of two. */
|
|
r = f;
|
|
scale = ldexp (1.0, p - e);
|
|
loerr = 1.0;
|
|
}
|
|
hierr = loerr;
|
|
|
|
/* Fixup. */
|
|
|
|
if (f == ldexp (1.0, p - 1))
|
|
{
|
|
/* Account for unequal gaps. */
|
|
hierr = ldexp (hierr, 1);
|
|
r = ldexp (r, 1);
|
|
scale = ldexp (scale, 1);
|
|
}
|
|
|
|
scale10 = ceil (scale / 10.0);
|
|
k = 0;
|
|
while (r < scale10)
|
|
{
|
|
--k;
|
|
r *= 10;
|
|
loerr *= 10;
|
|
hierr *= 10;
|
|
}
|
|
do
|
|
{
|
|
r2 = 2 * r;
|
|
|
|
while (r2 + hierr >= 2 * scale)
|
|
{
|
|
scale *= 10;
|
|
++k;
|
|
}
|
|
|
|
/* Perform any necessary adjustment of loerr and hierr to
|
|
take into account the formatting requirements. */
|
|
|
|
if (type == 'e')
|
|
cutoff += k; /* CutOffMode == "relative". */
|
|
/* Otherwise CutOffMode == "absolute". */
|
|
|
|
{ /* CutOffAdjust. */
|
|
int a = cutoff - k;
|
|
double y = scale;
|
|
while (a > 0)
|
|
{
|
|
y *= 10;
|
|
--a;
|
|
}
|
|
while (a < 0)
|
|
{
|
|
y = ceil (y / 10);
|
|
++a;
|
|
}
|
|
/* y == ceil (scale * pow (10.0, (double) (cutoff - k))) */
|
|
if (y > loerr)
|
|
loerr = y;
|
|
if (y > hierr)
|
|
{
|
|
hierr = y;
|
|
roundup = 1;
|
|
}
|
|
} /* End CutOffAdjust. */
|
|
} while (r2 + hierr >= 2 * scale);
|
|
|
|
/* End Fixup. */
|
|
|
|
/* First digit. */
|
|
--k;
|
|
r10 = r * 10;
|
|
digit = '0' + (unsigned int) floor (r10 / scale);
|
|
r = fmod (r10, scale);
|
|
loerr *= 10;
|
|
hierr *= 10;
|
|
|
|
low = 2 * r < loerr;
|
|
if (roundup)
|
|
high = 2 * r >= (2 * scale) - hierr;
|
|
else
|
|
high = 2 * r > (2 * scale) - hierr;
|
|
|
|
if (low || high || k == cutoff)
|
|
{
|
|
if ((high && !low) || (2 * r > scale))
|
|
++digit;
|
|
}
|
|
|
|
if (type == 'e')
|
|
{
|
|
/* Exponential notation. */
|
|
|
|
int expt = k; /* Base-10 exponent. */
|
|
|
|
/* Find the magnitude of the exponent. */
|
|
j = 1;
|
|
do
|
|
j *= 10;
|
|
while (j <= expt);
|
|
|
|
/* Write the first digit. */
|
|
put (digit);
|
|
|
|
if (low || high || k == cutoff)
|
|
{
|
|
if (prec > 0 || info->alt)
|
|
put (*decimal);
|
|
}
|
|
else
|
|
{
|
|
put (*decimal);
|
|
|
|
/* First post-decimal digit. */
|
|
--k;
|
|
r10 = r * 10;
|
|
digit = '0' + (unsigned int) floor (r10 / scale);
|
|
r = fmod (r10, scale);
|
|
loerr *= 10;
|
|
hierr *= 10;
|
|
|
|
low = 2 * r < loerr;
|
|
if (roundup)
|
|
high = 2 * r >= (2 * scale) - hierr;
|
|
else
|
|
high = 2 * r > (2 * scale) - hierr;
|
|
|
|
if (low || high || k == cutoff)
|
|
{
|
|
if ((high && !low) || (2 * r > scale))
|
|
++digit;
|
|
put (digit);
|
|
}
|
|
else
|
|
{
|
|
put (digit);
|
|
|
|
/* Remaining digits. */
|
|
while (1)
|
|
{
|
|
--k;
|
|
r10 = r * 10;
|
|
digit = '0' + (unsigned int) floor (r10 / scale);
|
|
r = fmod (r10, scale);
|
|
loerr *= 10;
|
|
hierr *= 10;
|
|
|
|
low = 2 * r < loerr;
|
|
if (roundup)
|
|
high = 2 * r >= (2 * scale) - hierr;
|
|
else
|
|
high = 2 * r > (2 * scale) - hierr;
|
|
|
|
if (low || high || k == cutoff)
|
|
{
|
|
if ((high && !low) || (2 * r > scale))
|
|
++digit;
|
|
put (digit);
|
|
break;
|
|
}
|
|
|
|
put (digit);
|
|
}
|
|
}
|
|
}
|
|
|
|
if (tolower (info->spec) != 'g' || info->alt)
|
|
/* Pad with zeros. */
|
|
while (k-- >= cutoff)
|
|
put ('0');
|
|
|
|
/* Write the exponent. */
|
|
put (isupper (info->spec) ? 'E' : 'e');
|
|
put (expt < 0 ? '-' : '+');
|
|
expt = abs (expt);
|
|
if (expt < 10)
|
|
/* Exponent always has at least two digits. */
|
|
put ('0');
|
|
do
|
|
{
|
|
j /= 10;
|
|
put ('0' + (expt / j));
|
|
expt %= j;
|
|
}
|
|
while (j > 1);
|
|
}
|
|
else
|
|
{
|
|
/* Decimal fraction notation. */
|
|
|
|
if (k < 0)
|
|
{
|
|
put ('0');
|
|
if (prec > 0 || info->alt)
|
|
put (*decimal);
|
|
}
|
|
|
|
/* Write leading fractional zeros. */
|
|
j = 0;
|
|
while (--j > k)
|
|
put ('0');
|
|
|
|
if (low || high || k == cutoff)
|
|
put (digit);
|
|
else
|
|
while (1)
|
|
{
|
|
put (digit);
|
|
|
|
--k;
|
|
digit = '0' + (unsigned int) floor ((r * 10) / scale);
|
|
r = fmod (r * 10, scale);
|
|
loerr *= 10;
|
|
hierr *= 10;
|
|
|
|
low = 2 * r < loerr;
|
|
if (roundup)
|
|
high = 2 * r >= (2 * scale) - hierr;
|
|
else
|
|
high = 2 * r > (2 * scale) - hierr;
|
|
|
|
if (low || high || k == cutoff)
|
|
{
|
|
if ((high && !low) || (2 * r > scale))
|
|
++digit;
|
|
put (digit);
|
|
break;
|
|
}
|
|
|
|
if (k == -1)
|
|
put (*decimal);
|
|
}
|
|
while (k > 0)
|
|
{
|
|
put ('0');
|
|
--k;
|
|
}
|
|
if (k == 0 && (prec > 0 || info->alt))
|
|
{
|
|
put (*decimal);
|
|
while (prec-- > 0)
|
|
put ('0');
|
|
}
|
|
}
|
|
}
|
|
|
|
#undef put
|
|
|
|
/* The number is all converted in BUF.
|
|
Now write it with sign and appropriate padding. */
|
|
|
|
if (is_neg || info->showsign || info->space)
|
|
--width;
|
|
|
|
width -= bp - buf;
|
|
|
|
if (!info->left && info->pad == ' ')
|
|
/* Pad with spaces on the left. */
|
|
while (width-- > 0)
|
|
outchar (' ');
|
|
|
|
/* Write the sign. */
|
|
if (is_neg)
|
|
outchar ('-');
|
|
else if (info->showsign)
|
|
outchar ('+');
|
|
else if (info->space)
|
|
outchar (' ');
|
|
|
|
if (!info->left && info->pad == '0')
|
|
/* Pad with zeros on the left. */
|
|
while (width-- > 0)
|
|
outchar ('0');
|
|
|
|
if (fwrite (buf, bp - buf, 1, s) != 1)
|
|
return -1;
|
|
done += bp - buf;
|
|
|
|
if (info->left)
|
|
/* Pad with spaces on the right. */
|
|
while (width-- > 0)
|
|
outchar (' ');
|
|
|
|
return done;
|
|
}
|