This is the mail archive of the libc-hacker@sources.redhat.com mailing list for the glibc project.
Note that libc-hacker is a closed list. You may look at the archives of this list, but subscription and posting are not open.
Index Nav: | [Date Index] [Subject Index] [Author Index] [Thread Index] | |
---|---|---|
Message Nav: | [Date Prev] [Date Next] | [Thread Prev] [Thread Next] |
Other format: | [Raw text] |
Hi! As shown by the testcase below, ecvt/ecvt_r doesn't work on denormals (well, it works for double denormals on IA-32 unless -mfpmath=sse, as there if you are lucky f is computed in long double precision and not spilled until done with it). The problem is that a denormal needs to be multiplied by a power of 10 bigger than what can fit into the floating type (double resp. long double) if the result of the multiplication should be bigger or equal than 1.0. So current CVS code overflows f into +Inf and efcvt_r gives incorrect result. If ecvt/ecvt_r was some function where we would care a lot about speed, it would be perhaps good to also use some precomputed powers of ten (1e-4096, 1e-2048, ... 1e1, 1e2, 1e4, 1e8, 1e16, 1e32, 1e64, ... 1e4096), so that we don't for close to smallest normalized doubles multiply f by 10 ~ 300 times and for really small long doubles ~ 4900 times. But given that it is a legacy function, I'm not sure if it is ok to increase .rodata by 16 * sizeof (double) and 24 * sizeof (long double) and code size because of that. What do you think? 2004-12-21 Jakub Jelinek <jakub@redhat.com> * misc/efgcvt_r.c (FLOAT_MIN_10_EXP, FLOAT_MIN_10_NORM): Define. (ecvt_r): Special case denormals. * misc/qefgcvt_r.c (FLOAT_MIN_10_EXP, FLOAT_MIN_10_NORM): Define. * misc/tst-efgcvt.c: Include float.h. (ecvt_tests): Add 2 new tests. --- libc/misc/tst-efgcvt.c.jj 2001-07-06 06:55:36.000000000 +0200 +++ libc/misc/tst-efgcvt.c 2004-12-21 18:41:02.890224202 +0100 @@ -1,4 +1,4 @@ -/* Copyright (C) 1998, 1999, 2000 Free Software Foundation, Inc. +/* Copyright (C) 1998, 1999, 2000, 2004 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 @@ -20,6 +20,7 @@ # define _GNU_SOURCE 1 #endif +#include <float.h> #include <math.h> #include <stdio.h> #include <stdlib.h> @@ -59,6 +60,10 @@ static testcase ecvt_tests[] = { 123.01, -4, 3, "" }, { 126.71, -4, 3, "" }, { 0.0, 4, 1, "0000" }, +#if DBL_MANT_DIG == 53 + { 0x1p-1074, 3, -323, "494" }, + { -0x1p-1074, 3, -323, "494" }, +#endif /* -1.0 is end marker. */ { -1.0, 0, 0, "" } }; --- libc/misc/qefgcvt_r.c.jj 2004-05-07 14:32:16.000000000 +0200 +++ libc/misc/qefgcvt_r.c 2004-12-21 18:18:31.261313136 +0100 @@ -24,6 +24,7 @@ #define FUNC_PREFIX q #define FLOAT_FMT_FLAG "L" #define FLOAT_NAME_EXT l +#define FLOAT_MIN_10_EXP LDBL_MIN_10_EXP #if LDBL_MANT_DIG == 64 # define NDIGIT_MAX 21 #elif LDBL_MANT_DIG == 53 @@ -40,5 +41,16 @@ # error "NDIGIT_MAX must be precomputed" # define NDIGIT_MAX (lrint (ceil (M_LN2 / M_LN10 * LDBL_MANT_DIG + 1.0))) #endif +#if LDBL_MIN_10_EXP == -37 +# define FLOAT_MIN_10_NORM 1.0e-37L +#elif LDBL_MIN_10_EXP == -307 +# define FLOAT_MIN_10_NORM 1.0e-307L +#elif LDBL_MIN_10_EXP == -4931 +# define FLOAT_MIN_10_NORM 1.0e-4931L +#else +/* libc can't depend on libm. */ +# error "FLOAT_MIN_10_NORM must be precomputed" +# define FLOAT_MIN_10_NORM exp10l (LDBL_MIN_10_EXP) +#endif #include "efgcvt_r.c" --- libc/misc/efgcvt_r.c.jj 2004-08-04 14:42:22.000000000 +0200 +++ libc/misc/efgcvt_r.c 2004-12-21 18:31:52.486396147 +0100 @@ -1,5 +1,5 @@ /* Compatibility functions for floating point formatting, reentrant versions. - Copyright (C) 1995,96,97,98,99,2000,01,02 Free Software Foundation, Inc. + Copyright (C) 1995,96,97,98,99,2000,01,02,04 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 @@ -31,6 +31,7 @@ # define FUNC_PREFIX # define FLOAT_FMT_FLAG # define FLOAT_NAME_EXT +# define FLOAT_MIN_10_EXP DBL_MIN_10_EXP # if DBL_MANT_DIG == 53 # define NDIGIT_MAX 17 # elif DBL_MANT_DIG == 24 @@ -43,6 +44,17 @@ # error "NDIGIT_MAX must be precomputed" # define NDIGIT_MAX (lrint (ceil (M_LN2 / M_LN10 * DBL_MANT_DIG + 1.0))) # endif +# if DBL_MIN_10_EXP == -37 +# define FLOAT_MIN_10_NORM 1.0e-37 +# elif DBL_MIN_10_EXP == -307 +# define FLOAT_MIN_10_NORM 1.0e-307 +# elif DBL_MIN_10_EXP == -4931 +# define FLOAT_MIN_10_NORM 1.0e-4931 +# else +/* libc can't depend on libm. */ +# error "FLOAT_MIN_10_NORM must be precomputed" +# define FLOAT_MIN_10_NORM exp10 (DBL_MIN_10_EXP) +# endif #endif #define APPEND(a, b) APPEND2 (a, b) @@ -171,6 +183,17 @@ APPEND (FUNC_PREFIX, ecvt_r) (value, ndi d = -value; else d = value; + /* For denormalized numbers the d < 1.0 case below won't work, + as f can overflow to +Inf. */ + if (d < FLOAT_MIN_10_NORM) + { + value /= FLOAT_MIN_10_NORM; + if (value < 0.0) + d = -value; + else + d = value; + exponent += FLOAT_MIN_10_EXP; + } if (d < 1.0) { do Jakub
Index Nav: | [Date Index] [Subject Index] [Author Index] [Thread Index] | |
---|---|---|
Message Nav: | [Date Prev] [Date Next] | [Thread Prev] [Thread Next] |