CODEWITHSUNDEEP
cfloating-point-precision
David Ranieri
David Ranieri

Reputation: 41046

Are multiples of 0.25 exactly representable as double?

I have the following code for finding quartiles:

#include <stdio.h>
#include <stdlib.h>

typedef struct {
    double qrt[3];
    double *value;
    int count;
} t_data;

static void set_qrt(t_data *data, int qrt)
{
    int n, e;
    double d;

    d = qrt * 0.25 * data->count + 0.5;
    n = (int)d;
    e = n != d;
    data->qrt[qrt - 1] = data->value[n - 1];
    if (e) {
        data->qrt[qrt - 1] += data->value[n];
        data->qrt[qrt - 1] *= 0.5;
    }
}

static void set_qrts(t_data *data)
{
    set_qrt(data, 2);
    if (data->count > 1) {
        set_qrt(data, 1);
        set_qrt(data, 3);
    } else {
        data->qrt[0] = 0.0;
        data->qrt[2] = 0.0;
    }
}

static int comp(const void *pa, const void *pb)
{
    const double a = *(const double *)pa;
    const double b = *(const double *)pb;

    return (a > b) ? 1 : (a < b) ? -1 : 0;
}

int main(void)
{
    double values[] = {3.7, 8.9, 7.1, 5.4, 1.2, 6.8, 4.3, 2.7};
    t_data data;

    data.value = values;
    data.count = (int)(sizeof(values) / sizeof(double));
    qsort(data.value, data.count, sizeof(double), comp);
    set_qrts(&data);
    printf("Q1 = %.1f\nQ2 = %.1f\nQ3 = %.1f\n", data.qrt[0], data.qrt[1], data.qrt[2]);
}

Is

d = qrt * 0.25 * data->count + 0.5;
n = (int)d;
e = n != d;

guaranteed to work as expected? (e == isinteger(d))

Upvotes: 6

Views: 1841

Answers (4)

Olof Forshell
Olof Forshell

Reputation: 3274

The double precision floating point format has 53 bits in its manitissa of which one is implicit. This means that it can represent all positive and negative integers in the range 2^0 to 2^53-1.

0 (zero) is a special case which has its own format.

When it comes to a 0.25 spacing the range is straight-forwardly calculated to be 2^-2 to 2^51-0.25. This means that quite a few but by no means all multiples of 0.25 are exactly representable in the double precision format, just as a quite a few but not all integers are exactly representable.

So if you have an exactly representable spacing of 2^x the representable range is 2^x to 2^(53+x)-2^x.

Upvotes: 1

Lee Daniel Crocker
Lee Daniel Crocker

Reputation: 13196

The values 0.5 and 0.25 themselves will be exact. The intermediate values of your calculation may or may not be, depending on their range. IEEE doubles have a 52-bit mantissa, so they will exactly represent to the 0.25 numbers that need 50 bits or fewer in the mantissa, which is about 15 decimal digits.

So if you add 0.25 to 100000000000000 (10^14), you'll get 100000000000000.25. But if you add 0.25 to 10000000000000000 (10^16), you'll lose the fraction.

Upvotes: 3

Barış Akkurt
Barış Akkurt

Reputation: 2256

dasblinkenlight is absolutely correct. Double/float and integer types are stored differently according to IEEE754. Watch this for an easy tutorial if you are curious about it.

Upvotes: 1

Sergey Kalinichenko
Sergey Kalinichenko

Reputation: 726929

Numbers 0.5, 0.25, 0.125 and so on represent negative powers of two, and therefore are representable exactly in IEEE 754 types. Using these numbers does not result in representation errors.

Upvotes: 8

Related Questions