%g should not print 0 as infe-308 (#185)
This commit is contained in:
parent
4335268a93
commit
e526e4f850
@ -490,21 +490,27 @@ static size_t _etoa(out_fct_type out, char *buffer, size_t idx, size_t maxlen, d
|
|||||||
} conv;
|
} conv;
|
||||||
|
|
||||||
conv.F = value;
|
conv.F = value;
|
||||||
int exp2 = (int) ((conv.U >> 52U) & 0x07FFU) - 1023; // effectively log2
|
int expval;
|
||||||
conv.U = (conv.U & ((1ULL << 52U) - 1U)) | (1023ULL << 52U); // drop the exponent so conv.F is now in [1,2)
|
if (conv.U) {
|
||||||
// now approximate log10 from the log2 integer part and an expansion of ln around 1.5
|
int exp2 = (int) ((conv.U >> 52U) & 0x07FFU) - 1023; // effectively log2
|
||||||
int expval = (int) (0.1760912590558 + exp2 * 0.301029995663981 + (conv.F - 1.5) * 0.289529654602168);
|
conv.U = (conv.U & ((1ULL << 52U) - 1U)) | (1023ULL << 52U); // drop the exponent so conv.F is now in [1,2)
|
||||||
// now we want to compute 10^expval but we want to be sure it won't overflow
|
// now approximate log10 from the log2 integer part and an expansion of ln around 1.5
|
||||||
exp2 = (int) (expval * 3.321928094887362 + 0.5);
|
expval = (int) (0.1760912590558 + exp2 * 0.301029995663981 + (conv.F - 1.5) * 0.289529654602168);
|
||||||
const double z = expval * 2.302585092994046 - exp2 * 0.6931471805599453;
|
// now we want to compute 10^expval but we want to be sure it won't overflow
|
||||||
const double z2 = z * z;
|
exp2 = (int) (expval * 3.321928094887362 + 0.5);
|
||||||
conv.U = (uint64_t) (exp2 + 1023) << 52U;
|
const double z = expval * 2.302585092994046 - exp2 * 0.6931471805599453;
|
||||||
// compute exp(z) using continued fractions, see https://en.wikipedia.org/wiki/Exponential_function#Continued_fractions_for_ex
|
const double z2 = z * z;
|
||||||
conv.F *= 1 + 2 * z / (2 - z + (z2 / (6 + (z2 / (10 + z2 / 14)))));
|
conv.U = (uint64_t) (exp2 + 1023) << 52U;
|
||||||
// correct for rounding errors
|
// compute exp(z) using continued fractions, see https://en.wikipedia.org/wiki/Exponential_function#Continued_fractions_for_ex
|
||||||
if (value < conv.F) {
|
conv.F *= 1 + 2 * z / (2 - z + (z2 / (6 + (z2 / (10 + z2 / 14)))));
|
||||||
expval--;
|
// correct for rounding errors
|
||||||
conv.F /= 10;
|
if (value < conv.F) {
|
||||||
|
expval--;
|
||||||
|
conv.F /= 10;
|
||||||
|
}
|
||||||
|
} else {
|
||||||
|
expval = 0;
|
||||||
|
conv.F = 1;
|
||||||
}
|
}
|
||||||
|
|
||||||
// the exponent format is "%+03d" and largest value is "307", so set aside 4-5 characters
|
// the exponent format is "%+03d" and largest value is "307", so set aside 4-5 characters
|
||||||
@ -513,7 +519,7 @@ static size_t _etoa(out_fct_type out, char *buffer, size_t idx, size_t maxlen, d
|
|||||||
// in "%g" mode, "prec" is the number of *significant figures* not decimals
|
// in "%g" mode, "prec" is the number of *significant figures* not decimals
|
||||||
if (flags & FLAGS_ADAPT_EXP) {
|
if (flags & FLAGS_ADAPT_EXP) {
|
||||||
// do we want to fall-back to "%f" mode?
|
// do we want to fall-back to "%f" mode?
|
||||||
if ((value >= 1e-4) && (value < 1e6)) {
|
if ((value == 0.0) || ((value >= 1e-4) && (value < 1e6))) {
|
||||||
if ((int) prec > expval) {
|
if ((int) prec > expval) {
|
||||||
prec = (unsigned) ((int) prec - expval - 1);
|
prec = (unsigned) ((int) prec - expval - 1);
|
||||||
} else {
|
} else {
|
||||||
|
Loading…
Reference in New Issue
Block a user