Introducción al Cálculo Numérico. Errores

Introducción

El Cálculo Numérico, o Análisis Numérico, como también se le denomina, es la rama de las Matemáticas que estudia los métodos numéricos de resolución de problemas, es decir, los métodos que permiten obtener una solución aproximada (a veces exacta) del problema considerado tras realizar un número finito de operaciones lógicas y algebraicas elementales. Los métodos numéricos son imprescindibles para resolver aquellos problemas para los cuales no existen métodos directos de resolución, o bien para aquellos otros en los que existe un método directo pero no es viable, como ocurre con la regla de Cramer para la resolución de sistemas lineales, que teóricamente es válida para resolver cualquiera de ellos pero que en la práctica sólo es aplicable cuando la dimensión es muy pequeña. El coste operacional de un método es el número de operaciones que requiere su aplicación. Conocidos dicho número y la velocidad de cálculo del ordenador, podemos saber el tiempo que tardaría en ser resuelto. Por ejemplo, resolver un sistema lineal de orden poco mayor que $20$, con los ordenadores actuales y aplicando la regla de Cramer, supondría más de un millón de años, porque cuando el valor de un determinante de orden $n $ se calcula como $n$ determinantes de orden $n-1$ el número de operaciones que hay que hacer es, aproximadamente, del orden de MATH, valor que crece muy rápidamente cuando $n$ aumenta.

Los métodos numéricos han sido utilizados desde la antigüedad en la búsqueda de soluciones aproximadas para los problemas procedentes de todas las ciencias y de las ingenierías. De hecho, hasta hace unos dos siglos la mayor parte de la matemática que se hacía era de tipo numérico y muchas de las demostraciones eran de tipo constructivo. No obstante, el nombre de Análisis Numérico aparece por primera vez en 1948, en California, donde se fundó el Instituto de Análisis Numérico, y su desarrollo ha sido paralelo al del ordenador, habiéndose producido un gran crecimiento en la segunda mitad del siglo XX, que continúa en la actualidad. Así pues, el Análisis Numérico tal y como se concibe hoy es una materia reciente. El ordenador es una herramienta imprescindible para esta disciplina, ya que los métodos numéricos requieren un gran número de cálculos y manipulaciones de datos.

Noción de algoritmo

Como hemos indicado, un método numérico es un procedimiento para aproximar la solución de un problema. Un algoritmo es una descripción, sin ambigüedad, de un proceso -o sucesión finita de pasos- que se puede seguir para obtener la solución aproximada de dicho problema. Incluye desde la entrada de datos hasta la impresión de resultados. Puede establecerse más de un algoritmo para resolver un mismo problema. Un algoritmo se puede esquematizar como sigue: MATH

Los dos ejemplos extremadamente simples que siguen ilustran este esquema.

Representación de números en el ordenador

Los ordenadores operan en un sistema que no suele ser el decimal (base $10$), sino que lo hacen en binario (base $2$) o en otro cuya base sea una potencia de $2$ (por ejemplo $16$). Sin embargo la comunicación con el ordenador se hace siempre en base $10$. Los datos numéricos que damos al ordenador están escritos en el sistema decimal y los resultados numéricos que el ordenador devuelve también. Para pasar un número de la base decimal a la binaria debemos dividir entre $2$ y anotar el resto, volver a dividir el cociente anterior entre $2$ y anotar el resto, y así sucesivamente hasta que el cociente sea cero. Los restos escritos en orden inverso al que se han obtenido proporcionan la expresión del número en binario. Por ejemplo, para escribir el número $53$ del sistema decimal en el sistema binario se realizan las operaciones que se indican a continuación:

Resto
$53$ $=$ $26$$\times 2+1$ $\rightarrow $ $1$
$26$ $=$ $13$$\times 2+0$ $\rightarrow $ $0$
$13$ $=$ $6$$\times 2+1$ $\rightarrow $ $1$
$6$ $=$ $3$$\times 2+0$ $\rightarrow $ $0$
$3$ $=$ $1$$\times 2+1$ $\rightarrow $ $1$
$1$ $=$ $0$$\times 2+1\ $ $\rightarrow $ $1$

Los restos dan $110101$ como expresión binaria de $53$, lo que podemos notar como $53_{10}=110101_{2}$. El paso de la representación binario a decimal es muy fácil; por ejemplo, MATH

Números que tienen una representación finita en el sistema decimal pueden tener una representación infinita en binario y viceversa. Por ejemplo, MATH, donde la parte señalada se repite periódicamente. Al introducir el número $0.7$ en un ordenador que opere en binario tiene que ser aproximado redondeándolo por defecto o por exceso a un número cuya representación sea admisible por el ordenador (un número máquina, es decir, una cadena finita de ceros y unos).

La notación científica (o coma flotante) para representar los números en el ordenador consta de signo, mantisa y exponente; es decir, un número diferente de cero queda representado de la forma MATH donde la mantisa $m$ es un número comprendido entre $0$ y $1$ cuyo primer dígito es no nulo; $b$ es la base ($10$ cuando estamos en el sistema decimal), que no es necesario dar; y $e$ es el exponente. Pensemos, por ejemplo, en $-0.23015*10^{23}$ o $0.95469*10^{-5}$. Para representar en binario el signo, la mantisa y el exponente de un número se usa una cadena de bit $.$ Un bit es un elemento físico que admite dos estados que, por tanto, podemos hacer corresponder al $0$ y al $1.$ Una cadena de bit se puede hacer corresponder a una de $0$ y $1$, y así obtener números en binario. Con un bit también podemos representar los dos signos $+$ y $-$. Matemáticamente, la cadena de bit que se usa para representar a un número equivale a una sucesión finita de ceros y unos. El primero de ellos se emplea para representar el signo del número, la mayor parte de los siguientes para la mantisa y los restantes para el exponente. No es imprescindible reservar una posición para el signo del exponente, ya que podemos realizar una traslación interna de los exponentes y, así, representar exponentes tanto positivos como negativos de, aproximadamente, igual valor. Por ejemplo, con siete cifras en binario se puede conseguir representar los números del $0$ al $127$; restando a todos ellos $64$, tendríamos como exponentes todos los números comprendidos entre $-64$ y $+63$.

Una situación normal es disponer de $64$ bit, de los que uno se dedica al signo, $52$ a la mantisa y $11$ al exponente. Los exponenetes posibles irán desde el $0$ hasta el MATH, de tal forma que, si se contempla un desplazamiento automático de $1024$, tendríamos como posibles exponentes los valores enteros comprendidos entre $-1024$ y $1023$, y el rango de números positivos con que podríamos trabajar va, aproximadamente, de MATH hasta MATH. La precisión es de MATH, es decir, de unas $16$ cifras decimales. Si tras una operación se obtiene un valor superior a $2^{1023}$ se produce un desbordamiento (overflow). Un valor inferior a $2^{-1024}$ es tomado como $0$ (underflow)

Sólo un conjunto finito de números puede ser representado en el ordenador de forma exacta. En nuestro ejemplo, el siguiente a MATH es MATH Cualquier otro número intermedio entre estos dos tiene que ser aproximado por uno de ellos (redondeo). Si se aproxima por el más cercano se dice que el redondeo es simétrico; si se cortan todos los dígitos de la representación que sobrepasan la dimensión de la mantisa el redondeo es por corte. Por ejemplo, en binario, con $52$ dígitos, el redondeo simétrico de MATH es $2^{-1}+2^{-52}$, mientras que por corte es redondeado a $2^{-1}=0.5$.

A partir de ahora suponemos que todos los números están escritos en sistema decimal y no lo indicaremos.

El número máximo de decimales que soporta el ordenador en la representación interna de un número real se denomina precisión. Supongamos que la precisión es $k\in \QTR{Bbb}{N}$; entonces todo número real se puede representar de la forma MATHDe hecho, supuesto que el redondeo sea simétrico, ese mismo número máquina representa a todos los números del intervalo MATHPor ejemplo, si la precisión es $9$, el número MATH representa indistintamente a todos los números del intervalo MATHLa diferencia entre cualquiera de estos números y MATH es el error de redondeo correspondiente. Lógicamente, los errores de redondeo se producen tanto en la entrada y salida de los datos, como en todas las operaciones intermedias. Por ejemplo, si la precisión del ordenador fuese de $5$ cifras y hubiese que multiplicar los números $a=0.42183$ y $b=0.38124$, entonces el valor MATH tiene que aproximarse por $0.16082$ si redondea de forma simétrica , o por $0.16081$ si se redondea por corte.

Reorganizando la memoria del ordenador se puede aumentar la precisión en los redondeos. De hecho, desde los primeros ordenadores se pudo trabajar en doble precisión. Programas como Mathematica van mucho más allá y son capaces de aumentar la precisión sin más límite que la memoria del ordenador.

Algunos conceptos elementales

A continuación se comentan brevemente algunos conceptos básicos.

Error absoluto y error relativo

Si $x$ es la solución exacta de un problema y $x^{*}$ es una aproximación de la misma, entonces $a_x=x-x^{*}$ es el error absoluto cometido. Si no tenemos más información, el error absoluto no dice nada por sí solo; por ejemplo, el error absoluto en una pesada fue de $300$ gramos. Cuando $x\neq 0$ se denomina error relativo a MATH. Como, en general, no se conoce $x$, el error relativo se suele aproximar por el valor MATH. El error relativo tiene que ser muy próximo a cero para que la aproximación sea buena; por ejemplo, si el peso exacto debía ser $500$ gramos, se habría cometido un error muy importante y el error relativo sería MATH. Sin embargo, el error relativo sería $0.00003$ y el error absoluto despreciable si el peso exacto fuese $10.000$ Kg.

En la práctica no se suelen conocer los errores absoluto y relativo, sino cotas de los mismos, $\delta _x$ y $\varepsilon _x$, de tal forma que se verifican las desigualdades MATH y MATH $\varepsilon _x$. También se suele notar MATH y MATH, o a la inversa MATH y MATH.

Propagación del error relativo en las operaciones elementales

Supuesto que se conozcan cotas para los errores absoluto y relativo de unas cantidades $x$, $y$, $\ldots $ es interesante conocer una cota del resultado obtenido al efectuar una operación elemental ($+,-,*,/$) con las mismas. Es inmediato comprobar que MATH ya que, si $x=x^{*}+a_x$ e $y=y^{*}+a_y$, entonces MATH

En cuanto al error relativo en la suma o la resta, de MATH e MATH se deduce que MATH; por tanto, supuesto que $x+y\neq 0$, la división por $x+y$ conduce a la igualdad MATH de tal forma que, si tanto $x$ como $y$ son de igual signo, se tiene MATH es decir, MATH Pero si $x$ e $y$ tienen signos opuestos y valores muy similares el error relativo se puede disparar. Es el efecto de la cancelación que se produce al restar cantidades casi iguales que comentaremos después.

En cuanto al producto y al cociente, supondremos que el producto de $r_x$ por $r_y$ y los cuadrados de ellos son cantidades despreciables.

Para la multiplicación se tiene que MATH por tanto, el error relativo en la multiplicación está acotado (aproximadamente) por MATH

De forma análoga, para el cociente se llega a la misma cota: MATH

En lo que se refiere a la evaluación en un punto $x$ de la función MATH, que suponemos derivable con continuidad, del desarrollo de Taylor se deduce que MATHDe forma análoga, si $Y$ es una función de varias variables,MATH, digamos MATH, una cota del error absoluto es MATH

Un desarrollo teórico más extenso sobre propagación de los errores en los operaciones escapa del interés de este libro. Una introducción suele encontrarse en muchos libros de Análisis Numérico; en general, una acotación del error cometido después de miles o millones de operaciones suele ser un valor muy pesimista, ya que supondría que el error es siempre de igual signo y se acumula, hecho que en la práctica rara vez ocurre.

Dígitos decimales significativos

Si $d$ es el mayor entero para el cual MATHse dice que $x^{\ast }$ aproxima a $x$ con $d$ dígitos decimales significativos. Por ejemplo, para MATH se tiene que $e^{\ast }=2.7182$ es una aproximación con cuatro dígitos significativos ya que MATH

Error por truncamiento (o discretización) y por cancelación

Cuando en un problema se dispone de una sucesión de valores, MATH, convergente a la solución del mismo y se usa un método numérico que emplee un valor $x_n$ como aproximación de la solución, el error cometido se llama error de truncamiento o truncatura. Por ejemplo, si consideramos la aproximación MATH, cometemos el error de truncadura MATH Análogamente, a veces se tiene fórmulas procedentes de despreciar el último término de un desarrollo finito de Taylor. El error cometido es por truncamiento. Por ejemplo, como MATH, entonces al aproximar la derivada de $f$ en $a$ por MATH supuesto que no se cometan errores al evaluar la función $f$ en $a+h$ y en $a$, ni en la resta, ni en la división, el error cometido es MATH y corresponde al truncamiento efectuado al despreciar el último término del desarrollo de Taylor.

Pero el ejemplo anterior, correspondiente a la aproximación de la derivada, es apropiado para mostrar la aparición de otros errores. Cuando $h$ es pequeño las cantidades MATH y $f\left( a\right) $ son muy similares, cada vez tendrán mayor número de cifras coincidentes y, por tanto, al efectuar la resta se pierden cifras significativas. Es lo que se denomina error por cancelación. Siempre que se restan cantidades muy próximas se produce un error por cancelación, que puede pasar a a ser un error muy grande cuando va multiplicado por un valor grande, como en el ejemplo anterior, en el que el factor $\frac 1h$ puede ser grande. Aparecen en la fórmula de derivación numérica dos errores con comportamiento diferente respecto del valor de $h$: el de truncamiento, que decrece cuando $h$ disminuye (prácticamente proporcional a $h$), y el de cancelación y redondeo en el numerador (digamos $\varepsilon $) que va dividido por $h$, es decir, inversamente proporcional a $h$. No por disminuir $h$ más y más el error va a ser menor. Al contrario, llegará un momento en el que aumentará rápidamente debido al término MATH. Pero llega otro momento en el que los valores de MATH y $f\left( a\right) $ son idénticos para el ordenador, de tal forma que, aunque disminuya $h$, la aproximación de la derivada será siempre $0$ y el error será MATH.

Para algunos problemas podemos organizar los cálculos de forma que se evite el error por cancelación. Por ejemplo, para todo $x\geq 0$, la expresión MATH es equivalente a MATH, cuyo valor es próximo a $\frac 12\sqrt{x}$; sin embargo, para $x$ suficientemente grande la primera expresión da cero (operando en modo científico).

Orden de aproximación MATH

Muchas veces disponemos de una expresión que nos da, para cualquier $h>0$, un valor $A\left( h\right) $ como aproximación del valor exacto $p$ y se sabe que existe una constante $C>0$ tal que MATH Entonces se dice que el error absoluto MATH es un infinitésimo de orden $n$, y se expresa MATH. Para valores positivos y pequeños de $h$ podemos considerar MATH

Orden de exactitud $n$

En ocasiones la expresión $A\left( h\right) $ es una combinación de datos relativos a una función $f$ (por ejemplo, una combinación de valores de $f$ y de algunas de sus derivadas), mientras que el valor $p$ correponde a una manipulación de $f$ que sólo puede realizarse de forma exacta para ciertas funciones (al menos las de tipo polinómico). En tal caso, los polinomios constituyen unas funciones test para medir la bondad de la aproximación, y se dice que $A\left( h\right) $ tiene orden de precisión $n$ si es exacta para las funciones $1$, $x$, $\ldots $, $x^n$ y comete error distinto de cero para la función $x^{n+1}.$ Por ejemplo, si usamos MATH para aproximar el valor MATH, cuando MATH se obtiene MATH y para MATH resulta que MATH, mientras que si MATH, se verifica que MATH Por tanto el orden de exactitud es $1$.

Condicionamiento

Podemos identificar un problema como una entrada de datos, $x$, que tras sufrir un proceso, $P$, da una salida $y$, lo que se nota $y=Px$. La entrada $x$ puede ser un número, un vector, una matriz, etc. Lo mismo es válido para $y$. ?`Cómo se mide la sensibilidad de $y$ frente a pequeñas variaciones de $x$, supuesto que el proceso $P$ se realiza con precisión infinita? En el caso más habitual en el que tanto $x$ como $y$ son distintos de cero, las variaciones se consideran en términos relativos. El factor que las relaciona es el condicionamiento. Veamos algún ejemplo.

Supongamos que el proceso consiste en evaluar una función real, $f$, en un punto, $x$, de modo que MATH. Si $f$ es derivable, entonces MATH Como suponemos que $x\neq 0\neq y$, podemos escribir MATH Esto sugiere definir el condicionamiento de $f$ en $x$, MATH, por MATH Mide la perturbación relativa de $y$ frente a una perturbación relativa de $x$.

Cuando $f$ es una función de $\QTR{Bbb}{R}^m$ en $\QTR{Bbb}{R}^n$ el condicionamiento de cada una de las componentes $y_\nu $ respecto de cada una de las variables $x_\mu $, MATH, es MATH

Para analizar la sensibilidad de $y$ frente a una perturbación del vector $x$ hemos de considerar normas vectoriales y matriciales compatibles. Por ejemplo, si MATH representa la matriz jacobiana de $f$, empleando la norma del máximo concluimos que MATH de forma que MATH es el condicionamiento en este caso.

Uno de los problemas en los que con más frecuencia hay que estimar el condicionamiento es el correspondiente a la resolución de sistemas lineales. Como se verá en aquel tema, el condicionamiento de una matriz viene dado por el producto de la norma de la matriz por la norma de su inversa, esto es MATH. Un sistema puede ser de dimensión muy pequeña y tener un mal condicionamiento. Por ejemplo, consideremos la matriz MATH Su determinante vale $-1$. Su inversa es MATH Así, con cualquier norma matricial $\left\| A\right\| $ será grande y MATH muy grande. Concretamente, si sumamos los valores absolutos de los elementos de cada fila y nos quedamos con el mayor valor de los obtenidos (norma del máximo para matrices), se tiene MATH y MATH con lo cual MATH (un valor enorme). Los condicionamientos bajos al trabajar con matrices son los que valen poco más de la unidad.

En todos los casos, si el condicionamiento (o número de condición) es alto se dice que el problema es mal condicionado y debemos cuidar mucho la precisión de los datos porque un pequeño error inicial en los mismos (por ejemplo, un redondeo) puede dar lugar a grandes diferencias en los resultados.

Estabilidad

La estabilidad se refiere a la forma en que repercuten los errores iniciales y los ocurridos durante el proceso en los resultados finales. Así, cuando estos errores se van propagando y ampliando de forma muy considerable se dice que el proceso es inestable, mientras que cuando los posibles errores iniciales y los que surgen en el proceso van amortiguándose se dice que el proceso es estable. El proceso de resolución de un problema puede ser bien condicionado y ser inestable. La inestabilidad afecta a todo el proceso; por tanto, es posible que cambiando el algoritmo de resolución pueda evitarse. Por ejemplo, es fácil demostrar por inducción que la sucesión de valores MATH puede generarse indistintamente a partir de los siguientes algoritmos:

  1. $s_0=1$, MATH $n\geq 1$.

  2. $s_0=1$, $s_1=\frac 12$, MATH, $n\geq 2$.

Sin embargo, con el segundo (operando con 6 cifras de precisión) el decimosexto término es $s_{15}=-113$, frente al valor MATH.

Análogamente, la sucesión MATH puede generarse a partir del algoritmo $s_0=1$, $s_1=\frac 13$, MATH, $n\geq 2$, que también es inestable.

Combinando estabilidad y condicionamiento sólo estamos seguros de que los resultados son correctos cuando aplicamos un algoritmo estable a un problema bien condicionado. En cualquiera de los otros tres casos no estamos seguros de que los resultados sean fiables.

Análisis a posteriori del error

Muchas veces no sabemos si el problema que estamos resolviendo está bien o mal condicionado y si el algoritmo usado es estable o no. Lo único que sabemos es que a partir de unos datos produce unos resultados. Si los resultados fuesen precisos seguramente el problema estaría bien condicionado y el método de resolución sería estable. Uno de los procedimientos más simples para estimar la veracidad de los resultados es el método de permutación-perturbación, que consiste en resolver el mismo problema varias veces, permutando las operaciones (o parte del proceso) que sea posible y redondeando los datos iniciales y los resultados intermedios cada vez, por exceso o por defecto, de forma aleatoria. Entonces la media de los resultados, $m$, dividida por la desviación típica de los mismos, $d$, sirve como estimador estadístico, de forma que $\log \frac md$ tiene que tener valores muy similares (que disten menos de $0.5$) cuando se aplica a dos tandas de perturbaciones. En la práctica, si después de efectuar tres o cuatro permutaciones-perturbaciones los resultados difieren poco en términos relativos es casi seguro que son fiables.

A veces el condicionamiento es tan malo que no es necesario recurrir a este estimador, sino que al aplicar varias perturbaciones los resultados son tan dispares que resulta evidente.

Por ejemplo, hemos aplicado a cada elemento de la matriz anterior, $A$, una perturbación aleatoria comprendida entre $-0.05$ y $0.05$ (sumando a cada elemento de $A$ el valor -0.05 Random[]) y hemos resuelto el sistema $Ax=b$, siendo MATH La solución exacta es MATH. En cada una de las cinco perturbaciones realizadas al menos una de las componentes del vector solución ha variado en más de $100$ veces, incluso ha cambiado de signo. Concretamente, si la matriz perturbada es MATH entonces MATH.

Ni que decir tiene que en un caso como el anterior cualquier error pequeño en los coeficientes del sistema puede llevarnos a una solución completamente errónea.

Bibliografía

A. Aubanell, A. Benseny, A. Delshams, Útiles Básicos de Cálculo Numérico, Labor, Barcelona, 1993.

W. Gautschi, Numerical Analysis. An Introduction, Birkhäuser,1997.

J. Martínez, A, Martínez y A. Lirola, Errores de redondeo: problemas y soluciones, Revista de Informática y Automática, Julio de 1985.

J. H Mathews and K. D. Fink, Métodos Numéricos con Matlab, Prentice Hall, 1999