In-depth: Comparing floating point numbers, 2012 edition

[In this reprinted #altdevblogaday in-depth piece, Valve Software programmer Bruce Dawson continues his series on the floating point format by looking at comparing floating point numbers.] We've finally reached the point in this series that I've been waiting for. In this post, I am going to share the most crucial piece of floating-point math knowledge that I have. Here it is:

You just won't believe how vastly, hugely, mind-bogglingly hard it is. I mean, you may think it's difficult to calculate when trains from Chicago and Los Angeles will collide, but that's just peanuts to floating-point math. Seriously. Each time I think that I've wrapped my head around the subtleties and implications of floating-point math, I find that I'm wrong and that there is some extra confounding factor that I had failed to consider. So, the lesson to remember is that floating-point math is always more complex than you think it is. Keep that in mind through the rest of the post where we talk about the promised topic of comparing floats, and understand that this post gives some suggestions on techniques, but no silver bullets. Previously on this channel… This is the fifth chapter in what is currently a four chapter series. The previous posts include:

Comparing for equality Floating point math is not exact. Simple values like 0.1 cannot be precisely represented using binary floating point numbers, and the limited precision of floating point numbers means that slight changes in the order of operations or the precision of intermediates can change the result. That means that comparing two floats to see if they are equal is usually not what you want. GCC even has a warning for this: "warning: comparing floating point with == or != is unsafe". Here's one example of the inexactness that can creep in:

  float f = 0.1f;
  float sum;
  sum = 0;

  for (int i = 0; i < 10; ++i)
    sum += f;
  float product = f * 10;
  printf("sum = %1.15f, mul = %1.15f, mul2 = %1.15fn",
      sum, product, f * 10);

This code shows two ancient and mystical techniques for calculating 'one'. I call these techniques "iterative-adding" and "multiplying". And, just to show how messy this stuff is I do the multiplication twice. That's three separate calculations of the same thing. Naturally we get three different results, and only one of them is the unity we were seeking:

Disclaimer: the results you get will depend on your compiler and your compiler settings, which actually helps make the point. So what happened, and which one is correct? What do you mean 'correct'? Before we can continue I need to make clear the difference between 0.1, float(0.1), and double(0.1). In C/C++ 0.1 and double(0.1) are the same thing, but when I say "0.1" in text I mean the exact base-10 number, whereas float(0.1) and double(0.1) are rounded versions of 0.1. And, to be clear, float(0.1) and double(0.1) don't have the same value, because float(0.1) has fewer binary digits, and therefore has more error. Here are the values for 0.1, float(0.1), and double(0.1):