20

I am trying to determine the double machine epsilon in Java, using the definition of it being the smallest representable double value x such that 1.0 + x != 1.0, just as in C/C++. According to wikipedia, this machine epsilon is equal to 2^-52 (with 52 being the number of double mantissa bits - 1).

My implementation uses the Math.ulp() function:

double eps = Math.ulp(1.0);
System.out.println("eps = " + eps);
System.out.println("eps == 2^-52? " + (eps == Math.pow(2, -52)));

and the results are what I expected:

eps = 2.220446049250313E-16
eps == 2^-52? true

So far, so good. However, if I check that the given eps is indeed the smallest x such that 1.0 + x != 1.0, there seems to be a smaller one, aka the previous double value according to Math.nextAfter():

double epsPred = Math.nextAfter(eps, Double.NEGATIVE_INFINITY);
System.out.println("epsPred = " + epsPred);
System.out.println("epsPred < eps? " + (epsPred < eps));
System.out.println("1.0 + epsPred == 1.0? " + (1.0 + epsPred == 1.0));

Which yields:

epsPred = 2.2204460492503128E-16
epsPred < eps? true
1.0 + epsPred == 1.0? false

As we see, we have a smaller than machine epsilon such which, added to 1, yields not 1, in contradiction to the definition.

So what is wrong with the commonly accepted value for machine epsilon according to this definition? Or did I miss something? I suspect another esoteric aspect of floating-point maths, but I can't see where I went wrong...

EDIT: Thanks to the commenters, I finally got it. I actually used the wrong definition! eps = Math.ulp(1.0) computes the distance to the smallest representable double > 1.0, but -- and that's the point -- that eps is not the smallest x with 1.0 + x != 1.0, but rather about twice that value: Adding 1.0 + Math.nextAfter(eps/2) is rounded up to 1.0 + eps.

2
  • 1
    Have you tried with strictfp?
    – davmac
    Commented Feb 26, 2015 at 13:14
  • Yes, strictfp didn't help here.
    – Franz D.
    Commented Feb 26, 2015 at 15:49

3 Answers 3

21

using the definition of it being the smallest representable double value x such that 1.0 + x != 1.0, just as in C/C++

This has never been the definition, not in Java and not in C and not in C++.

The definition is that the machine epsilon is the distance between one and the smallest float/double larger than one.

Your “definition” is wrong by a factor of nearly 2.

Also, the absence of strictfp only allows a larger exponent range and should not have any impact on the empirical measurement of epsilon, since that is computed from 1.0 and its successor, each of which and the difference of which can be represented with the standard exponent range.

7
  • Is 1.0 + x != 1.0 && (1.0 + x) - 1.0 == x an adequate definition? Or "1.0 minus (the smallest double value greater than 1.0)"?
    – Random832
    Commented Feb 26, 2015 at 14:08
  • 2
    @Random832 If you want to “compute” it, use Math.nextAfter(1.0, Double.POSITIVE_INFINITY) - 1.0. The first proposal is a property that is true for many numbers including 1.0 (this is true: 1.0 + 1.0 != 1.0 && (1.0 + 1.0) - 1.0 == 1.0) Commented Feb 26, 2015 at 14:15
  • So the proposed definition (1.0 + x != 1.0) is invalid because of rounding? Is that the case?
    – Brian J
    Commented Feb 26, 2015 at 15:34
  • @BrianJ Yes. 1.0 + x is generally inexact for values of x of the magnitude we are concerned with here. By contrast, though the definition by subtraction leaves this fact implicit, the distance between 0x1.0000000000001 and 0x1.0000000000000 can be represented as a double and 0x1.0000000000001 - 0x1.0000000000000, as a subtraction of close floating-point numbers, is exact and computes exactly this value. Commented Feb 26, 2015 at 15:38
  • 1
    @FranzD. The wikipedia page claims that the C standard uses that definition. It does not. I actually cite the C standard in my blog post. Please do check the reference for yourself: open-std.org/jtc1/sc22/wg14/www/docs/n1256.pdf (the sentence to look for is given in the blog post). The wrong definition gives the wrong result, as anyone can check for themselves (and you have). The correct definition is used in glibc. What else do you need? Commented Feb 26, 2015 at 16:06
7

I'm not sure your experimental method / theory is sound. The documentation for the Math class states:

For a given floating-point format, an ulp of a specific real number value is the distance between the two floating-point values bracketing that numerical value

The documentation for the ulp method says:

An ulp of a double value is the positive distance between this floating-point value and the double value next larger in magnitude

So, if you want the smallest eps value such that 1.0 + eps != 1.0, your eps really should generally be less than Math.ulp(1.0), since at least for any value greater than Math.ulp(1.0) / 2, the result will be rounded up.

I think the smallest such value will be given by Math.nextAfter(eps/2, 1.0).

8
  • The mathematical ulp(x) is defined for any real x, so for most values of x, it doesn't matter whether the definition has ≤ or <. It is not generally agreed on what ulp(x) means when x is a power of two (the case for which it would matter whether it's ≤ < or < ≤). There is a quite thorough treatment in ens-lyon.fr/LIP/Pub/Rapports/RR/RR2005/RR2005-09.pdf Commented Feb 26, 2015 at 13:38
  • @davmac: I'm also not sure that my approach is sound, that's why I asked :) But I don't understand: Why isn't ulp(1.0), "the positive distance between 1.0 and the double next larger in magnitude", the smallest value such that 1.0 + eps != 1.0? There shouldn't be any smaller value with that property, as then ulp(1.0) isn't the distance to the next larger double. But, as @Pascal said, there doesn't seem to be a general definition of ulp(1.0), so maybe my problem lies here. Still I don't get it...
    – Franz D.
    Commented Feb 26, 2015 at 16:02
  • @FranzD. The value given by the “definition” as “smallest representable double value x such that 1.0 + x != 1.0” is not a power of two, so it is not the result of the ulp function for any definition of the mathematical ulp, nor for the Java method named ulp. Commented Feb 26, 2015 at 16:10
  • @Pascal: Huh? Who said ulp(x) returns only powers of two? I only said that the input x = 1.0 is a power of two, so ulp's results may be dubious.
    – Franz D.
    Commented Feb 26, 2015 at 16:23
  • 1
    @FranzD. I said that ulp returns powers of two. You do not have to believe me. Did you read ens-lyon.fr/LIP/Pub/Rapports/RR/RR2005/RR2005-09.pdf ? Commented Feb 26, 2015 at 16:30
0

In fairness to the OP, the definition he gives for machine epsilon is the same one I learned in Numerical Methods lo these many years ago. According to Wikipedia (https://en.wikipedia.org/wiki/Machine_epsilon), as of the date of this comment, both definitions are in widespread use. The two definitions may be either equivalent, or may differ by a factor of two, depending on whether the computer is using rounding or truncation.

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.