Stable and Unstable Algorithms - Introduction to Scientific Computing

In the previous section, we saw that the condition number measures how sensitive a problem is to input perturbations. But sensitivity of the problem is only half the story — the algorithm we choose matters just as much.

Definitions¶

Example: Two Algorithms for the Same Function¶

Example 1 (Stable vs. Unstable Evaluation of

\sqrt{x+1} - \sqrt{x}

)

Consider $f(x) = \sqrt{x+1} - \sqrt{x}$ . The condition number satisfies $\kappa \leq 1/2$ — this is a well-conditioned problem.

Algorithm 1 (unstable): Evaluate directly using 4-digit decimal arithmetic:

f(1984) = \sqrt{1985} - \sqrt{1984} = 0.4455 \times 10^2 - 0.4454 \times 10^2 = 0.0001 \times 10^2 = 0.1000 \times 10^{-1}

(1)

The true value is $0.1122 \times 10^{-1}$ . We lost two significant digits — 10% relative error from a well-conditioned problem! The subtraction of nearly equal numbers is ill-conditioned, with $\kappa = (|a|+|b|)/|a-b| \approx 0.8909 \times 10^2 / 0.0001 \times 10^2 \approx 10^4$ , causing catastrophic cancellation.

Algorithm 2 (stable): Rationalize the expression:

\sqrt{x+1} - \sqrt{x} = \frac{(\sqrt{x+1} - \sqrt{x})(\sqrt{x+1} + \sqrt{x})}{\sqrt{x+1} + \sqrt{x}} = \frac{1}{\sqrt{x+1} + \sqrt{x}}

(2)

Now:

f(1984) = \frac{1}{\sqrt{1985} + \sqrt{1984}} = \frac{1}{0.4455 \times 10^2 + 0.4454 \times 10^2} = \frac{1}{0.8909 \times 10^2} = 0.1122 \times 10^{-1} \quad \checkmark

(3)

Same problem, different algorithm, full accuracy. Algorithm 1 is unstable because it introduces an ill-conditioned subtraction step; Algorithm 2 avoids the subtraction and is stable.

Example: The Finite Difference Revisited¶

The same pattern — well-conditioned problem, unstable algorithm — appears in a computation we already know well.

Example 2 (Computing

f'(x)

— Well-Conditioned Problem, Unstable Algorithm)

Consider computing $f'(x)$ for $f(x) = \sin(x)$ at $x = 1$ . The problem’s condition number is:

\kappa_{\text{problem}} = \left|\frac{x\,f''(x)}{f'(x)}\right| = \left|\frac{-\sin(1)}{\cos(1)}\right| = |\tan(1)| \approx 1.56

(4)

The problem is well-conditioned. But the forward difference algorithm computes $f'(x)$ via the subtraction $f(x+h) - f(x)$ , whose condition number is:

\kappa_{\text{subtraction}} \approx \frac{2|f(x)|}{h|f'(x)|} \sim \frac{1}{h} \to \infty \quad \text{as } h \to 0

(5)

The algorithm introduces an ill-conditioned subtraction step that the problem itself does not require. A different method (e.g., automatic differentiation) avoids this subtraction entirely and is stable.

Recognizing this pattern — is the problem sensitive, or is the algorithm choosing a sensitive path? — is one of the central skills in numerical analysis.

Floating-Point Operations to Watch For¶

Not every source of error is an ill-conditioning issue. The three operations below can each cause trouble, but for different reasons:

Subtracting nearly equal numbers ( $x - y$ when $x \approx y$ ): Subtraction of nearly equal numbers is genuinely ill-conditioned — the condition number of subtraction is $\kappa = (|a|+|b|)/|a-b| \to \infty$ as $a \to b$ . The key insight is that an algorithm may choose to subtract nearly equal numbers as an intermediate step, even when the overall problem doesn’t require it. If the subtraction can be avoided by algebraic reformulation, the algorithm becomes stable.
Dividing by a small number ( $x / y$ when $|y| \ll 1$ ): This is not ill-conditioned — the condition number of $f(y) = x/y$ with respect to $y$ is $\kappa = 1$ . The issue is that small absolute errors in $y$ become large absolute errors in the result. If $y$ itself was computed with some absolute error $\delta$ , then $x/(y+\delta) - x/y \approx -x\delta/y^2$ , which is large when $y$ is small. Remedy: rewrite using logarithms or rescale the computation.
Adding numbers of vastly different magnitude ( $x + y$ when $|x| \gg |y|$ ): This is also not an ill-conditioning issue. Subtraction is not involved, so there is no cancellation. The problem is that floating-point numbers have finite precision: if $y$ is smaller than the spacing between consecutive floats near $x$ , then $\text{fl}(x + y) = x$ and $y$ is simply lost. Remedy: sort and sum smallest first, or use compensated summation (Kahan summation).

Looking Ahead¶

We will formalize these ideas using forward and backward error when we study linear systems. There, we will see that backward stability — the precise version of “solves a nearby problem exactly” — is the gold standard for numerical algorithms. For now, the key takeaway is: distinguish between the problem (conditioning) and the algorithm (stability).