2

Important Edit: The original question was about getting the density of both doubles and fractions. As i get the answer for doubles and not for fractions, I'm changing the topic to close this question. The other half of the original question is here

New question

I want to find the density of doubles between 2 given numbers but I can't think of a good way. So I'm looking for a closed-form expressions doublesIn(a,b). Or some code that does the work in a reasonable time.

With doubles i should use some formula with mantissa and exponent I'm not aware of. I already have a code using nextafter and it's awfully slow close to [-1,1] (below 1e6 is very slow)

.

Any ideas? Thanks in advance! :)

PS: If you want to know, I'm coding some math stuff for myself and I want to find how useful would be to replace double with a fraction (long,long or similar) on certain algorithms (like Gaussian elimination, newton's method for finding roots, etc), and for that I want to have some measures.

Pernoctador
  • 108
  • 8
  • 1
    Try: http://www.exploringbinary.com/the-spacing-of-decimal-floating-point-numbers/ and http://www.exploringbinary.com/tag/floating-point/ – Richard Critten Dec 31 '17 at 20:42
  • 2
    The straightforward way is to use [std::nextafter](http://en.cppreference.com/w/cpp/numeric/math/nextafter). – Incomputable Dec 31 '17 at 20:44
  • 1
    All representible floating point values are rational (excluding infinities, NaNs, etc). – Peter Dec 31 '17 at 20:46
  • There are [std::nextafter and std::nexttoward](http://en.cppreference.com/w/cpp/numeric/math/nextafter) functions. It might take a while to go through all the values between the two `double`s. – Ron Dec 31 '17 at 20:48
  • Can you access the double's bit pattern as a long integer? – Patricia Shanahan Dec 31 '17 at 22:53
  • Patricia, it seems your question is answer here: https://stackoverflow.com/questions/18754174/separate-a-double-into-its-sign-exponent-and-mantissa?rq=1, though i haven't fully check it. Do you have an idea how to make use of it? – Pernoctador Dec 31 '17 at 23:15

1 Answers1

1

In what follows, including the program, I am assuming double is represented by IEEE 754 64-bit binary floating point. That is the most likely case, but not guaranteed by the C++ standard.

You can count doubles in a range in constant time, because you can count unsigned integers in a range in constant time by subtracting the start from the end and adjusting for whether the range is open or closed.

The doubles in a finite non-negative range have bit patterns that form a consecutive sequence of integers. For example, the range [1.0,2.0] contains one double for each integer in the range [0x3ff0_0000_0000_0000, 0x4000_0000_0000_0000].

Finite non-positive ranges of doubles behave the same way except the unsigned bit patterns increase in value as the doubles become more negative.

If your range includes both positive and negative numbers, split it at zero, so that you deal with one non-negative range and another non-positive range.

Most of the complications arise when you want to get the count exactly right. In that case, you need to adjust for whether the range is open or closed, and to count zero exactly once.

For your purpose, being off by one or two in a few hundred million may not matter much.

Here is a simple program that demonstrates the idea. It has received little error checking, so use at your own risk.

#include <iostream>
#include <cmath>
using namespace std;

uint64_t count(double start, double end);

void testit(uint64_t expected, double start, double end) {
    cout << hex << "Should be " << expected << ": " << count(start, end)
            << endl;
}

double increment(double data, int count) {
    int i;
    for (i = 0; i < count; i++) {
        data = nextafter(data, INFINITY);
    }
    return data;
}

double decrement(double data, int count) {
    int i;
    for (i = 0; i < count; i++) {
        data = nextafter(data, -INFINITY);
    }
    return data;
}

int main() {
    testit((uint64_t) 1 << 52, 1.0, 2.0);
    testit(5, 3.0, increment(3.0, 5));
    testit(2, decrement(0, 1), increment(0, 1));
    testit((uint64_t) 1 << 52, -2.0, -1.0);
    testit(1, -0.0, increment(0, 1));
    testit(10, decrement(0,10), -0.0);
    return 0;
}

// Return the bit pattern representing a double as
// a 64-bit unsigned integer.
uint64_t toInteger(double data) {
    return *reinterpret_cast<uint64_t *>(&data);
}

// Count the doubles in a range, assuming double
// is IEEE 754 64-bit binary.
// Counts [start,end), including start but excluding end
uint64_t count(double start, double end) {
    if (!(isfinite(start) && isfinite(end) && start <= end)) {
        // Insert real error handling here
        cerr << "error" << endl;
        return 0;
    }
    if (start < 0) {
        if (end < 0) {
            return count(fabs(end), fabs(start));
        } else if (end == 0) {
            return count(0, fabs(start));
        } else {
            return count(start, 0) + count(0, end);
        }
    }
    if (start == -0.0) {
        start = 0.0;
    }
    return toInteger(end) - toInteger(start);
}
Patricia Shanahan
  • 25,849
  • 4
  • 38
  • 75
  • Thanks for the answer. As you say, i don't need the exact amount of doubles and fractions in a range, only a good approximation to compare which can represent more numbers. The range also could only be non-negative, because density is the same between (2,3) and (-3,-2). However I don't see how I can get with this information a closed-form expression. And i guess counting each one is a bit slow even if it's linear, because there is a lot of representable numbers between the smallest positive integers. – Pernoctador Jan 02 '18 at 02:13
  • @Pernoctador Think about getting the number of integers in a range. Apply that to the result of `unsigned long L = *(unsigned long*) &d;` where `d` is your double. – Patricia Shanahan Jan 02 '18 at 02:50
  • Awesome! I have tested and I didn't find anything wrong. Thank you so much. Now i need to get a reasonable code or math expression for fractions – Pernoctador Jan 02 '18 at 21:16