What's the closest double to 1.0, that isn't 1.0?

Question

Is there a way to programmatically get the double that is closest to 1.0, but isn't actually 1.0?

One hacky way to do this would be to memcpy the double to a same-sized integer, and then subtract one. The way IEEE754 floating-point formats work, this would end up decreasing the exponent by one while changing the fractional part from all zeros (1.000000000000) to all ones (1.111111111111). However there exist machines where integers are stored little-endian while floating-point is stored big-endian, so that won't always work.

Calandra

You can't assume that +1 is the same distance (from 1.0) as -1. The interleaving of base 10 and base 2 floating point representations means that the gaps are uneven. – Aliphatic

@Richard: you're right. It is very unlikely that subtracting one ULP will get the, er, "nextbefore" value, because I guess the exponent would have to be adjusted as well. nextafter() is the only proper way of achieving what he wants. – Denman

FYI have a read of this blog (not mine): exploringbinary.com/… – Aliphatic

@RudyVelthuis: It does work on every IEEE754 binary floating point format. – Northway

@EdgarBonet: have you tried it? – Denman

@RudyVelthuis: No, I just believe what IEEE 754 says. If you don't, you should try: that's how you build confidence in your understanding of the spec. – Northway

Ok, then tell me: what "does work on every IEEE754 floating point format"? It is simply not true that if you decrement the significand that you get the "firstbefore()" value, especially not for 1.0, which has a significand which is a power of two. That means that 1.0000... binary is decrement to 0.111111.... and to normalize it, you must shift it to the left: 1.11111... which requires you to decrement the exponent. And then you are 2 ulp away from 1.0. So no, subtracting one from the integral value does NOT give you what is asked here. – Denman

@Edgar: You might want to show your confidence in your reading skills and quote that part of the IEEE 754 spec that says that "it does work on every IEEE754 binary floating point". – Denman

@RudyVelthuis: The actual standard is behind a paywall, and I will gladly read it for you if you buy me a copy. The data format, however, is available from many secondary sources. It is really common knowledge, so you have no excuse for telling me wrong before checking the facts. I'll just give you a hint: 1.0 = 0–01111111111–000...000, 1−ε/2 = 0–01111111110–111...111. I'll let you do the integral subtraction (it's not that hard). – Northway

@RudyVelthuis: For all positive, non-∞ FP values, the result of comparing the actual bit representation, as integer, of two FP numbers is the same as the result of comparing their FP values. Consequently, the closest two values to any positive floating-point value, other than DBL_MIN and DBL_MAX, can be computed by writing the FP value to memory, doing an integer increment/decrement, and then reading it back. At least, in IEEE754. Trust me, I've done this. Trouble is, you can't write PORTABLE code to do that because on some architectures, the endian-ness is different, between integers and FP. – Calandra

If you don't believe me, run the code: ideone.com/kewPwu – Calandra

Answer 1

Since C++11, you may use nextafter to get next representable value in given direction:

std::nextafter(1., 0.); // 0.99999999999999989
std::nextafter(1., 2.); // 1.0000000000000002

Demo

Reisinger

This is also a nice way to increment a double to the next representable integer: std::ceil(std::nextafter(1., std::numeric_limits<double>::max())). – Bettyebettzel

The next question's going to be "how is this implemented in the stdlib" :P – Selfdelusion

My man page for nextafter says "CONFORMING TO C99, POSIX.1-2001, POSIX.1-2008. This function is defined in IEC 559 (and the appendix with recommended functions in IEEE 754/IEEE 854)." So you may be able to get it via C even if you're 12 years behind C++11 :) – Shade

After reading @LightnessRacesinOrbit's comment I got curious. This is how glibc implements nextafter, and this is how musl implements it in case anyone else wants to see how it's done. Basically: raw bit twiddling. – Permeate

@Cornstalks: I am not surprised it's down to bit twiddling, the only other option would be having CPU support. – Thyrotoxicosis

Bit twiddling is the only way to do it properly, IMO. You could do a lot of test tries, trying to slowly approach it, but that could be very slow. – Denman

@Permeate For positive normal doubles it's just the bit pattern, cast to int, plus one. Which is pretty cool, actually. – Bewley

Answer 2

In C and C++, the following gives the closest value to 1.0:

#include <limits.h>

double closest_to_1 = 1.0 - DBL_EPSILON/FLT_RADIX;

Note however that in later versions of C++, limits.h is deprecated in favour of climits. But then, if you are using C++ specific code anyway, you can use

#include <limits>

typedef std::numeric_limits<double> lim_dbl;
double closest_to_1 = 1.0 - lim_dbl::epsilon()/lim_dbl::radix;

And as Jarod42 writes in his answer, since C99 or C++11 you can also use nextafter:

#include <math.h>

double closest_to_1 = nextafter(1.0, 0.0);

Of course in C++ you can (and for later C++ versions should) include cmath and use std::nextafter instead.

Hispid

Answer 3

In C, you can use this:

#include <float.h>
...
double value = 1.0+DBL_EPSILON;

DBL_EPSILON is the difference between 1 and the least value greater than 1 that is representable.

You'll need to print it to several digits in order to see the actual value.

On my platform, printf("%.16lf",1.0+DBL_EPSILON) gives 1.0000000000000002.

Worms

So that wont work for some other values than 1. as 1'000'000 Demo – Reisinger

@Jarod42: You're right, but OP asks specifically about 1.0. BTW, it also gives the closest value greater than 1, and not the absolute closest value to 1 (which is possibly smaller than 1). So I agree that this is a partial answer, but I thought it could nevertheless contribute. – Worms

@LưuVĩnhPhúc: I give precision about the restriction of the answer, and the closest in the other direction. – Reisinger

This does not give the closest double to 1.0, as (assuming base 2) the double right before 1.0 is only half as far away as the double right after 1.0 (which is the one you calculate). – Hispid

@celtschk: You're right, I've explained that in the comment above. – Worms

The closest value to one is indeed smaller than one, as @barakmanos sugests, since double values between 0.5 and 1.0 have twice the precision of values between 1.0 and 2.0. – Besmear

@HelloGoodbye: barakmanos is the same user who wrote this answer. – Worms

The question is very explicit: “closest double to 1.0, that isn't 1.0”, and this answer is just wrong. – Northway

How does this work at all since 1.0 is not representable with absolute accuracy? Or is 1.0 representable in memory with absolute accuracy? - Isn't that the point of the question, that 1.0 is not representable in memory with absolute accuracy? Does this even make sense? What? – Pit

@user3728501: Sorry, I do not understand your question. How does what work? – Worms

@barakmanos double value = 1.0 + DBL_EPSILON; - DBL_EPSILON has some binary value in memory, probably something like 0 00000000000 0000000000000000000000000000000000000000000000000001 (0 sign, 0 exponent, all fraction set to zero except last bit? - maybe this isn't quite right but you see the principle of the idea)... but essentially my question can be answered by considering is 1.0 representable exactly in floating point format?... At first I thought "no", however I now think "yes", and I think that binary representation I just gave is 1.0? - In which case what is the representation – Pit

of DBL_EPSILON in memory as a 64 bit float? – Pit

@user3728501: Why don't you give it a try and see exactly how it is represented in memory? Here is a piece of code that breaks strict aliasing rule but should still give you the answer: double x = 1.0; printf("%X",*(uint64_t*)&x);. – Worms

Could do that I suppose, but like everyone else in the world, we look for the fastest solution - ie a simple yes/no answer to the question is 1.0 absolutely representable in memory will do. Having thought about it further I conclude the answer is yes. – Pit

Answer 4

In C++ you can also use this

1 + std::numeric_limits<double>::epsilon()

Silencer

Like barak manos' answer, this will not work for any value other than 1. – Diecious

@Diecious tehnically for typical binary floating-point implementations it will work for any value between 1 and 2-epsilon. But, yes, you're right that it's only guaranteed to apply to 1. – Sophey

Technically, it doesn't work for 1, since the closest number to 1 is the number right before 1, not the one right after it. double's precision between 0.5 and 1 is twice as high as its precision between 1 and 2, hence the number right before 1 ends up closer to 1. – Besmear