Explain why catastrophic cancellation happens

After my own research, the following picture emerges as the most frequently used example of catastrophic cancellation (It is indeed used in my class).

Could anyone explain why the plot takes that shape? (i.e., the widely fluctuating jagged line to the sides, and the vertical drop near the middle)

Minor point: why does the demonstrating function divides by $x^2$? Isn't $1 - \cos(x)$ enough to cause catastrophic cancellation?

Solution 1:

When $\cos(x)$ is computed for a reasonably small $x$, a number like 0.99999999999994598649857etc is obtained; the IEEE double precision standard defines a way to round this number to roughly 16 decimals. Thus out comes something like 0.99999999999995. I didn't count the digits, but you get the idea (and anyway, rounding is done in base 2). This does not seem catastrophic. But now subtract 1 from both numbers, and you see that a huge relative error has been committed to that result. The division by $x^2$ roughly gives you this relative error, because when $x$ is small, then $\cos(x) \approx 1 - \frac12 x^2$ (according to Taylor series). So when $x$ is very small, $\cos(x)$ is rounded to 1 (no digits after the comma), which explains the plateau in your plot.

Follow up on OP's request: To formalize this, let $fl$, for "floating point arithmetics", be the function that rounds a real number to the nearest machine number ($=$ a number representable as double precision IEEE number, of which there are only finitely many). Although maybe not strictly accurately speaking (there may be optimization; there may be hidden extended precision computation; and the rounding is not always to the nearest machine number), the machine computes not $(1 - \cos(x)) / x^2$ but $$ fl( fl( 1 - {\color{red}{fl}}(\cos(x)) ) / fl(x^2) ) , $$ supposing we start with $x$ such that $x = fl(x)$. The issue originates with the application of the $fl$ that is marked in red, because $$ {\color{red}{fl}}(\cos(x)) = \cos(x) (1 + \delta) $$ where $|\delta| \leq \text{eps} \approx 10^{-16}$, and $\text{eps}$ is the smallest positive number such that $$ fl(1 + \text{eps}) \neq 1 . $$ Therefore, $$ 1 - {\color{red}{fl}}(\cos(x)) = 1 - \cos(x) + \delta \cos(x) \approx \tfrac12 x^2 + \delta $$ when $x$ is not too small (otherwise ${\color{red}{fl}}(\cos(x)) = 1$). When $x$ is of order $\sqrt{\text{eps}} \approx 10^{-8}$ (as in your picture), this result $\tfrac12 x^2$ is contaminated by an absolute error of $\delta$, or a relative error of $2 \delta / x^2$. So the final result $(1 - \cos(x)) / x^2$ is contaminated by an absolute error of $\delta / x^2 \approx 1$ when $x \approx \sqrt{\text{eps}}$, which, again, is what you see in the picture.

ps. There have been infamous consequences of finite precision arithmetics.