Skip to content

Commit

Permalink
Fixes #114, #115: Special-casing powers close to the extremal ones.
Browse files Browse the repository at this point in the history
* Special-casing the power-of-10 function for 10^308 and 10^-308.
* Special-casing `get_normalized_components()` to not strive for better precision for near-extremal powers and settle for the "rougher" way to obtain the integral and decimal parts.
* Added a test-case of a near-extremal-power which fails due to the re-normalization even if the powering is correct.
  • Loading branch information
eyalroz committed Feb 5, 2022
1 parent 5b76bd6 commit 50cd7a8
Show file tree
Hide file tree
Showing 2 changed files with 22 additions and 8 deletions.
17 changes: 14 additions & 3 deletions src/printf/printf.c
Original file line number Diff line number Diff line change
Expand Up @@ -650,11 +650,20 @@ static struct scaling_factor update_normalization(struct scaling_factor sf, doub
return result;
}

static struct double_components get_normalized_components(bool negative, printf_size_t precision, double non_normalized, struct scaling_factor normalization)
static struct double_components get_normalized_components(bool negative, printf_size_t precision, double non_normalized, struct scaling_factor normalization, int floored_exp10)
{
struct double_components components;
components.is_negative = negative;
components.integral = (int_fast64_t) apply_scaling(non_normalized, normalization);
double scaled = apply_scaling(non_normalized, normalization);

bool close_to_representation_extremum = ( (-floored_exp10 + (int) precision) >= DBL_MAX_10_EXP - 1 );
if (close_to_representation_extremum) {
// We can't have a normalization factor which also accounts for the precision, i.e. moves
// some decimal digits into the mantissa, since it's unrepresentable, or nearly unrepresentable.
// So, we'll give up early on getting extra precision...
return get_components(negative ? -scaled : scaled, precision);
}
components.integral = (int_fast64_t) scaled;
double remainder = non_normalized - unapply_scaling((double) components.integral, normalization);
double prec_power_of_10 = powers_of_10[precision];
struct scaling_factor account_for_precision = update_normalization(normalization, prec_power_of_10);
Expand Down Expand Up @@ -819,6 +828,8 @@ static double log10_of_positive(double positive_number)

static double pow10_of_int(int floored_exp10)
{
if (floored_exp10 == 308) { return 1e308; }
if (floored_exp10 == -308) { return 1e-308; }
// Compute 10^(floored_exp10) but (try to) make sure that doesn't overflow
double_with_bit_access dwba;
int exp2 = bastardized_floor(floored_exp10 * 3.321928094887362 + 0.5);
Expand Down Expand Up @@ -887,7 +898,7 @@ static void print_exponential_number(output_gadget_t* output, double number, pri
struct double_components decimal_part_components =
should_skip_normalization ?
get_components(negative ? -abs_number : abs_number, precision) :
get_normalized_components(negative, precision, abs_number, normalization);
get_normalized_components(negative, precision, abs_number, normalization, floored_exp10);

// Account for roll-over, e.g. rounding from 9.99 to 100.0 - which effects
// the exponent and may require additional tweaking of the parts
Expand Down
13 changes: 8 additions & 5 deletions test/test_suite.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -888,11 +888,11 @@ TEST_CASE("tiny floating-point values", "[]" ) {
PRINTING_CHECK("1.380651569e-23", ==, sprintf_, buffer, "%.9e", 1.380651569e-23);
PRINTING_CHECK("1.3806515690e-23", ==, sprintf_, buffer, "%.10e", 1.380651569e-23);
PRINTING_CHECK("1.38065156900e-23", ==, sprintf_, buffer, "%.11e", 1.380651569e-23);
// PRINTING_CHECK("1.380651569000e-23", ==, sprintf_, buffer, "%.12e", 1.380651569e-23);
// PRINTING_CHECK("1.3806515690000e-23", ==, sprintf_, buffer, "%.13e", 1.380651569e-23);
// PRINTING_CHECK("1.38065156900000e-23", ==, sprintf_, buffer, "%.14e", 1.380651569e-23);
// PRINTING_CHECK("1.380651569000000e-23", ==, sprintf_, buffer, "%.15e", 1.380651569e-23);
// PRINTING_CHECK("1.3806515690000000e-23", ==, sprintf_, buffer, "%.16e", 1.380651569e-23);
PRINTING_CHECK("1.380651569000e-23", ==, sprintf_, buffer, "%.12e", 1.380651569e-23);
PRINTING_CHECK("1.3806515690000e-23", ==, sprintf_, buffer, "%.13e", 1.380651569e-23);
PRINTING_CHECK("1.38065156900000e-23", ==, sprintf_, buffer, "%.14e", 1.380651569e-23);
PRINTING_CHECK("1.380651569000000e-23", ==, sprintf_, buffer, "%.15e", 1.380651569e-23);
PRINTING_CHECK("1.3806515690000000e-23", ==, sprintf_, buffer, "%.16e", 1.380651569e-23);
#endif
}

Expand Down Expand Up @@ -1012,6 +1012,9 @@ TEST_CASE("floating-point specifiers, precision and flags", "[]" ) {
PRINTING_CHECK("-4e+04", ==, sprintf_, buffer, "%.1g", -40661.5);
PRINTING_CHECK("-4.e+04", ==, sprintf_, buffer, "%#.1g", -40661.5);
PRINTING_CHECK("100.", ==, sprintf_, buffer, "%#.3g", 99.998580932617187500);
// Ensure high precision does make us exceed the representation threshold
// internally (e.g. with normalization factors and such
PRINTING_CHECK("1.2345678901e-308", ==, sprintf_, buffer, "%.10e", 1.2345678901e-308);
// Rounding-focused checks
PRINTING_CHECK("4.895512e+04", ==, sprintf_, buffer, "%e", 48955.125);
PRINTING_CHECK("9.2524e+04", ==, sprintf_, buffer, "%.4e", 92523.5);
Expand Down

0 comments on commit 50cd7a8

Please sign in to comment.