Where it falls apart: catastrophic cancellation

Run the column further and the answers get worse, then silly, then exactly zero. The machine, not the mathematics.

The loop in the last lesson stopped at $h = 10^{-7}$, and it stopped there on purpose. The algebra says the quotient is $6 + h$, so smaller steps should give better answers forever. Take the same loop out to $h = 10^{-18}$ and something happens that the algebra gives no warning about at all.

Three separate things are visible there. Around $10^{-8}$ the answers stop getting better. Below $10^{-12}$ they get actively worse, and by $10^{-14}$ the first decimal place is wrong. At $10^{-16}$ the answer is exactly zero, reported with no complaint of any kind. If you had trusted the loop, you would have concluded that the ball was standing still.

A double holds about sixteen digits and no more

Nearly every language stores a number with a decimal point in it as a double: 64 bits, of which 52 are given over to the digits of the value and the rest carry its sign and its size. There are 53 digits in play rather than 52, because the leading digit of a binary number is always 1 and so never has to be stored. Fifty-three binary digits works out at roughly 15 to 17 significant decimal digits. Not 15 to 17 decimal places, note; 15 to 17 significant digits, counted from the first one that is not a zero.

Anything beyond that is not stored inaccurately. It is not stored. The value is rounded to the nearest number the format can represent and the rest is gone, permanently, before your next line of code runs.

The second line is the famous one, and it is worth understanding rather than memorising: 0.1 and 0.2 are not exactly representable in binary any more than a third is exactly representable in decimal, so their stored versions are slightly off and their sum is slightly off in a way that does not happen to match the stored version of 0.3.

The third line is the one that matters here. Nine plus a ten-thousandth of a trillionth needs seventeen significant digits to write down. The format has sixteen. So the sum is rounded, and it rounds back to exactly 9.

The subtraction is where the digits die

Now watch what that does inside the difference quotient. With $a = 3$ and $h = 10^{-8}$, the rounding happens at a + h, before the squaring ever runs. The sum you asked for is 3.00000001, and the nearest double to it is 3.000000009999999939, so a little has gone already. Squaring that stored value gives 9.000000059999999635, where the true $(3 + 10^{-8})^2$ is 9.0000000600000001.

a + h       wanted   3.000000010000000000
            stored   3.000000009999999939   rounded already
f(a + h)    stored   9.000000059999999635   the stored a + h, squared
f(a)        stored   9.000000000000000000   exact, nothing to round
                     --------------------
difference           0.000000059999999635   eight good digits, then noise

Sixteen digits went into the subtraction and eight came out. The leading nine and the seven zeros after it were identical in both numbers, so they cancelled to nothing, and the digits that would have replaced them at the bottom end were thrown away when the values were stored. Dividing by h scales the survivors back up: 0.000000059999999635 over $10^{-8}$ is 5.999999963535, which is the row in the table. The exact answer is 6.00000001, so the scaling has enlarged the noise along with the answer.

This is called catastrophic cancellation: subtracting two numbers that are nearly equal, so that the significant digits cancel and what is left is short, and made mostly of whatever rounding error was already there. The word catastrophic is doing real work. The subtraction itself is computed perfectly; the loss happened earlier, when the operands were rounded, and the subtraction merely exposes it.

At $k = 16$ the line reads a + h = 3.0. The step was absorbed entirely: 3 + 1e-16 rounds to 3, so the top of the fraction is a number minus itself, which is exactly zero, and zero divided by anything is zero. The program is not confused and does not warn you. It has computed exactly what you asked, and the answer is worthless.

Two errors, pulling opposite ways

There are two sources of error in that column and they behave in opposite directions, which is why the total has a floor rather than heading to nothing.

The rule of thumb is worth remembering, because it explains the shape of the column exactly. The smallest relative step a double can resolve is about $2.2 \times 10^{-16}$. The best h is roughly the square root of that, about $1.5 \times 10^{-8}$, and at that step you get about half the digits: eight of them. That is the row where the answers stopped improving.

So the experiment has a ceiling. A program can take h down to about $10^{-8}$ and no further, which means no program can ever finish the journey to $h = 0$. If the value 6 is going to be established, something other than a loop has to establish it, and that something is the subject of the next lesson.