Thu 15 Jun 2017 11:33:46 AM UTC, comment #3:
I did a bit of investigation into this, based on some data Jaap sent me offline.
The number Jaap is dealing with is in fact 0.9349999999999994981791929 (at least that's how GSL calculates it and printf ("%0.22f") displays it.)
I think the problem is in data-out.c - the function should_round_up is a bit too naive:
/* Returns true if the magnitude represented by R should be
rounded up when chopped off at DECIMALS decimal places, false
if it should be rounded down. */
static bool
should_round_up (const struct rounder *r, int decimals)
{
int digit = r->string[r->integer_digits + decimals + 1];
assert (digit >= '0' && digit <= '9');
return digit >= '5';
}
It simply looks as the next number after the desired decimal places - which in this case is '4' so concludes that rounding down is appropriate.
I suggest we change the implementation of should_round_up to iterate from the right hand side and carry the previous digit. Something like:
static bool
should_round_up (const struct rounder *r, int decimals)
{
int x;
int digit = 0;
int carry_digit = '0';
for (x = strlen (r->string) - 1;
x >= r->integer_digits + decimals + 1;
--x)
{
digit = r->string[x];
assert (digit >= '0' && digit <= '9');
if (carry_digit >= '5')
digit++;
carry_digit = digit;
}
return digit >= '5';
}
The comment for this function would need to be adjusted too, as would numerous test results.
|