The Approximation Comes Due

Chapter 1 said that a double stores a very close approximation of the number you wrote rather than always the number itself, and promised that a later chapter would deal with the consequences. This is that chapter. Nothing about double changes here; what changes is that you stop trusting it in the one place where trusting it quietly breaks programs. Here is the demonstration nearly everyone meets eventually, usually while drafting a bug report against their compiler.

#include <stdio.h>

int main(void)
{
    double sum = 0.1 + 0.2;

    printf("0.1 + 0.2 shown as %%.2f is %.2f\n", sum);
    printf("0.1 + 0.2 shown as %%.17g is %.17g\n", sum);
    printf("0.3 alone shown as %%.17g is %.17g\n", 0.3);
    printf("sum == 0.3 is %d\n", sum == 0.3);

    return 0;
}

0.1 + 0.2 shown as %.2f is 0.30
0.1 + 0.2 shown as %.17g is 0.30000000000000004
0.3 alone shown as %.17g is 0.29999999999999999
sum == 0.3 is 0

Read it from the bottom. The last line is the one that matters: two values every calculator agrees are the same are not the same to C. The line above shows why. The literal 0.3 is stored as a value slightly below three tenths, 0.1 + 0.2 computes to one slightly above it, those are different bit patterns, and == correctly reports that they differ. The first line explains how you can write a lot of C without ever noticing, since %.2f rounds to two decimal places and both of these round to 0.30. The error was there all along and the output format was hiding it. %.17g is the specifier that stops hiding it, seventeen significant digits being enough to tell any two double values apart.

Why Binary Fractions Miss

This is neither a compiler bug nor sloppiness. It is the same thing that happens to you in decimal. Write one third in decimal and you get 0.3333 running on forever; stop at any digit and you have an approximation, and no number of digits ever makes it exact, because a decimal expansion can only spell fractions whose denominators are built out of 2 and 5. One third is not one of them. A double spells numbers as binary fractions. The bits after the point are worth one half, one quarter, one eighth and so on, and the format gives you 52 of them alongside a sign bit and an 11 bit exponent saying where the point sits. Binary can only spell fractions whose denominators are powers of two, and one tenth is not one of those either, so 0.1 in binary is 0.0001100110011 repeating forever exactly as one third repeats in decimal. The 52 bits fill up, the rest is dropped, and what gets stored is the nearest double to one tenth rather than one tenth itself. Add two such approximations and you get a third approximation, with no reason at all for it to be the same approximation you get by rounding 0.3 directly. None of this is specific to C: it is IEEE 754 binary floating point, which nearly every language and nearly every processor uses, so that same 0.30000000000000004 comes out of Python, JavaScript and Java.

What Is Exact, and Why double Is the Default

A great deal is exact, and knowing which parts keeps this from turning into superstition.

#include <stdio.h>

int main(void)
{
    float narrow = 0.1f;

    printf("0.5 is %.17g and 87.5 is %.17g\n", 0.5, 87.5);
    printf("0.1 is %.17g\n", 0.1);
    printf("0.1 stored in a float is %.17g\n", (double)narrow);

    return 0;
}

0.5 is 0.5 and 87.5 is 87.5
0.1 is 0.10000000000000001
0.1 stored in a float is 0.10000000149011612

0.5 is one half, and 0.25 is one quarter, both powers of two and both stored perfectly. 87.5 is 87 plus one half, so its whole part and its fraction are each representable and the value is exact, which is why the %.1f averages and prices earlier in this course never surprised anybody. Whole numbers are exact too, up to 2 to the 53rd, which is 9007199254740992: below that every integer gets its own double, and above it they begin to share, so adding 1.0 to 9007199254740992.0 gives back 9007199254740992 unchanged. A double counting items or cents is not approximating anything until the count grows astronomical. The second line of output is the exception that started the lesson, and the third is that same failure in a float, roughly a hundred million times larger.

That third line also settles which type to reach for. A float is 4 bytes against a double's 8, and its 23 fraction bits give it roughly 7 significant decimal digits where a double manages about 15. Use double by default. Reach for float when memory or bandwidth genuinely demands it, in a large array or a graphics buffer, which is almost never the case in a program written while learning. One asymmetry is worth stating once: printf promotes a float argument to double so %f prints either, while scanf is handed an address and can promote nothing, so it needs %f for a float * and %lf for a double *.

Compare With a Tolerance

The rule follows directly and it has no exceptions in this course. Never compare computed floating point values with ==. Ask instead whether they are close enough, which means subtracting them, taking the size of that difference regardless of sign, and testing it against a small number you chose deliberately. fabs returns the absolute value of a double and lives in <math.h>, C's mathematics header, which you have not needed until now. Nothing extra is required to use it here, since the compile command this course runs links the maths library for you.

#include <math.h>
#include <stdio.h>

static const double TOLERANCE = 1e-9;

static int nearly_equal(double a, double b, double tolerance)
{
    return fabs(a - b) < tolerance;
}

int main(void)
{
    double sum = 0.1 + 0.2;

    printf("the gap is %.17g\n", fabs(sum - 0.3));
    printf("nearly_equal(sum, 0.3, TOLERANCE) is %d\n", nearly_equal(sum, 0.3, TOLERANCE));

    return 0;
}

the gap is 5.5511151231257827e-17
nearly_equal(sum, 0.3, TOLERANCE) is 1

The gap is about 5.6e-17, vastly smaller than the 1e-9 we have agreed to ignore, so the tolerant comparison says yes where == said no. nearly_equal is one line of body and worth writing once per program rather than spelling fabs(a - b) < 1e-9 at every call site, because the name records what the comparison means. 1e-9 is scientific notation for 0.000000001, and TOLERANCE is a const double rather than a bare literal so the number has a name and exactly one place to change. Now the honest part, because a fixed tolerance is a tool with a working range rather than a law. 1e-9 is a good answer for values near 1 and a poor one at both ends of the scale. Up near 1e16 consecutive doubles are already several units apart, so 1e16 and 1e16 + 4.0 differ by 4 and a 1e-9 tolerance calls them different even though they are about as close as doubles get at that size. Down at 1e-12 the opposite happens: 1e-12 and 3e-12 differ by a factor of three, their gap of roughly 2e-12 sits comfortably inside 1e-9, and the same test calls them equal. The general answer is a relative comparison that scales the tolerance to the size of the values involved, which has real subtleties of its own around zero and is worth meeting when you need it. At this course's scale the working rule is simpler: choose a tolerance suited to the magnitudes your program actually handles, give it a name, and remember that it is a choice you made rather than a constant of nature.

The Error Accumulates

One rounding is invisible. Repeating it is how the invisible becomes a bug.

#include <math.h>
#include <stdio.h>

int main(void)
{
    double total = 0.0;

    for (int i = 0; i < 10; ++i)
    {
        total += 0.1;
    }

    printf("ten additions of 0.1 give %.17g\n", total);
    printf("the same value with %%.2f is %.2f\n", total);
    printf("total == 1.0 is %d, and the gap is %.17g\n", total == 1.0, fabs(total - 1.0));

    return 0;
}

ten additions of 0.1 give 0.99999999999999989
the same value with %.2f is 1.00
total == 1.0 is 0, and the gap is 1.1102230246251565e-16

Ten additions of a value that was already slightly wrong give a total wrong by about 1.1e-16, and the direction is not predictable from the input: this total landed just below 1.0 while 0.1 + 0.2 landed just above 0.3. Printed with %.2f it reads 1.00, which is exactly the trap, because the output looked right and the comparison still failed. Accumulate a few million values and the drift grows large enough to show in the printed answer too. A loop written as while (total != 1.0) over a total built by adding fractions is a loop that may never finish, and nearly_equal is the fix there as well.

Infinity, NaN, and What Is Not Undefined

Chapter 2 was blunt that dividing an integer by zero is undefined behaviour. Dividing a double by zero is a different thing entirely, and the difference is worth being exact about.

#include <stdio.h>

int main(void)
{
    double zero = 0.0;
    double huge = 1.0 / zero;
    double undefined_ratio = zero / zero;

    printf("1.0 / 0.0 is %f\n", huge);
    printf("0.0 / 0.0 is %f\n", undefined_ratio);
    printf("undefined_ratio == undefined_ratio is %d\n", undefined_ratio == undefined_ratio);

    return 0;
}

1.0 / 0.0 is inf
0.0 / 0.0 is nan
undefined_ratio == undefined_ratio is 0

IEEE 754 defines both results and the double format reserves bit patterns for them. Dividing a nonzero value by zero gives an infinity, which prints as inf and compares as you would hope, being greater than every finite value. Dividing zero by zero has no defensible answer, so the result is a NaN, short for "not a number", which prints as nan and spreads: arithmetic involving a NaN produces a NaN, so one bad division early on turns an entire calculation into nan rather than into a merely wrong number. The last line is the property people trip over. A NaN does not equal itself. IEEE 754 defines every comparison involving a NaN as false, so x == x is 0 when x is a NaN, and isnan(x) from <math.h> is how you actually ask the question. Hold the contrast in one sentence: 1 / 0 on integers is undefined behaviour and the standard then promises nothing about your program at all, while 1.0 / 0.0 is a defined operation with a defined result you can print, compare and test for. Neither is something to do on purpose, and only one of them is a hole in the language.

Key Takeaways

  • A double stores the nearest representable binary fraction, and most decimal fractions are not representable. 0.1 in binary repeats forever exactly as one third repeats in decimal, so 0.1 + 0.2 is 0.30000000000000004 while the literal 0.3 is 0.29999999999999999, and == between them is 0. This is IEEE 754, not a C defect. %.17g shows enough digits to distinguish any two double values, while formats like %.2f round the error out of sight, which is how the bug hides for a long time before it bites.
  • Powers of two are exact, so 0.5, 0.25 and 87.5 are stored perfectly, and whole numbers are exact up to 2 to the 53rd, which is 9007199254740992.
  • Never compare computed floating point values with ==. Write fabs(a - b) < tolerance with fabs from <math.h>, wrapped in a named function such as nearly_equal, with the tolerance held in a named const double. That tolerance is a choice with a working range. Near 1e16 two adjacent doubles differ by more than 1e-9, and near 1e-12 values a factor of three apart differ by less than it. Relative comparison is the general answer; a tolerance chosen to suit the magnitudes you handle is this course's answer.
  • Rounding error accumulates: ten additions of 0.1 give 0.99999999999999989, which prints as 1.00 under %.2f and still fails == 1.0. Never write a loop condition as != against a floating point total.
  • Use double by default. A float is 4 bytes with roughly 7 significant digits, for memory or bandwidth pressure. printf promotes a float to double so %f serves both, but scanf needs %f for a float * and %lf for a double *.
  • Floating point division by zero is defined: 1.0 / 0.0 is inf and 0.0 / 0.0 is nan, unlike integer 1 / 0, which is undefined behaviour. A NaN spreads through arithmetic and does not equal itself, so test it with isnan rather than with ==.