Skip to content

English · Español

03 — Orden de suma, cancelación y Kahan

🇪🇸 La segunda clase de fallos numéricos: cuando sumas muchos números o restas dos casi iguales, los bits "bajos" desaparecen. La cancelación catastrófica y la suma de Kahan son dos caras del mismo problema: la mantisa es finita. Ejemplo concreto: sumar las 600 probabilidades por forma verbal del vocabulario (§A13) y verificar que la suma da 1.0 con suficiente precisión.

El planteamiento

La suma ingenua en coma flotante tiene error proporcional a la longitud de la suma. Si sumas N números cada uno de magnitud μ, el resultado tiene error absoluto en el peor caso O(N · ε · μ) donde ε es el épsilon de máquina. Para N = 10⁶ sumas fp32, eso es 10⁶ × 1.2 × 10⁻⁷ × μ ≈ 0.12 · μ — un 12% de error relativo, sobre una suma que debería ser exacta a 7 dígitos.

En el universo de §A13, el patrón de suma más común es "normalizar un vector de probabilidades" o "calcular la función de partición". Un vector de longitud 600 de probabilidades, cada una de magnitud ~1/600, sumado ingenuamente, acumula 600 × ε × (1/600) ≈ ε. Tolerable pero no despreciable. Un acumulador de longitud 10⁶ sobre las pérdidas de muchos pasos de entrenamiento no es tolerable.

Esta página hace dos cosas:

  1. Cuantifica los modos de fallo (cancelación catastrófica, crecimiento del error en sumas).
  2. Te da los dos arreglos estándar (Kahan, Neumaier).

Modo de fallo 1 — cancelación catastrófica

Resta dos números en coma flotante casi iguales. El resultado es pequeño, pero los dígitos líderes se cancelan, y lo que queda lo determinan los bits de orden bajo de la mantisa — que son precisamente los bits con el mayor error relativo de redondeo.

Ejemplo concreto. Digamos que dos probabilidades para formas verbales adyacentes en nuestra distribución están muy próximas:

p_a = 0.16739871        # P(token = "worked", in fp32)
p_b = 0.16739864        # P(token = "walked", in fp32)

Ambos se almacenan con ~7 dígitos decimales de precisión; los valores representables reales difieren de los decimales declarados en ~p × ε.

Calcula p_a - p_b simbólicamente: 7e-8. Calcúlalo en fp32: los 6 dígitos líderes se cancelan, y lo que queda es la diferencia del último ~1 bit de la mantisa de cada operando. El resultado es correcto dentro de ~1 ULP (unit in the last place), pero ese ULP es ahora la magnitud entera de la respuesta. El error relativo del resultado es del 100%.

Esto es cancelación catastrófica. No da NaN, no desborda, no imprime un aviso. Simplemente te da un número que parece plausible y es ruido total.

Patrones que esconden cancelación catastrófica

La formulación estilo C es (big_number) - (almost_the_same_big_number). Dónde se esconde esto en aprendizaje automático (machine learning):

  • log(1 + x) para x diminuto. 1 + x ≈ 1; la suma descarta los bits significativos de x; entonces log(1) es exactamente 0. Usa log1p(x), que se calcula mediante una serie de Taylor y preserva los bits.
  • exp(x) - 1 para x diminuto. Misma forma. Usa expm1(x).
  • sqrt(1 + x²) - 1 (distancia pitagórica). Numéricamente inestable para x pequeño. Usa el reordenamiento algebraico x² / (sqrt(1+x²) + 1).
  • Depuración sustractiva del softmax. softmax(x)[i] - softmax(y)[i] para dos vectores de logits casi idénticos parecerá ruidoso sobre probabilidades pequeñas. Calcula log-probabilidades y resta esas en su lugar.
  • Varianza vía E[X²] - E[X]² (la fórmula "ingenua"). Si E[X] es grande, E[X²] y E[X]² son ambos grandes y casi iguales; restar cancela. Usa el algoritmo online de Welford.

Una heurística útil: si tu fórmula tiene una resta de dos términos que esperas que sean próximos, pregúntate si puedes transformarla algebraicamente en una suma.

Modo de fallo 2 — deriva en sumas largas

La suma ingenua acumula un error de redondeo por cada suma. Tras N sumas:

\[ \mathrm{error}_\mathrm{naive} \leq N \cdot \varepsilon \cdot \max_i |s_i| \]

Donde s_i es la suma parcial corriente. En el peor caso (todos los sumandos positivos), |s_i| crece linealmente con i, dando error total O(N² · ε). En el caso medio (signos aleatorios), el error es O(√N · ε) (paseo aleatorio). Para sumas de todos positivos de términos de magnitud similar — el caso de uso típico de aprendizaje automático — la cota relevante es O(N · ε).

Para N = 10⁶ y ε = 1.19 × 10⁻⁷ (fp32), la cota es ~0.12. Doce por ciento. El caso fp64 (ε ≈ 2.2 × 10⁻¹⁶) acota en ~2.2 × 10⁻¹⁰ — invisible — por lo cual "usa fp64 para acumuladores" es un mantra común.

La solución — suma de Kahan

La idea: rastrea el error de redondeo de cada suma y reinyéctalo en el siguiente paso.

sum   = 0          # running total
c     = 0          # running compensation (the "lost" bits from prior rounds)
for x in seq:
    y    = x - c       # correct the input by the accumulated error
    t    = sum + y     # this addition discards low-order bits of y
    c    = (t - sum) - y   # what got discarded
    sum  = t

La magia es la línea c = (t - sum) - y. En aritmética exacta, esto es cero. En aritmética fp, t = sum + y redondeado, así que t - sum (calculado en fp) no es exactamente y — es y menos lo que se descartó. Así que (t - sum) - y es el negativo de lo que se descartó. En la siguiente iteración, eso vuelve a sumarse vía x - c.

El error total se vuelve O(ε) independiente de N. Específicamente:

\[ \mathrm{error}_\mathrm{Kahan} \leq (2 \varepsilon + O(N \varepsilon^2)) \cdot \sum |s_i| \]

Para N = 10⁶ fp32, la cota es ~2.4 × 10⁻⁷ de la suma — cinco órdenes de magnitud mejor que ingenuo.

Variante de Neumaier

El paso de compensación de Kahan asume |sum| ≥ |y|. Si |y| > |sum| (el nuevo término es enorme), Kahan pierde precisión. Neumaier (1974) añadió una guarda:

sum   = 0
c     = 0
for x in seq:
    t = sum + x
    if abs(sum) >= abs(x):
        c += (sum - t) + x       # lost low bits of x
    else:
        c += (x - t) + sum       # lost low bits of sum
    sum = t
sum += c                         # add accumulated compensation at the end

Para sumas de aprendizaje automático con términos de magnitud similar, Kahan y Neumaier se comportan idénticamente. Para magnitud mixta (suma de términos 1e10 y 1e-10), Neumaier gana.

Un ejemplo trabajado — sumar 600 probabilidades por forma verbal

Vocabulario de §A13: 600 formas. Supón que la probabilidad predicha de tu modelo para cada forma, en fp32, suma (matemáticamente) exactamente 1.0. Cada p_i ≈ 1/600 ≈ 0.001667.

Suma ingenua:

total = 0.0
for p in probabilities:    # length 600
    total += p

Error en el peor caso: 600 × ε × 1/600 = ε ≈ 1.2 × 10⁻⁷. Así que total aterriza en [1 - 1.2e-7, 1 + 1.2e-7]. Bien para casi cualquier uso.

Pero escala esto. Supón que el bucle de entrenamiento de Borja acumula pérdidas por batch sobre 100,000 batches × 64 ejemplos cada uno = 6.4 × 10⁶ sumandos:

running_loss = 0.0
for batch in train_loader:
    running_loss += loss_of_batch

Ahora N = 100,000 (o 6.4 × 10⁶ si sumas por ejemplo). El error crece a N × ε × μ, que para μ ≈ 2.0 (CE típica) y N = 10⁵ da ~0.02. Dos por ciento de error en el promedio corriente reportado. El "loss promedio" que Borja ve en su log de entrenamiento está mal por un 2% solo por la suma ingenua.

La solución:

total = 0.0
c     = 0.0
for p in probabilities:
    y     = p - c
    t     = total + y
    c     = (t - total) - y
    total = t

O, equivalentemente en NumPy: usa np.float64 para el acumulador y np.float32 para los datos. Subir el casting te da precisión equivalente a Kahan gratis, al coste de 8 bytes por acumulador (despreciable). El lab 02 mide ambos enfoques.

Cuándo usar realmente Kahan vs acumuladores fp64

Situación Mejor solución
Bucle interno, crítico en rendimiento, sin kernel consciente de Kahan Acumulador fp64
Solo tienes hardware fp32 y el bucle interno es grande Kahan
Estás calculando varianza / media Algoritmo online de Welford
Sumando gradientes a través de dispositivos Reducción en árbol (incorporada en NCCL, FSDP)
Calculando log-likelihood sobre N tokens logsumexp es una suma estable; no se necesita truco extra

Para los laboratorios de la Fase 2, implementarás Kahan en NumPy y compararás con el ingenuo con acumulador fp64. La conclusión: para aprendizaje automático, "usa fp64 para acumuladores" es la solución correcta más barata en casi todos los casos. Reserva Kahan para cuando fp64 no esté disponible (algunos caminos de inferencia embebida).

Determinismo

Una consecuencia sutil de la no asociatividad: el orden de la suma cambia el resultado, incluso con el mismo conjunto de entradas. Para reproducibilidad:

  • Fija el orden. np.sum(arr) usa suma por pares (ya mucho mejor que ingenuo); su orden depende de la forma, eje y stride.
  • Evita reducciones paralelas en sitios donde la reproducibilidad bit-exacta importe. NumPy es mono-hilo por defecto; PyTorch no (Fase 25).
  • Tests de CI que comprueben loss == expected_loss hasta 1e-6 serán flaky si el orden de reducción varía entre ejecuciones (p. ej., dependiendo de OMP_NUM_THREADS). Usa tolerancia relativa o acumuladores fp64.

scripts/seed_everything (Fase 0) pone los thread counts a 1 por esto durante los tests deterministas.

Problemas de práctica

Soluciones en solutions/03-summation-ref.md (al abrir la fase).

  1. Suma el vector fp32 [1.0] + [1e-9] * 10⁶ ingenuamente. ¿Qué esperas? ¿Qué da Kahan? ¿Qué da el acumulador fp64? ¿Por qué?
  2. Calcula var([1e10 + 1, 1e10 + 2, 1e10 + 3]) de dos formas: con la fórmula ingenua E[X²] - E[X]² y con el algoritmo online de Welford. Compara con la respuesta exacta (var = 2/3).
  3. Calcula log(1 + 1e-9) y log1p(1e-9). Predice el error relativo de cada uno; luego verifica.
  4. Muestra que la suma de Kahan de una lista de longitud N de números de igual magnitud tiene error total O(ε), no O(N·ε). (Pista: acota el término de compensación.)
  5. Dadas las 600 probabilidades del vocabulario de §A13, cada una ~1/600 + δ_i con Σ δ_i = 0 (la distribución está normalizada), ¿cuál es el error esperado de la suma ingenua? ¿Cuándo importa?

Recapitulación en un párrafo

La suma fp ingenua acumula un ε de error por cada suma; sobre N términos crece a O(N·ε), que es ~12% para N = 10⁶ en fp32 — inaceptable para acumuladores de aprendizaje automático. La cancelación catastrófica en restas de valores casi iguales destruye el resultado en silencio. Las dos soluciones son (1) suma de Kahan (rastrea y reinyecta el error de redondeo en cada paso, dando O(ε) total) y (2) log1p / expm1 / algoritmo de Welford y reformulaciones similares para expresiones propensas a cancelación. Para aprendizaje automático, la solución correcta más barata es "usa fp64 para acumuladores"; recurre a Kahan cuando fp64 no esté disponible.

Lo que esta página NO cubre

  • Reducciones paralelas en árbol (postergadas a la Fase 35 distribuida).
  • Suma por pares (la opción predeterminada de NumPy; un punto medio entre ingenua y Kahan).
  • Productos escalares compensados. Fuera de alcance.
  • Aritmética exacta (fractions, Decimal). Usadas como oráculos en tests solamente.

Siguiente: theory/04-precision-zoo.md.