c++numerical-stability

Numerical precision of comparison operators


Short question: is the behavior of comparison operators defined at max precision? Meaning, if I have two numbers (x and y) that are identical up to precision, should I always expect the same answer when I do x < y?

I realize this might seem like a dumb question, but let me elaborate. I am working with double and I have an array of numbers like:

0:  62536.5477752959
1:  62536.4840613718
2:  62536.4576412381
3:  62522.8487197062
4:  62536.5473896233
5:  62536.5467941254
6:  62527.3508907998
7:  62536.5477752959
8:  62517.5900098039
9:  62536.5477752959

Notice that entries 0, 7 and 9 have the same value.

When I do (something like this) :

int low = 0, high = 0;
for(int i = 0; i < N; ++i) {
  if(x[i] < x[low])
    low = i;
  if(x[i] > x[high]) 
    high = i;
 }
 cout << "low: " << low << " high: " << high << endl;

I sometimes get: low: 8 high: 0, sometimes: low: 8 high: 7 I would have expected always the lowest index value.

Any ideas?

[edit missing braces.]


Solution

  • Yes it is, assuming IEEE754 for your floating point types. Any two double values x and y say are such that exactly one of x < y, x == y, or x > y holds, with the exception of some edge cases like +Inf, -Inf, or NaN.

    Confusion starts to arise when you use decimal notation to represent floating point values; e.g. there is no such double as 62536.5477752959 (or any other one on your list for that matter).

    The numbers that you present have been truncated by your debugger / standard outputter, they are not the ones being used in the actual algorithm that you present. Be assured that the same decimal number always produces the same double, there is no arbitrary choice being made here: IEEE754 mandates that the closest double is picked.

    For further reading, see Is floating point math broken?

    Finally, replace int i with int i = 0. Currently the behaviour of your program is undefined.