Check if double are given with a specific accuracy

Question

I have a vector of double. Those numbers are given with a precision such as 0.001 meaning that the numbers are given with a thousandth of units accuracy. For example I'm expecting to have 123456789.012 but not 123456789.01231. In practice, because of floating-point arithmetic the actual numbers are more like 123456789.01199999452 which is the valid and closest representation of 123456789.012.

I'd like to robustly test if my numbers are actually given with the accuracy there are supposed to be. For example I'd like to warn if I have the 123456789.01231 which does not match with the accuracy 0.001. Another valid option to me is to find the accuracy of a given number.

What I've done so far is to test if (x - offset)/accuracy is an integer. The offset protects against integer overflow. In practice, because of floating-point arithmetic, it cannot be an integer so I added an arbitrary threshold. It works but it is not a robust solution in my opinion.

dx = (x[i]-offset)/accuracy;
ix = std::round(dx);
if (std::abs(dx - ix) > 1e-5)
  throw exception("...");

No string is involved in this process. The numbers come from binary files or from runtime computations. The accuracy is given by the user and the numbers are expected to be compliant with this accuracy. In my previous example of 123456789.012 I do know that this number does not actually exist but 123456789.01199999452 is valid because it is its best representation with double. But 123456789.01231 (actually 123456789.01230999827) is not valid because it is not the best representation of the 3 digits rounded value.

Answer 1

The problem stated in the question is suggestive of an attempt to use floating-point in ways it was not designed for, but there is insufficient information given to suggest alternatives. So, with some reluctance, here is a solution for the problem as stated.

/*  IsNearestValue(x, A) returns true iff x is the double number nearest to the
    nearest multiple of A', where A' is the unit fraction nearest A.  (All
    negative powers of 10, such as 10**-3, are unit fractions.)

    For example, IsNearestValue(x, .001) returns true iff x is the result of
    converting some number with three decimal digits after the decimal point
    to double.

    This routine presumes x/A <= 2**53.
*/
bool IsNearestValue(double x, double A)
{
    //  Calculate 1/A'.  The result has no rounding error.
    double reciprocal = std::round(1/A);

    /*  Find what multiple of A' x must be near.  This produces an exact
        result.  That is, t is an integer such that t * A' = x, with
        real-number arithmetic, not floating-point arithmetic.
    */
    double t = std::round(x * reciprocal);

    //  Calculate the double nearest t/A'.
    t /= reciprocal;

    //  Return true iff x is the double nearest t/A'.
    return x == t;
}

The concept here is pretty simple. First, we correct the problem that A is given as a double, but none of the desired numbers (.1, .01, .001) can be represented in a double. However, their reciprocals can, so we take the reciprocal and round to get exactly the reciprocal of the desired number.

Then we multiply x by the reciprocal and round to the nearest integer. Then we divide that integer by the reciprocal, and that gives us the double that x must be if it is indeed the double nearest some multiply of the desired A.

I am not sure the constraint x/A ≤ 2⁵³ is necessary, but I have not tried to prove it is not, so I am leaving that bound unless there is further request.

Answer 2

It seems you want to ensure that your numbers are as close as possible to a set of numbers, which are an arithmetic progression from 0 to infinity with common difference of the accuracy value. In this case, this is basically a quantization task with further evaluation of the quantization error.

#include <cmath>
#include <iostream>

bool is_number_accurate(double n, double accuracy)
{
    bool accurate = false;

    std::cout << "n: " << n << ", ";

    n = std::abs(n);
    if (n > accuracy)
    {
        // Epsilon to account for precision loss from the calculations below,
        // as well as floating point representation error.
        const double epsilon = 0.00001;

        // Remainder from quantizing n to accuracy precision
        // (i.e. quantization error)
        double rem = std::fmod(n, accuracy);

        // Normalize to 1.0. Here 0 means n exactly matches the quantized value,
        // and slightly greater than 0 means that n is slightly greater than the
        // quantized value.
        // Values close to 1.0 also mean n is close to match the quantized value,
        // but is slightly less than it.
        double error = rem / accuracy;
        std::cout << "error: " << error << ", ";

        accurate = error < epsilon || (error <= 1.0 && error > (1.0 - epsilon));
    }
    else
    {
        accurate = n == 0.0 || n == accuracy;
    }

    std::cout << "accurate: " << accurate << std::endl;

    return accurate;
}

int main()
{
    is_number_accurate(0.0, 0.001);
    is_number_accurate(0.0000000000001, 0.001);
    is_number_accurate(0.001, 0.001);
    is_number_accurate(0.01, 0.001);
    is_number_accurate(123456789.012, 0.001);
    is_number_accurate(123456789.0121, 0.001);
    is_number_accurate(123456789.0125, 0.001);
    is_number_accurate(123456789.0129, 0.001);
    is_number_accurate(123456789.013, 0.001);
    is_number_accurate(123456789.01231, 0.001);
}

The output is:

n: 0, accurate: 1
n: 1e-13, accurate: 0
n: 0.001, accurate: 1
n: 0.01, error: 0, accurate: 1
n: 1.23457e+08, error: 0.999992, accurate: 1
n: 1.23457e+08, error: 0.0999936, accurate: 0
n: 1.23457e+08, error: 0.5, accurate: 0
n: 1.23457e+08, error: 0.899992, accurate: 0
n: 1.23457e+08, error: 0.999994, accurate: 1
n: 1.23457e+08, error: 0.309996, accurate: 0