Lo que todo científico —o programador— debería saber sobre la aritmética de coma flotante

Los cálculos que realizamos en nuestros ordenadores no son exactamente los mismos cálculos que realizamos con nuestra cabeza. En nuestra cabeza, podemos dividir 1 entre 3, y sabemos que es exactamente un tercio (1/3), o bien aproximadamente 0,333333…, indicando los puntos suspensivos una repetición del periodo.

Si estamos haciendo un cálculo con números en coma flotante, sabemos que 1 – 0,9 – 0,1 es igual a cero (después de todo, 0,9 más 0,1 debe ser 1), y sabemos que eso es equivalente a la operación 1 – 9/10 – 1/10.

Pero si usamos un ordenador para calcularlo, por ejemplo, usando Python…

mac:~ user$ python
Python 2.3.5 (#1, Mar 20 2005, 20:38:20)
[GCC 3.3 20030304 (Apple Computer, Inc. build 1809)] on darwin
Type “help”, “copyright”, “credits” or “license” for more information.
>>> print 1 – 0.9 – 0.1
-2.77555756156e-17

¿Por qué el resultado no es 0?

La respuesta, básicamente, es que los números que son periódicos puros en una base no tienen por qué serlo en otra… y en este caso, por ejemplo, 1/10, en binario, es un número periódico, por lo que al truncarlo siempre tendremos un error que dependerá de la precisión numérica que estemos utilizando.

En un interesantísimo artículo, titulado What Every Computer Scientist Should Know About Floating-Point Arithmetic, y que ahora forma parte de la Numerical Computation Guide de Sun, David Goldberg se adentra en estos aspectos, y nos explica cómo establecer comparaciones entre números en coma flotante, y cómo tener en cuenta la propagación de errores de redondeo que hace que se produzcan resultados como el anteriormente expuesto.

Por cierto, el artículo lo encontré en la bitácora NSBlog de Mike Ash, una bitácora muy recomendable que encontré gracias a la no menos recomendable Lista Enlazada de John Gruber…

[Actualización 1] El artículo de la guía de cálculo numérico de Sun habla de cómo se describen los números en coma flotante en los diferentes lenguajes de programación, y habla de conceptos como exponentes y mantisas, pero sin especificar sus tamaños. Existe un estándar, el IEEE-754, que describe justamente los tamaños y forma de uso de los diferentes formatos de coma flotante (simple precisión, doble precisión, cuádruple precisión…), y podéis encontrarlo en los enlaces.

Si se calcula el valor absoluto del logaritmo en base 2 del error relativo cometido al calcular 1 – 0,9 – 0,1, obtendríamos el tamaño de la mantisa de un número en coma flotante de doble precisión (52):

>>> from math import *; print abs(log(2.77555756156e-17/0.1)/log(2))
51.6780719051

Recordemos que el logaritmo en base n de un número x se puede calcular siempre como el logaritmo en cualquier base de x, entre el logaritmo en la misma base de n:

logn(x) = log(x)/log(n).

Enlaces

• Sun’s Numerical Computation Guide: What Every Computer Scientist Should Know About Floating-Point Arithmetic
• W. Kahan: Lecture Notes on the Status of IEEE Standard 754 for Binary Floating-Point ArithmeticVia NSBlog < Daring Fireball Linked List