Numerical Stability and Conditioning

Why you should care (before any formulas)

You can write down a correct mathematical expression and still get nonsense on a real computer.

• •A softmax that should produce probabilities returns NaN because $e^{1000}$ overflowed.
• •A variance that should be nonnegative becomes slightly negative, and then a square root yields NaN.
• •The same sum of numbers gives different answers depending on the order you add them.

These aren’t “bugs” in deep learning libraries—they’re consequences of finite-precision arithmetic.

Numerical analysis separates two ideas:

1. 1)Conditioning (problem sensitivity): if we slightly perturb the input, how much can the true answer change?
2. 2)Stability (algorithm robustness): given rounding errors and finite precision, does our computed answer stay close to the true answer of a nearby input?

A well-conditioned problem can be solved poorly by an unstable algorithm. A well-designed stable algorithm cannot fully rescue an ill-conditioned problem—but it can avoid making things worse.

Floating-point reality: rounding + limited dynamic range

Computers represent real numbers in floating-point format (float16/float32/float64). Two consequences matter constantly:

• •Rounding error: most real numbers cannot be represented exactly, so operations are rounded.
• •Limited dynamic range: numbers that are too large overflow to ∞; too small underflow toward 0 (sometimes with subnormals).

A standard model of floating-point arithmetic is:

where:

• •$\circ$ is $+$, $-$, $\times$, or $/$,
• •$u$ is the unit roundoff, often roughly half of machine epsilon.

We’ll use machine epsilon (often denoted $\varepsilon_{\text{machine}}$) as the “scale” of relative rounding error.

Typical values (approximate):

Two important takeaways:

• •float16 has much worse precision and a much smaller max value than float32.
• •bfloat16 has low precision like float16 but a large exponent range like float32.

Conditioning: “how sensitive is the math problem?”

Consider computing $f(x)$. If a tiny perturbation $\Delta x$ can cause a large change in $f(x)$, the problem is ill-conditioned near $x$.

A classic measure is the relative condition number:

Interpretation:

• •If $\kappa_f(x)$ is large, then a small relative error in $x$ can produce a large relative error in $f(x)$.

Example intuition:

• •$f(x) = x$ has $\kappa \approx 1$ (well-conditioned).
• •$f(x) = 1/x$ near $x=0$ is ill-conditioned (tiny changes in $x$ blow up the output).

In deep learning, conditioning shows up as:

• •nearly singular matrices (e.g., Hessians, Jacobians)
• •poorly scaled features/activations
• •exploding/vanishing gradients

Stability: “does the algorithm amplify rounding errors?”

An algorithm is (backward) stable if the computed result equals the exact result for a slightly perturbed input:

This is the gold standard: rounding happens anyway; stability means your algorithm behaves like it just saw a tiny bit of noise.

Key mental model:

• •Conditioning is about the landscape you’re standing on.
• •Stability is about whether you’re walking carefully or running with untied shoelaces.

Supporting static diagram (conceptual): float spacing vs magnitude

A core stability intuition is that floats are denser near 0 and sparser at large magnitudes.

value magnitude:   0 ---- 1 ---- 10 ---- 1e3 ---- 1e6 ---- 1e9
float spacing:     tiny  small  bigger   larger   huge     enormous

At large magnitudes, the next representable number can be far away, so adding a small number can do nothing (it rounds away). This drives reduction-order issues and “small updates disappear” problems.

This section builds a concrete catalog of failure modes you will repeatedly see in ML code. The goal is to recognize them quickly, then reach for the right stable pattern.

1) Overflow and underflow

What it is

• •Overflow: a value exceeds the largest representable finite float → becomes ∞.
• •Underflow: a value is too small in magnitude → becomes 0 (or subnormal, then 0).

Why ML triggers it

Exponentials, products over many terms, and repeated multiplications are common.

• •Softmax uses $e^{z_i}$.
• •Likelihoods multiply probabilities.
• •Deep nets multiply Jacobians across layers.

Classic example: softmax overflow

Naive softmax:

If some $z_i$ is large (say 1000 in float32), $e^{z_i}$ overflows to ∞, and you get ∞/∞ → NaN.

2) Catastrophic cancellation

What it is

Subtracting two nearly equal numbers destroys significant digits.

If $a \approx b$, then $a - b$ is small, but the floating-point representation of $a$ and $b$ each has rounding error. The subtraction exposes that error.

Why ML triggers it

• •variance computations
• •normalization
• •loss differences
• •some distance computations

Example: variance via “E[x²] − (E[x])²”

The identity is correct:

But numerically, if $x$ has large mean and small variance, both terms are large and close, so subtraction cancels digits. You may get a small negative result.

3) Reduction-order (summation) error

What it is

Floating-point addition is not associative:

because each intermediate sum is rounded.

Why ML triggers it

Training is full of reductions:

• •summing losses over a batch
• •summing gradients
• •dot products and matrix multiplies

On GPU/TPU, parallel reductions change order depending on hardware and scheduling. That can change results—usually slightly, sometimes significantly.

A helpful intuition

If you add a tiny number to a huge number, the tiny number may not change the sum at all (it gets rounded away). So order matters: adding small numbers together first can preserve more information.

4) Division by tiny numbers / unstable normalization

Normalization often has the form:

If $\sigma$ is near 0, this blows up. Many algorithms add an $\varepsilon$ term:

The $\varepsilon$ here is a stability constant, not machine epsilon (though it’s motivated by it). In float16, $\varepsilon$ often must be larger.

5) Gradients that vanish due to underflow or saturation

Even if your forward pass is stable, gradients can underflow or become 0 in low precision.

Examples:

• •sigmoid/tanh saturation: for large |x|, derivatives are near 0.
• •multiplying many small factors across layers.

This interacts with floating-point: values smaller than ~1e−8 may become subnormals or 0 in float16.

6) Infs/NaNs that propagate through a computational graph

Given prerequisites about computational graphs, here’s the key operational point:

• •One NaN in forward → many NaNs in backward.
• •One overflow in an activation (e.g., exp) → NaNs in loss → NaN gradients.

This is why stable local formulas (softmax, log-softmax, normalization) matter so much: they sit early in the graph and affect everything.

Interactive canvas specification (explicitly maps to 3 failure modes)

To address these issues systematically, your interactive canvas should have three mini-labs, each with:

• •a toggle for precision: float16 / float32 / float64
• •side-by-side traces: naive vs stable
• •a plot of intermediate values (including detection of inf/nan)

Mini-lab A: Softmax overflow

Controls:

• •slider: max logit magnitude (e.g., from 0 to 2000)
• •toggle: stable shift (subtract max) on/off
• •precision selector

Display:

• •naive exp(logits), sum exp, softmax output
• •stable exp(logits − max), sum exp, softmax output
• •highlight when overflow happens

Mini-lab B: Cancellation in variance

Controls:

• •slider: mean magnitude μ (e.g., 0 to 1e6)
• •slider: true std σ (e.g., 1e−2 to 10)
• •toggle: one-pass (E[x²]−E[x]²) vs Welford
• •precision selector

Display:

• •computed variance (naive vs stable)
• •difference from ground truth
• •highlight negative variance events

Mini-lab C: Reduction-order error

Controls:

• •choose data distribution: (many small + few large), alternating signs, random
• •toggle: left-to-right sum vs pairwise sum vs Kahan
• •precision selector

Display:

• •final sum
• •absolute/relative error vs high-precision reference
• •step-by-step accumulation trace to show where small terms get lost

This canvas ties directly to: overflow in softmax, cancellation in variance, and reduction-order error—the three most common “why did my model blow up?” numerics in practice.

Now we turn the failure modes into stable, reusable techniques.

1) Log-sum-exp (LSE): the workhorse for stable exponentials

Motivation

You want to compute:

Directly computing $e^{z_i}$ is risky if any $z_i$ is large.

Trick: factor out the maximum

Let $m = \max_i z_i$. Then:

Taking log:

Key benefits:

• •all exponents $z_i - m \le 0$, so $e^{z_i - m} \in (0, 1]$
• •prevents overflow
• •reduces risk of underflow dominating the sum (though very negative values may still underflow safely to 0)

Connection to softmax and log-softmax

Stable softmax:

Stable log-softmax:

In practice you compute log-softmax directly (stable) and then NLL loss, rather than softmax then log.

2) Stable variance: Welford’s algorithm

Motivation

Variance is everywhere: batch norm, layer norm, feature scaling, monitoring.

The naive formula is prone to cancellation:

Welford’s one-pass stable update

Given samples $x_1, x_2, \dots, x_n$, maintain:

• •running mean $\mu_k$
• •running sum of squares of deviations $M2_k$

Initialize:

• •$\mu_1 = x_1$
• •$M2_1 = 0$

For $k \ge 2$:

Then:

Why it helps

It avoids subtracting two nearly equal large numbers. It tracks deviations relative to an evolving mean, keeping intermediate quantities closer to the scale of the variance.

3) Stable summation: pairwise and Kahan

Motivation

Summing many numbers is ubiquitous: losses, gradient accumulation, dot products.

Pairwise summation

Instead of summing linearly, sum in a tree (balanced) structure. This tends to reduce error growth from $O(nu)$ to closer to $O(\log n \cdot u)$ in many practical settings.

Kahan compensated summation

Maintain a compensation term $c$ for lost low-order bits.

Algorithm (for summing values $x_i$):

• •$s = 0$, $c = 0$
• •For each $x$:
• •$y = x - c$
• •$t = s + y$
• •$c = (t - s) - y$
• •$s = t$

Interpretation: $c$ estimates what was lost to rounding.

Practical note for ML

• •Many BLAS/GEMM kernels already use strategies to reduce error, but you still see reduction issues in custom kernels, metrics, and loss aggregation.
• •Mixed precision often accumulates in float32 even if operands are float16 (a major stability improvement).

4) Normalization and scaling as stability tools

Normalization isn’t only about optimization speed; it’s also about preventing overflow/underflow and controlling conditioning.

Feature scaling

If one feature has magnitude ~1e6 and another ~1e−3, dot products and gradients become ill-scaled. Scaling features can reduce condition numbers and improve numerical behavior.

Activation scaling and initialization

Careful initialization (e.g., Xavier/He) aims to keep variance of activations/gradients roughly stable across depth. That’s partly a conditioning story (keeping Jacobian singular values reasonable) and partly a numerical range story (avoid overflow/underflow).

5) Conditioning vs stability: a quadrant you should internalize

A useful mental diagram:

Deep learning often has locally ill-conditioned regions (e.g., during early training or near sharp minima). Your job is to avoid compounding that with unstable computations.

6) Mixed precision stability patterns

Modern training often uses float16/bfloat16 for speed, with float32 master weights or accumulators.

Common stability patterns:

• •Loss scaling: multiply loss by a scale factor to prevent gradient underflow in float16; then unscale before the optimizer step.
• •FP32 accumulation: accumulate sums/products in float32 (especially reductions).
• •Keep sensitive ops in FP32: softmax/log-softmax, normalization stats, some reductions.

These are not “cheats”—they’re explicit numerics engineering.

This section connects the numerics toolbox to everyday deep learning components.

1) Softmax + cross-entropy: do not compute them separately (naively)

In classification, the loss is often:

Naively:

1) compute softmax

2) take log

This is numerically fragile because:

• •softmax can overflow
• •probabilities can underflow to 0, then log(0) → −∞

Stable approach:

• •compute log-softmax via log-sum-exp
• •index the correct class

This avoids both overflow and log(0) underflow.

2) BatchNorm/LayerNorm: stable mean/variance and safe eps

Normalization layers compute statistics over batches/features.

Key stability decisions:

• •Use stable variance computation (often a two-pass or Welford-like approach internally).
• •Choose $\varepsilon$ appropriate to dtype and scale.

Important nuance:

• •$\varepsilon$ is not a universal constant; it depends on typical activation magnitudes and precision.
• •In float16, using too-small $\varepsilon$ can fail to prevent division blow-ups.

3) Attention and Transformers: “scaled dot-product” is a numerics fix

Attention uses logits:

Then softmax over $j$.

The $1/\sqrt{d}$ scaling is partly about gradients/optimization, but also about numerics:

• •without scaling, dot products can grow with $d$, producing large logits
• •large logits make softmax overflow more likely and cause saturation (one-hot attention)

FlashAttention and similar kernels implement careful stabilization (blockwise softmax with max-subtraction) explicitly to keep values in safe ranges.

4) Log-likelihoods, HMMs, CRFs: stay in log space

Whenever you multiply many probabilities, underflow is inevitable:

Use logs:

And whenever you need to sum probabilities in log space, use log-sum-exp.

This pattern generalizes far beyond softmax.

5) Gradient clipping and normalization as overflow control

Exploding gradients can cause:

• •weight updates that overflow
• •activations that overflow next step

Gradient clipping (by norm) can be viewed as enforcing a bound on magnitude to keep values in range.

6) Debugging numerics systematically (a checklist)

When training becomes NaN/Inf:

1. 1)Locate first NaN/Inf in the forward pass (hooks / anomaly detection).
2. 2)Check common culprits:

• •softmax/log-softmax
• •normalization variance / division by small numbers
• •exp, log, sqrt
• •mixed precision overflow/underflow

3. 3)Try interventions:

• •replace softmax+log with log-softmax
• •increase normalization eps
• •enable FP32 accumulation
• •use loss scaling
• •clamp logits or activations (carefully)

A key mindset: you’re not only debugging code—you’re debugging a numeric process under rounding.

7) Connecting conditioning to optimization behavior

Conditioning appears as “how hard the optimizer must work.”

• •If the Hessian has a wide spread of eigenvalues, gradient descent zig-zags.
• •Normalization, adaptive optimizers, and preconditioning methods try to reduce effective condition number.

While this lesson focuses on floating-point issues, it’s important to see that:

• •better conditioning → often fewer extreme values → fewer numeric failures
• •stable numerics → fewer NaNs → more reliable optimization signals