round double with N significant decimal digits in an overflow-safe manner

Viewed 116

I want a overflow-safe function that round a double like std::round in addition it can handle the number of significant decimal digts.

f.e.

round(-17.747, 2) -> -17.75
round(-9.97729, 2) -> -9.98
round(-5.62448, 2) -> -5.62
round(std::numeric_limits<double>::max(), 10) ...

My first attempt was

double round(double value, int precision)
{
    double factor=pow(10.0, precision);
    return floor(value*factor+0.5)/factor;
}

but this can easily overflow.

Assuming IEEE, it is possible to decrease the possibility of overflows, like this.

double round(double value, int precision)
{
    // assuming IEEE 754 with 64 bit representation
    // the number of significant digits varies between 15 and 17

    precision=std::min(17, precision);
    double factor=pow(10.0, precision);
    return floor(value*factor+0.5)/factor;
}

But this still can overflow.

Even this performance disaster does not work.

double round(double value, int precision)
{
    std::stringstream ss;
    ss << std::setprecision(precision) << value;
    std::string::size_type sz;
    return std::stod(ss.str(), &sz);
}
round(std::numeric_limits<double>::max(), 2.0) // throws std::out_of_range

Note:

  • I'm aware of setprecision, but i need rounding not only for displaying purpose. So that is not a solution.
  • Unlike this post here How to round a number to n decimal places in Java , my question is especially on overflow safety and in C++ (the anwser in the topic above are Java-specific or do not handle overflows)
1 Answers

I haven't heavily tested this code:

/* expects x in (-1, 1) */
double round_precision2(double x, int precision2) {
    double iptr, factor = std::exp2(precision2);
    double y = (x < 0) ? -x : x;
    std::modf(y * factor + .5, &iptr);
    return iptr/factor * ((x < 0) ? -1 : 1);
}

double round_precision(double x, int precision) {
    int bits = precision * M_LN10 / M_LN2;
            /* std::log2(std::pow(10., precision)); */
    double iptr, frac = std::modf(x, &iptr);
    return iptr + round_precision2(frac, bits);
}

The idea is to avoid overflow by only operating on the fractional part of the number.

We compute the number of binary bits to achieve the desired precision. You should be able to put a bound on them with the limits you describe in your question. Next, we extract the fractional and integer parts of the number. Then we add the integer part back to the rounded fractional part.

To compute the rounded fractional part, we compute the binary factor. Then we extract the integer part of the rounded number resulting from multiplying fractional part by the factor. Then we return the fraction by dividing the integral part by the factor.

Related