CalculoNumericoV4

Vista previa del material en texto

Calculo Numérico - UNJu - Facultad de Ingeniería �San Salvador de Jujuy
RAÍCES DE ECUACIONES NO LINEALES
Uno de los problemas básicos del análisis numérico es el llamado problema de búsqueda de raíces, que consiste
en encontrar los valores de la variable x que satisfacen la ecuación f(x) = 0 . A una solución de este problema
se le llama un cero de f o una raíz de f(x) = 0.
Podemos clasi�car las ecuaciones de acuerdo al tipo de función que es f(x). A grandes rasgos se tiene:
Ecuaciones

Lineales
No lineales
 Algebraicas
{
Racionales
Irracionales
Trascendentes
Una ecuación lineal en la variable x es una ecuación que puede escribirse en la forma ax+ b = 0, donde a y
b son constantes que generalmente llamamos parámetros y a 6= 0.
Algunos ejemplos de ecuaciones no lineales son:
� Ecuación algebraica racional 9x3x−1 = 2 +
3
3x−1
� Ecuación algebraica irracional, la variable x sometida a la operación de radicación x
√
2.1−0.5x
(1−x)
√
1.1−0.5x −3.69 =
0 0 < x < 1
� Ecuación trascendente, incluye funciones trigonométricas, exponenciales, logarítmicas y otras menos fa-
miliares tg(x) = tgh(2x)
La razón principal para resolver ecuaciones no lineales por medio de métodos computacionales radica en la di�-
cultad de encontrar una solución por métodos convencionales. Por su parte, excepto para muy pocos problemas,
la solución analítica de las ecuaciones polinomiales existe sólo hasta el orden cuatro, pero no existen métodos
generales para arribar a las soluciones en forma exacta para órdenes superiores. Por lo tanto, las raíces de esas
ecuaciones no lineales se obtienen mediante los métodos del análisis numérico.
Los métodos usuales para la obtención de una aproximación numérica a una solución o raíz de f(x) = 0
consisten en procesos iterativos en los que se parte de un valor inicial x0 de la raíz buscada α y, se usa cierta
relación de recurrencia para generar una secuencia de aproximaciones sucesivas x1, x2,..., xn,... que convergen
al límite xn, generalmente por métodos analíticos, muchas veces con ayuda de grá�cos.
El problema se plantea de la siguiente manera: Dada f : R→ R(o bien f : [a, b] �R) se quiere encontrar α tal
que f(α) = 0. El cálculo aproximado de raíces puede dividirse en dos etapas:
� Se separan las raíces, es decir se busca un subintervalo de [a, b] que contenga una y sólo una raíz de f.
Para asegurar la existencia de al menos una raíz en el intervalo propuesto se usa el Teorema de Bolzano.
Para asegurar que no hay más de una raíz se usa el Teorema de Rolle, es decir, se veri�ca que la derivada
primera no cambie de signo en dicho intervalo.
� Se aplica un método para aproximar la raíz aislada.
MÉTODO GRÁFICO
Además de la utilidad para determinar valores iniciales, también son útiles para visualizar las propiedades
de las funciones y el comportamiento de los métodos numéricos.
Consiste en gra�car la función y observar en donde cruza el eje x. Este punto, que representa el valor de x
para el cual f(x) = 0, proporciona una aproximación inicial de la raíz.
Ejemplo:
Empléense grá�cas para obtener una raíz aproximada de la función f(x) = e−x − x
Solución:
1
Calculo Numérico - UNJu - Facultad de Ingeniería �San Salvador de Jujuy
Figure 1:
Por medio del programa GnuPlot para generar grá�cas de funciones se obtiene la grá�ca de la función
f(x) = e−x − x. Ver Figura 1
Un vistazo a la grá�ca proporciona una estimación aproximada de la raíz de 0.57 que se acerca a la raíz
exacta de 0.56714328...., que se debe determinar con métodos numéricos. La validez de la estimación visual se
puede veri�car sustituyendo su valor en la ecuación original para obtener: f(0.57) = e−0.57− 0.57 = −0.0045 la
cual se acerca a cero.
MÉTODOS NUMÉRICOS PARA APROXIMACIÓN DE RAÍCES
Se dividen en dos categorías generales:
1. Métodos cerrados que usan intervalos: Bisección (intervalo medio), Regula Falsi (falsa posición).
Requieren un intervalo de x que contenga a la raíz, siempre son convergentes, pero la velocidad de con-
vergencia puede ser demasiada lenta.
2. Métodos abiertos: Iteración de punto �jo, Newton Raphson, Secante.
Requieren información únicamente de un punto, o de dos pero que no necesariamente encierran a la raíz,
para extrapolar una nueva aproximación a la raíz. La convergencia es más rápida pero existe también la
posibilidad de divergencia.
MÉTODO DE LA BISECCIÓN (DEL INTERVALO MEDIO)
Es el más simple, aunque también el más seguro y sólido para encontrar una raíz en un intervalo dado, donde
se sabe que existe dicha raíz. Se apoya en la idea geométrica del teorema de Bolzano: Dada f(x) = 0 si f es tal
que es monótona y continua en [a, b], f(a) y f(b) tienen signos distintos entonces existe, por lo menos un α, a<
α <b, tal que f(α) = 0. En general puede decirse que en el intervalo [a, b] existe un número impar de raíces.
El método requiere de dividir repetidamente a la mitad los subintervalos de [a, b] y, en cada paso, localizar
la mitad que contiene a la aproximación de la raíz xi. Ver Figura 2
2
Calculo Numérico - UNJu - Facultad de Ingeniería �San Salvador de Jujuy
Figure 2:
Este proceso se repite hasta que se veri�que algún criterio de parada.
Algoritmo:
1. Veri�car el teorema de Bolzano f(a) . f(b) < 0, si no se veri�ca el teorema no signi�ca que la función no
posea raíces y ante esto se sugiere buscar otro intervalo [a, b].
2. Establecer tolerancia de error ξ.
3. Calcular una nueva aproximación a la raíz en el punto medio del intervalo:
xi =
(a+b)
2
4. Calcular f(xi), si f(xi) = 0 entonces la raíz es igual a xi y el proceso �naliza.
5. Calcular f(a) . f(xi), si f(a) . f(xi) < 0 entonces la raíz se encuentra en [a,xi], se actualiza el extremo
superior haciendo b=xi.
6. Calcular f(a) . f(xi), si f(a) . f(xi) > 0 entonces la raíz se encuentra en [xi,b], se actualiza el extremo
inferior haciendo a=xi.
7. Veri�car criterio de parada |xi−xi−1xi | < ξ, si esta condición se veri�ca �naliza el proceso iterativo y se toma
a xi como la aproximación a la raíz. Por lo contrario si la condición de parada no se veri�ca se vuelve a
repetir el proceso desde el punto 3.
Al repetir este proceso, el tamaño del intervalo con la raíz se vuelve cada vez más pequeño. En cada paso, se
toma el punto medio del intervalo como la aproximación más actualizada de la raíz. Se genera una sucesión
x1 =
a+b
2 ∈ [a1, b1], x2 ∈ [a2, b2], x3 ∈ [a3, b3],. . . , donde cada intervalo [an, bn] mide la mitad del anterior.
b1 − a1 = b−a2
b2 − a2 = b2−a22 =
b−a
4
bn − an = ... = b−a2n
Además
a ≤ a1 ≤ a2 ≤ ... ≤ b
b ≥ b1 ≥ b2 ≥ ... ≥ a
Entonces an y bn son sucesiones monótonas y acotadas y en consecuencia convergen, es decir existen los
límites:
3
Calculo Numérico - UNJu - Facultad de Ingeniería �San Salvador de Jujuy
liman
n→∞ y
lim bn
n→∞
y como |bn − an| ≤ b−a2n → 0, se tiene
liman
n→∞ =
lim bn
n→∞ = α
En cada paso se veri�ca f(an) . f(bn) ≤ 0 y tomando límite ( f continua) resulta f
′′
(α) ≤ 0
Entonces α es la raíz buscada pues cumple, f(α)=0
En este método el cálculo de cotas de error es muy simple. Por cota de error entendemos un número que
acote superiormente, en módulo, el error máximo que podríamos llegar a cometer cuando nos quedamos con uno
de los puntos medios de los intervalos construidos mediante el algoritmo, en vez de con la solución del problema.
El error se puede acotar de la siguiente forma.
Tenemos:
xn =
(an−1+bn−1)
2
Entonces |α− xn| ≤ bn−1−an−12n =
b−a
2n
Criterio de parada y estimación del error
Como se planteó anteriormente el método termina cuando se alcance un error más bajo, por ejemplo, al
0.1%, pero esta estrategia resulta inconveniente, ya que la estimación del error en el se basa en el conocimiento
del valor verdadero de la raíz de la función. Éste no es el caso de una situación real, ya que no habría motivo
para utilizar el método si se conoce la raíz.
Por lo tanto, sin la necesidad del conocimiento previo de la raíz, se puede calcular el error relativo porcentual
ξade la siguiente manera:
ξa = |xi−xi−1xi| . 100%
Donde xi es la raíz en la iteración actual y xi−1 es el valor de la raíz en la iteración anterior. Se utiliza
el valor absoluto, ya que por lo general importa sólo la magnitud de ξasin considerar su signo. Cuando ξa es
menor que un valor previamente �jado ξ, termina el cálculo.
Inconvenientes del Método
El método de bisección tiene inconvenientes importantes. Converge muy lentamente (o sea, N, número de
iteraciones, puede ser muy grande antes que, x - xn , sea su�cientemente pequeño) y, una buena aproximación
intermedia puede ser desechada sin que nos demos cuenta. Además, hay que tener en cuenta que en el caso
de existir más de una raíz (siempre en número impar) en el intervalo; el método sólo encuentra una de ellas,
desechándose las otras. Entonces puede darse la situación paradójica que, se encuentre una raíz y sin embargo
no sea esta la solución más conveniente al problema.
Aplicaciones del Método
La bisección suele recomendarse para encontrar un valor aproximado de la raíz, y luego este valor se re�na
por medio de métodos más e�caces. La razón es que la mayor parte de los otros métodos requieren un valor
inicial cerca de una raíz; al carecer de dicho valor pueden fallar por completo.
Converge para cualquier f continua.
Ejemplo 1:
Use el método de bisección para determinar la raíz de f(x) = e−x − x con una exactitud de 10−3.
Solución:
1. Análisis grá�co (Si es posible): Si se analiza La Figura 1 se puede ver que la raíz se encuentra entre 0 y
1, por lo tanto se acota la raiz al intervalo [0,1].
2. Veri�cación inicial de convergencia: Se veri�ca que f(0) . f(1) < 0, por consiguiente la estimación inicial
de la raíz se sitúa en el punto medio de este intervalo.
3. Se establece criterio de parada del método: ξ = 10−3
4. Aplicación iterativa del método mediante software. Ver Figura 3
4
Calculo Numérico - UNJu - Facultad de Ingeniería �San Salvador de Jujuy
Figure 3:
Como se puede observar en la tabla en la iteración 10 se cumple que el error absoluto es menor que ξ ,
por lo tanto terminan las iteraciónes, y así obtenemos como aproximación a la raiz xi = 0.5673828
5. Solución a traves de lenguaje de programación (ForTran 77):
Program Biseccion
real eps ,x,x0 ,x1
integer i, maxitera
x0=0.0
x1=1.0
eps =0.001
maxitera =100
if (f(x0)*f(x1).gt.0) then
print *, 'No hay cambio de signo en el intervalo inicial '
goto 200
endif
dif =10.0
i=1
do while ((dif.gt.eps).and.(i.le.maxitera ))
x=(x0+x1 )/2.0
if((f(x0)*f(x)).le.0) then
x1=x
dif=error(x1,x0)
else
x0=x
dif=error(x0,x1)
endif
write (*,50)i,x0,x1,dif
i=i+1
end do
if (i.le.maxitera) then
if (dif.lt.eps) then
write (*,51) 'Raiz Aproximada=', x
write (*,51) 'f(x)=', f(x)
end if
else
print *, 'No converge en ',maxitera , ' iteraciones '
end if
50 format(I5, 1x, F16.8, F16.8, F16.8)
51 format(A, F16.8)
stop 'Fin del programa '
200 end
real function error(x1,x2)
real x1, x2
error=abs(x1 -x2)
return
5
Calculo Numérico - UNJu - Facultad de Ingeniería �San Salvador de Jujuy
end
real function f(x)
real x
f=exp(-x)-x
return
end
MÉTODO DE LA FALSA POSICIÓN (REGULA-FALSI)
Este método aprovecha la idea de unir los puntos con una línea recta. La intersección de esta línea con el eje
x proporciona una estimación de la raíz. El reemplazo de la curva por una línea recta da una �posición falsa�
de la raíz, de aquí el nombre de método de la regla falsa o en latín regula falsi. Ver Figura 4
Figure 4:
La recta L que une los puntos (a, f(a))con (b, f(b))tiene la ecuacion:
y − f(a) = f(b)−f(a)b−a (x1 − a)
Como x1es el valor de x que cumple y=0, se tiene:
x1 = a− f(a)(b−a)f(b)−f(a)
esta es la fórmula de la regla de la falsa posición.
Algoritmo:
1. Veri�car el teorema de Bolzano f(a) . f(b) < 0, si no se veri�ca el teorema no signi�ca que la función no
posea raíces y ante esto se sugiere buscar otro intervalo [a, b].
2. Establecer tolerancia de error ξ.
3. Calcular una nueva aproximación a la raíz con la fórmula de la falsa posición:
xi = a− f(a)(b−a)f(b)−f(a)
6
Calculo Numérico - UNJu - Facultad de Ingeniería �San Salvador de Jujuy
4. Calcular f(xi), si f(xi) = 0 entonces la raíz es igual a xi y el proceso �naliza.
5. Calcular f(a) . f(xi), si f(a) . f(xi) < 0 entonces la raíz se encuentra en [a,xi], se actualiza el extremo
superior haciendo b=xi.
6. Calcular f(a) . f(xi), si f(a) . f(xi) > 0 entonces la raíz se encuentra en [xi,b], se actualiza el extremo
inferior haciendo a=xi.
7. Veri�car criterio de parada |xi−xi−1xi | < ξ, si esta condición se veri�ca �naliza el proceso iterativo y se toma
a xi como la aproximación a la raíz. Por lo contrario si la condición de parada no se veri�ca se vuelve a
repetir el proceso desde el punto 3.
Ventajas y Desventajas del Método
Comparte con el método del intervalo medio la ventaja de converger en cualquier circunstancia, y su principal
desventaja es la de encontrar sólo un resultado en el caso de raíces múltiples. Sin embargo, este método es más
veloz que su análogo, el del intervalo medio o bisección. Adicionalmente, pueden aparecer extremos �jos, como
muestra la Figura 4, en donde uno de los extremos de la sucesión de intervalos no se mueve del punto original,
por lo que las aproximaciones convergen a la raíz exacta solamente por un lado. En cuyo caso, no siempre,
pueden presentarse situaciones de convergencia rápida.
Ejemplo 2:
Use el método de la falsa posición para determinar la raíz de f(x) = e−x − x con una exactitud de 10−3.
Solución:
1. Análisis grá�co (Si es posible): Si se analiza La Figura 1 se puede ver que la raíz se encuentra entre 0 y
1, por lo tanto se acota la raiz al intervalo [0,1].
2. Veri�cación inicial de convergencia: Se veri�ca que f(0) . f(1) < 0.
3. Se establece condición de parada del método: ξ = 10−3
4. Aplicación iterativa del método mediante software. Ver Figura 5
Figure 5:
Como se puede ver en la iteración 4 se cumple que el error absoluto es menor que ξ , terminan las iteraciónes,
y así obtenemos como aproximación a la raiz xi = 0.5672055
CANTIDAD DE ITERACIONES
Se puede predecir �a priori� el número de iteraciones que se deben realizar con el método de bisección o de
la falsa posición para obtener una aproximación con una presición deseada ξ.
Debido a que en cada iteración se reduce el error a la mitad, la fórmula general que relaciona el error deseado
ξ y el número de iteraciones n es:
ξ = b−a2n
Si ξ es el error deseado, de esta ecuación se despeja n
entonces n =
log( b−aξ )
log(2)
7
Calculo Numérico - UNJu - Facultad de Ingeniería �San Salvador de Jujuy
MÉTODO DE ITERACIÓNES DE PUNTO FIJO (ITERACIONES SUCESIVAS)
La idea es reemplazar la ecuación f(x) = 0 por otra de la forma x = g(x) de manera que la solución de
ésta sea la solución del problema original. Esta transformación se puede llevar a cabo mediante operaciones
algebraicas o simplemente agregando x a cada lado de la ecuación original. A una solución de esta ecuación se
le llama un punto �jo de la función g.
La ecuación x = g(x) proporciona una fórmula para predecir un nuevo valor de x en función del valor anterior
de x. Es decir, dada una aproximación inicial a la raíz, xi, la ecuación puede usarse para obtener una nueva
aproximación xi+1, expresada por las fórmulas iterativas.
xi = g(xi−1), xi+1 = g(xi)
La ventaja de este método consiste en su gran sencillez y �exibilidad para elegir la forma de g(x). Sin
embargo, es muy importante la formación de la función g(x) en la ecuación x = g(x); de las múltiples opciones
que pueden existir, ya que no siempre converge con cualquier forma elegida de g(x) y cuando lo hace, lo hace
más rápido que los métodos cerrados.
Un planteamiento grá�co consiste en separar la ecuación x = g(x) en dos partes, como:
y1 = x y2 = g(x)
estas funciones se pueden gra�car por separado. Los valores de x correspondientes a las intersecciones de
estas funciones representan las raíces de f(x) = 0. Ver Figura 6
Figure 6:
Convergencia
Sea x un punto �jo de una función continua y derivableg(x) tal que | g′(x) |< 1. Entonces el método de
iteraciones de punto �jo xi = g(xi−1) para i=0,1,2,3,... es localmente convergente en el vencidario del punto
�jo x de g(x), es decir la iteración xi = g(xi−1) converge a dicho punto �jo.
En otras palabras, la convergencia se da cuando la magnitud de la pendiente de g(x) es menor que la
pendiente de la recta f(x) = x. Ver Figura 7
8
Calculo Numérico - UNJu - Facultad de Ingeniería �San Salvador de Jujuy
Figure 7:
Algoritmo:
1. Plantear valor inicial x0 (cercano a la raíz).
2. Establecer tolerancia de error ξ.
3. Plantear fórmula de iteración de punto �jo x = g(x).
4. Calcular una nueva aproximación a la raíz con la fórmula:
xi = g(xi−1)
5. Veri�car criterio de parada |xi−xi−1xi | < ξ, si esta condición se veri�ca �naliza el proceso iterativo y se toma
a xi como la aproximación a la raíz. Por lo contrario si la condición de parada no se veri�ca se vuelve a
repetir el proceso desde el punto 4.
Ejemplo 3:
Use el método de iteraciones sucesivas para determinar la raíz de f(x) = e−x−x con una exactitud de 10−3.
Solución:
1. Analisis grá�co (Si es posible): Si se analiza La Figura 1 se puede ver que la raíz se encuentra cerca de
0.6 , por lo tanto se plantea el valor inicial x0 = 0.6.
2. Se establece condición de parada del método: ξ = 10−3
3. Se plantean 2 fórmulas o transformaciones de f(x) :
(a) g1(x) = e−x =⇒ xi = e−xi−1
(b) g2(x) = − log(x)log(e) =⇒ xi = −
log(xi−1)
log(e)
4. Aplicación iterativa del método para g1(x) mediante software. Ver Figura 8
9
Calculo Numérico - UNJu - Facultad de Ingeniería �San Salvador de Jujuy
Figure 8:
Como se puede ver en la iteración 8 se cumple que el error relativo es menor que ξ , terminan las iteraciónes
y obtenemos como aproximación a la raiz xi = 0.56749133. Por lo tanto se puede decir que con la fórmula
iterativa g1(x) el método converge a la raíz.
Se puede observar también que en cada iteración la distancia del valor actual de xi respecto del valor
anterior xi−1 se hace mas pequeña.
5. Aplicación iterativa del método para g2(x) mediante software. Ver Figura 9
Figure 9:
Para este ejemplo en particular en la iteración 8 se observa que se presenta un caso de indeterminación
debido al valor anterior xi−1 es cual resulta ser negativo . Por lo tanto se puede decir que con g2(x) el
método no converge a la raíz.
Se puede observar que en cada iteración la distancia del valor actual de xi respecto del valor anterior xi−1
se hace mas grande.
6. Solución con g1(x) a traves de lenguaje de programación (ForTran 77):
Program PuntoFijo
real eps , xi , xanterior
integer i, maxitera
xi=0.6
eps =0.001
maxitera =100
dif =10.0
i=1
10
Calculo Numérico - UNJu - Facultad de Ingeniería �San Salvador de Jujuy
do while ((dif.gt.eps).and.(i.le.maxitera ))
xanterior=xi
xi=g(xanterior)
dif=error(xi,xanterior)
write (*,50)i, xi, dif
i=i+1
end do
if (i.le.maxitera) then
if (dif.lt.eps) then
write (*,51) 'Raiz Aproximada=', xi
write (*,51) 'f(x)=', f(xi)
end if
else
print *, 'No converge en ',maxitera , ' iteraciones '
end if
50 format(I5, 1x, F16.8, F16 .8)
51 format(A, F16.8)
stop 'Fin del programa '
200 end
real function error(x1,x2)
real x1, x2
error=abs(x1-x2)
return
end
real function g(x)
real x
g=exp(-x)
return
end
real function f(x)
real x
f=exp(-x)-x
return
end
MÉTODO DE LA TANGENTE O DE NEWTON-RAPHSON
El método consiste en empezar con un valor de x0 (cercano a la raíz) y trazar la tangente en el punto
(x0, f(x0)). El punto donde esta tangente cruza al eje x se toma el mismo como la siguiente aproximación. Este
proceso se repite hasta que los valores de xi sucesivos están su�cientemente próximos o el valor de la función
está su�cientemente cerca de cero. Ver Figura 10.
11
Calculo Numérico - UNJu - Facultad de Ingeniería �San Salvador de Jujuy
Figure 10:
Hay por lo menos tres maneras usuales de introducir el método de Newton-Raphson:
1. Se puede derivar geométricamente o técnica grá�ca
2. El uso de la serie de Taylor
3. Obtener el método de Newton a partir de la técnica de iteración de punto �jo.
Vamos a desarrollar 2) y 3).
2) Recuérdese que la serie de Taylor se puede representar como:
f(xi+1) = f(xi) + f
′(xi)(xi+1 + xi) +
f(ξ)
2 (xi+1 − xi)
2
en donde ξ se encuentra en alguna parte del intervalo entre xi y xi+1 . Truncando la serie después de la
primera derivada, se obtiene una versión aproximada:
f(xi+1) ≈ f(xi) + f ′(xi)(xi+1 − xi)
En la intersección con el eje x, f(xi+1) debe ser igual a cero , o:
0 ≈ f(xi) + f ′(xi)(xi+1 − xi)
que se puede resolver para:
xi+1 = xi − f(xi)f ′(xi)
a la que se conoce como fórmula de Newton-Raphson. Observemos que para que la fórmula tenga
sentido f ′(xi) 6= 0.
3) Para resolver una ecuación de la forma f(x) = 0, supongamos que la ecuación f(x) = 0 tiene una solución
x∗ tal que f ′(x∗) 6= 0 .
Consideremos el esquema de punto �jo
xi+1 = g(xi)
con la g de la forma
12
Calculo Numérico - UNJu - Facultad de Ingeniería �San Salvador de Jujuy
g(x) = x− φ(x)f(x)
Donde φes una función arbitraria que se escogerá más adelante.
Si φ(x) está acotada, entonces g(x∗) = x∗ , y para que el procedimiento iterativo derivado de g sea cuadráti-
camente convergente, es su�ciente que g′(x∗) = 0 . Pero
g′(x) = 1− φ′(x)f(x)− f ′(x)φ(x) y
g′(x∗) = 1− f ′(x∗)φ(x∗)
Consecuentemente, g′(x∗) = 0 si y solo si φ(x∗) = 1f ´(x∗)
El proceso iterativo que de�ne esta elección es:
xi+1 = xi − f(xi)f ′(xi)
que es la fórmula de Newton�Raphson.
El método de Newton se emplea ampliamente porque, al menos en la vecindad próxima de una raíz, converge
más rápido que cualquiera de los métodos analizados hasta ahora.
Algoritmo:
1. Plantear valor inicial x0 (cercano a la raíz).
2. Establecer tolerancia de error ξ.
3. Obtener f ′(x).
4. Calcular una nueva aproximación a la raíz con la fórmula de Newton-Raphson:
xi+1 = xi − f(xi)f ′(xi)
5. Veri�car criterio de parada |xi−xi−1xi | < ξ, si esta condición se veri�ca �naliza el proceso iterativo y se toma
a xi como la aproximación a la raíz. Por lo contrario si la condición de parada no se veri�ca se vuelve a
repetir el proceso desde el punto 4.
Ejemplo 4:
Use el método de iteraciones sucesivas para determinar la raíz de f(x) = e−x−x con una exactitud de 10−3.
Solución:
1. Analisis grá�co (Si es posible): Si se analiza La Figura 1 se puede ver que la raíz se encuentra cerca de
0.6 , por lo tanto se plantea el valor inicial x0 = 0.6.
2. Se establece condición de parada del método: ξ = 10−3
3. Se encuentra f ′(x):
f ′(x) = −e−x − 1
4. Aplicación iterativa del método mediante software. Ver Figura 11
Figure 11:
13
Calculo Numérico - UNJu - Facultad de Ingeniería �San Salvador de Jujuy
De esta manera, el método converge rápidamente a la raíz.
Desventajas del método de Newton-Raphson
Aunque el método de Newton en general es muy e�ciente, hay situaciones en que se comporta de�cientemente.
Ejemplo:
Determínese la raíz positiva de f(x) = x10 − 1 usando el método de Newton�Raphson con un valor inicial
de x = 0.5.
Solución:
1. Se establece condición de parada del método: ξ = 10−3
2. Se encuentra f ′(x):
f ′(x) = 10x9
3. Aplicación iterativa del método mediante software. Ver Figura 12
Figure 12:
......
De esta forma, después de la primera predicción de�ciente, el método converge a la raíz 1, pero con una
velocidad muy lenta. Además de la convergencia lenta, debida a la naturaleza de la función, se pueden
originar otras di�cultades, como:
14
Calculo Numérico - UNJu - Facultad de Ingeniería �San Salvador de Jujuy
Figure 13:
La �gura 13 muestra el caso donde un punto de in�exión ocurre en la vecindad de una raíz. Nótese que las
iteraciones que empiezan en x0 divergen progresivamente de la raíz.
Figure 14:
En la �gura 14 se observa la tendencia del método de Newton�Raphson a oscilar alrededor de un punto
mínimo o máximo local.
15
Calculo Numérico - UNJu - Facultad de Ingeniería �San Salvador de Jujuy
Figure 15:
Se observa en la �gura 15 que un valor inicialcercano a una raíz puede saltar a una posición varias raíces
lejos. Esta tendencia de alejarse del área de interés se debe a que se encuentran pendientes cercanas a cero.
Una pendiente cero causa una división por cero en la fórmula de Newton- Raphson. Grá�camente, esto signi�ca
que la solución se dispara horizontalmente y jamás toca el ejex. La única solución en estos casos es la de tener
un valor inicial cercano a la raíz. Este conocimiento, de hecho, lo proporciona el conocimiento físico o quimico
del problema o mediante el uso de herramientas tales como las grá�cas que proporcionan mayor claridad en el
comportamiento de la solución.
16
Calculo Numérico - UNJu - Facultad de Ingeniería �San Salvador de Jujuy
INTERPOLACIÓN
El concepto de interpolación surge, por ejemplo, cuando disponemos de datos que provienen de mediciones
experimentales o estadísticos, puesto que queremos determinar la evolución general de estos datos con el ob-
jetivo de estimar/predecir los valores que no conocemos. Por ejemplo, esto ocurre si tenemos partes de una
imagen fotográ�ca y queremos reconstruir la imagen completa. En otras palabras, buscamos una función (lla-
mada función interpolante) que toma valores predeterminados en algunos puntos. Notemos que otra aplicación
de la interpolación es la aproximación de funciones dadas. Normalmente se utilizan funciones de un tipo pre-
determinado (polinomios, funciones trigonométricas, etc) dando lugar a diferentes métodos de interpolación.
Estudiaremos la interpolación polinómica.
INTERPOLACIÓN POLINOMIAL
Una de las más y bien conocidas clases de funciones reales de variable real es la clase de los polinomios
algebraicos, o sea, el conjunto de funciones de la forma
f(x) = a0 + a1x+ a2x
2 + a3x
3 + .....+ anx
n
donde n es un entero no negativo y a0, a1, a2..., an son constantes reales. Una razón primordial de su
importancia es que aproximan uniformemente funciones continuas; esto es, dada una función de�nida y continua
en un intervalo cerrado, existe un polinomio que está tan cerca de la función dada como se desee.
Teorema de Aproximación de Weierstrass:
Si f está de�nida y es continua en [a, b], dado ζ > 0, existe un polinomio P, de�nido en [a, b], con la
propiedad que
| f(x)�P (x) | < ζ para toda x ∈[a, b]. Ver Figura 16.
Figure 16:
Otro aspecto importante para considerar a los polinomios en la aproximación de funciones es que es sencillo
determinar la derivada y la integral inde�nida de cualquier polinomio y el resultado es otra vez un polinomio.
Por estas razones, los polinomios se usan con frecuencia para aproximar otras funciones quedamos se conoce o
se supone son continuas.
POLINOMIO DE LAGRANGE
Planteo del problema
Sea f(x) la función que se quiere interpolar y se supone conocida en un conjunto de puntos x0, x1, x2, ..., xn:
y0 = f(x0)
17
Calculo Numérico - UNJu - Facultad de Ingeniería �San Salvador de Jujuy
y1 = f(x1)
y2 = f(x2)
. .
. .
yn = f(xn)
La interpolación de Lagrange consiste en encontrar un polinomio de grado n, P (x) (polinomio de interpo-
lación de Lagrange), que pase por los puntos dados. Dicho polinomio cumple las condiciones:
P (x0) = y0
P (x1) = y1
P (x2) = y2
. .
. .
P (xn) = yn
Caso Lineal
Vamos a comenzar por plantearnos el caso de interpolar mediante una línea recta que une 2 puntos cua-
lesquiera. La ecuación de la recta que pasa por los puntos (x0, y0) y x1, y1 es la que presentamos a continuación:
y = P (x) = y0 + (y1 − y0) (x−x0)(x1=x0)
Vamos a tratar de reescribir la misma expresión tal cual lo hizo Lagrange:
P1(x) = y0
(x−x1)
(x0=x1)
+ y1 (x=x0)(x1=x0) = L0(x)y0 + L1(x)y1
Con L0(x) =
(x=x1)
(x0=x1)
y L1(x) =
(x=x0)
(x1=x0)
, llamados coe�cientes de Lagrange, que cumplen:
Cuando x = x0, L0(x0) = 1 mientras que L1(x0) = 0 =⇒ P1(x0) = y0
Cuando x = x1, L0(x1) = 0 mientras que L1(x1) = 1 =⇒ P1(x1) = y1
Caso General: Polinomio de grado n
Teorema 1: Si x0, x1, x2...., xn son (n+1) números diferentes y f es una función cuyos valores están dados
en estos puntos, entonces existe un único polinomio P de grado n con la propiedad de que :
f(xk) = P (xk) para cada k = 0, 1, 2, ....., n
Este polinomio está dado por:
P (x) = f(x0)L0(x) + f(x1)L1(x) + f(x2)L2(x) + ....+ f(xn)Ln(x) =
n∑
i=0
f(xi)Li(x)
Demostración:
Por una serie de n+1 puntos pasa un polinomio de grado n que, lo podemos expresar en función de sus raíces
y tiene la siguiente forma:
Pn(x) = a0(x − x1)(x − x2)...(x − xn) + a1(x − x0)(x − x2)...(x − xn) + a2(x − x0)(x − x1)(x − x3)...(x −
xn)+ ...+ ai(x−x0)(x−x1)...(x−xi−1)(x−xi+1)...(x−xn)+ ...+ an−1(x−x0)(x−x1)...(x−xn−2)(x−xn)+
an(x− x0)(x− x1)...(x− xn−2)(x− xn−1)
Los coe�cientes del polinomio se determinan haciendo cumplir las condiciones:
P (x0) = y0
18
Calculo Numérico - UNJu - Facultad de Ingeniería �San Salvador de Jujuy
P (x1) = y1
P (x2) = y2
. .
. .
P (xn) = yn
obtenemos:
a0 =
f(x0)
(x0=x1)(x0=x2)...(x0=xn)
a1 =
f(x1)
(x1=x0)(x1=x2)...(x1=xn)
a2 =
f(x2)
(x2−x0)(x2−x1)...(x2−xn)
en general:
ai =
f(xi)
(xi=x0)(xi=x1)...(xi=xi−1)(xi=xi+1)...(xi=xn)
= f(xi)
n∏
j = 0
j 6= i
1
(xi−xj)
Reemplazando en el polinomio Pn(x), obtenemos la fórmula del polinomio de Lagrange
Pn(x)=
n∑
i=0
f(xi)Li(x), con Li(x) =
n∏
j = 0
j 6= i
(x−xj)
(xi−xj)
El siguiente paso consiste en calcular un término residual o cota para el error involucrado en la aproximación
de una función mediante un polinomio interpolante. Esto se hace en el teorema siguiente:
Teorema 2:
Si x0, x1, x2....., xn son puntos distintos en [a, b] y si f es derivable hasta el orden (n+1) en [a, b], entonces,
para cada x en [a, b], existe un número ξ(x) en (a, b) tal que:
f(x) = P (x) + f
(n+1)(ξ(x))
(n+1)! (x-x0)(x− x1)...(x− xn)
donde P (x) es el polinomio interpolante. El segundo término corresponde a la fórmula del error. Esta
fórmula es un resultado teórico importante, su uso práctico está restringido a funciones cuyas derivadas tengan
cotas conocidas.
Ejemplo 1:
La tabla muestra los valores de una función en diversos puntos. Compararemos las aproximaciones a
f(1.5)obtenidas con varios polinomios de Lagrange.
x f(x)
1, 0 0,7651977
1,3 0,6200860
1,6 0,4554022
1,9 0,2818186
2,2 0,1103623
Aproximación por Interpolación lineal:
Como x = 1, 5 se encuentra entre 1,3 y 1,6, el polinomio lineal utilizará x0 = 1, 3 y x1 = 1, 6.
P1(x) =
1∑
i=0
f(xi)Li(x) = f(x0)L0(x) + f(x1)L1(x)
P1(1.5) = 0.6200860
(1.5−1.6)
(1.3−1.6) + 0.4554022
(1.5−1.3)
(1.6−1.3) = 0.5102968
Aproximación por polinomio de grado 2:
Suponemos que x0 = 1.3, x1 = 1.6 y x2 = 1.9
P2(x) =
2∑
i=0
f(xi)Li(x) = f(x0)L0(x) + f(x1)L1(x) + f(x2)L2(x)
19
Calculo Numérico - UNJu - Facultad de Ingeniería �San Salvador de Jujuy
P2(1.5) = 0.6200860
(1.5−1.6)
(1.3−1.6)
(1.5=1.9)
(1.3=1.9) + 0.4554022
(1.5−1.3)
(1.6−1.3)
(1.5=1.9)
(1.6−1.9) + 0.2818186
(1.5=1.6)
(1.9−1.6)
(1.5=1.3)
(1.9=1.3) = 0.5112857
Aproximación por polinomio de grado 3:
Suponemos que x0 = 1.3, x1 = 1.6, x2 = 1.9 y x3 = 2.2
P3(x) =
3∑
i=0
f(xi)Li(x) = f(x0)L0(x) + f(x1)L1(x) + f(x2)L2(x) + f(x3)L3(x)
P3(1.5) = 0.6200860
(1.5−1.6)
(1.3−1.6)
(1.5−1.9)
(1.3−1.9)
(1.5−2.2)
(1.3−2.2) + 0.4554022
(1.5−1.3)
(1.6−1.3)
(1.5−1.9)
(1.6−1.9)
(1.5−2.2)
(1.6−2.2)+
0.2818186 (1.5−1.3)(1.9−1.3)
(1.5−1.6)
(1.9−1.6)
(1.5−2.2)
(1.9−2.2) + 0.1103623
(1.5−1.3)
(2.2−1.3)
(1.5−1.6)
(2.2−1.6)
(1.5−1.9)
(2.2−1.9) = 0.5118302
Aproximación por polinomio de grado 4 (se utilizan todos los valores de la tabla):
x0 = 1.0, x1 = 1.3, x2 = 1.6, x3 = 1.9 y x4 = 2.2
P4(x) =
3∑
i=0
f(xi)Li(x) = f(x0)L0(x) + f(x1)L1(x) + f(x2)L2(x) + f(x3)L3(x) + f(x4)L4(x)
P4(1, 5) = 0, 7651977
(1,5−1,3)
(1,0−1,3)
(1,5−1,6)
(1,0−1,6)
(1,5−1,9)
(1,0−1,9)
(1,5−2,2)
(1,0−2,2) + 0, 6200860
(1,5−1,0)
(1,3−1,0)
(1,5−1,6)
(1,3−1,6)
(1,5−1,9)
(1,3−1,9)
(1,5−2,2)
(1,3−2,2)+
0, 4554022 (1,5−1,0)(1,6−1,0)
(1,5−1,3)
(1,6−1,3)
(1,5−1,9)
(1,6−1,9)
(1,5−2,2)
(1,6−2,2) + 0, 2818186
(1,5−1,0)
(1,9−1,0)
(1,5−1,3)
(1,9−1,3)
(1,5−1,6)
(1,9−1,6)
(1,5−2,2)
(1,9−2,2)+
0, 1103623 (1,5−1,0)(2,2−1,0)(1,5−1,3)
(2,2−1,3)
(1,5−1,6)
(2,2−1,6)
(1,5−1,9)
(2,2−1,9) = 0, 5118200
Ejemplo 2:
Determinar el polinomio interpolador de Lagrange para la función f(x) = sen(x) en los siguientes puntos
x0 =
π
4 , x1 =
π
3 , x2 =
π
2 . Estimar sen(1.2). ¾Serviría este polinomio para evaluar sen(3)?.
Determinación de valores yi:
x f(x)
π
4
√
2
2
π
3
√
3
2
π
2 1
Aproximación por polinomio de grado 2:
P2(x) =
√
2
2
(x=π3 )
(π4−
π
3 )
(x=π2 )
(π4−
π
2 )
+
√
3
2
(x=π4 )
(π3=
π
4 )
(x=π2 )
(π3=
π
2 )
+ 1
(x=π4 )
(π2=
π
4 )
(x=π3 )
(π2=
π
3 )
Realizamos operaciones algebraicas y agrupamos términos:
P2(x) =
√
2
2 (
48(x=π3 )(x=
π
2 )
π2 ) +
√
3
2 (
−72(x=π4 )(x=
π
2 )
π2 ) + 1(
24(x=π4 )(x=
π
3 )
π2 )
El polinomio interpolador con coe�cientes decimales es:
P2(x) = −0, 447100350x2 + 1, 426378617x− 0, 137374388
Puesto que disponemos de la funci�ón original y el polinomio que la interpola, podemos ver sus gr�á�cos y
compararlos. Ver Figura 17. Se puede observar también que si ampliamos el intervalo hasta sobrepasar el
intervalo [π4 ,
π
2 ], la similitud entre la función y el polinomio que la interpola se pierde pues s�ólo se aproximan
en el intervalo de interpolación. El proceso de interpolar por fuera de intervalo se llama extrapolación.
20
Calculo Numérico - UNJu - Facultad de Ingeniería �San Salvador de Jujuy
Figure 17:
Valor estimado de sen(1, 2):
P2(1, 2) =0, 930455446
Valor real de sen(1, 2):
sen(1, 2) = 0, 932039086
Desventajas del método
Las desventajas de la interpolación de Lagrange son las siguientes:
� La cantidad de cálculos necesaria para una interpolación es grande.
� La interpolación para otro valor de x necesita la misma cantidad de cálculos adicionales, ya que no se
pueden utilizar partes de la aplicación previa.
� Cuando el número de datos tiene que aumentar o disminuir, no se pueden utilizar los resultados de los
cálculos previos.
� Aumentar el número de datos en el intervalo no implica mejora en los resultados.
� La evaluación del error no es fácil.
Problemas con interpolación � Fenómeno de Runge
Supongamos que dado un intervalo [a, b] lo vamos subdividiendo en más y más puntos, más concretamente
tomemos:
xi = a+ ih para i = 0, 1, 2, 3, ..., n donde h =
(b−a)
n
y supongamos que construimos con estos puntos el polinomio de interpolación Pn(x) para una función dada
f, esto es, que Pn(xi) = f(xi), para estos n puntos.
¾Se tiene
lim
n→∞
pn(x) = f(x) cuando n tiende a in�nito?
La respuesta es negativa. En realidad, al aumentar el número de puntos se mejora la aproximación en
la parte central del intervalo, pero la diferencia entre la función y el polinomio interpolador puede aumentar
21
Calculo Numérico - UNJu - Facultad de Ingeniería �San Salvador de Jujuy
rápidamente en los extremos. No es bueno hacer demasiado extenso el intervalo de interpolación, ya que además
de aumentar el número de operaciones con la consecuente acumulación de errores, podemos aumentar la pérdida
de precisión en los extremos. Este fenómeno es conocido como fenómeno de Runge.
Ejemplo 3
Dada la función:
f(x) = 11+x2 en el intervalo [−5, 5]. Ver �gura 18.
Figure 18:
Si construimos el polinomio de interpolación Pn(x) en este intervalo, entonces seguro que no hay convergencia
en los puntos donde |x| > 3, 63. El problema parece pues que se tuerce, pero por otro lado se vuelve más
interesante. Resulta que para ciertas funciones, por ejemplo para f(x) = ex, si hay convergencia. Pero para
otras no. El problema con estas últimas es que sus derivadas van creciendo demasiado en el intervalo considerado,
esto es lo que sucede con esta función de Runge, que parecía al principio bastante inofensiva. En las �guras
siguiente se demuestra este comportamiento.
Demostración del fenómeno de Runge en interpolación por Lagrange.
Caso1: Polinomio de Grado 2
Figure 19:
Caso 2: Polinomio de Grado 4
22
Calculo Numérico - UNJu - Facultad de Ingeniería �San Salvador de Jujuy
Figure 20:
Caso 3: Polinomio de Grado 10
Figure 21:
Caso 4: Polinomio de Grado 20
Figure 22:
23
Calculo Numérico - UNJu - Facultad de Ingeniería �San Salvador de Jujuy
Caso 5: Polinomio de Grado 50
Figure 23:
MÉTODO DE DIFERENCIAS DIVIDIDAS
Los métodos para determinar la representación explícita de un polinomio interpolante a partir de datos
tabulados se conocen como métodos de diferencias divididas. Estos métodos se usaron más con propósitos de
cómputo antes de que el equipo de cómputo digital llegara a ser fácilmente disponible. Sin embargo, los métodos
pueden usarse también para derivar técnicas para aproximar las derivadas y las integrales de funciones, así como
para aproximar las soluciones de ecuaciones diferenciales.
Interpolación de Newton en puntos con separación no uniforme
Supongamos que Pn es el polinomio de Lagrange de grado n que coincide con la función f en los números
distintos x0, x1, x2....., xn. Las diferencias divididas de f con respecto a x0, x1, x2....., xn se pueden derivar
demostrando que Pn tiene la representación:
Pn(x) = a0+a1(x−x0)+a2(x−x0)(x−x1)+a3(x=x0)(x=x1)(x=x2)+...+an(x−x0)(x−x1)........(x−xn−1)
con constantes apropiadas a0, a1, a2, ..., an
Evaluando Pn en x0 : a0 = Pn(x0) = f(x0)
Evaluando en x1 : f(x0) + a1(x1 − x0) = Pn(x1) = f(x1)⇒ a1 = f(x1)−f(x0)x1−x0
Introducimos lo que se conoce como notación de diferencia dividida.
Diferencias divididas de orden cero de la función f:
f [x0] = f(x0), f [x1] = f(x1), f [x2] = f(x2), ........., f [xn] = f(xn)
Diferencias divididas de orden 1 de la función f:
f [x0, x1] =
f [x1]−f [x0]
x1−x0 , f [x1, x2] =
f [x2]−f [x1]
x2−x1 ,....,f [xi, xi+1] =
f [xi+1]−f [xi]
xi+1−xi
Diferencias divididas de orden 2 de la función f:
f [x0, x1, x2] =
f [x1,x2]−f [x0,x1]
x2−x0 , f [x2, x3, x4] =
f [x3,x4]−f [x2,x3]
x4−x2 ,....,f [xi, xi+1, xi+2] =
f [xi+1,xi+2]−f [xi,xi+1]
xi+2−xi
Cuando las (k − 1) diferencias divididas
f [xi, xi+1, xi+2, ..., xi+k−1] y f [xi+1, xi+2, ..., xi+k=1, xi+k]
han sido determinadas, la k-ésima diferencia dividida de f relativa a xi, xi+1, xi+2, ...., xi+k, está dada por:
f [xi, xi+1, ..., xi+k−1, xi+k] =
f [xi+1,xi+2,...,xi+k]−f [xi,xi+1,...,xi+k−1]
xi+k−xi
24
Calculo Numérico - UNJu - Facultad de Ingeniería �San Salvador de Jujuy
Los coe�cientes a1, a2, a3, ..., an, se pueden expresar en términos de las diferencias divididas
a1 =
f(x1)−f(x0)
x1−x0 = f [x0, x1], a2 = f [x0, x1, x2], a3 = f [x0, x1, x2, x3], a4 = f [x0, x1, x2, x3, x4]....
Sustituyendo en el polinomio interpolante
Pn(x) = a0+a1(x−x0)+a2(x−x0)(x−x1)+a3(x=x0)(x=x1)(x=x2)+...+an(x−x0)(x−x1)........(x−xn−1)
Pn(x) = f [x0] + f [x0, x1](x − x0) + f [x0, x1, x2](x − x0)(x − x1) + f [x0, x1, x2, x3](x=x0)(x=x1)(x=x2) +
...+ f [x0, x1, ..., xn](x=x0)(x=x1)(x=x2)........(x=xn=1)
Pn(x) =
n∑
k=0
f [x0, x1, ..., xk](x− x0)(x=x1)...(x− xk−1)
En forma mas concisa:
Pn(x) =
n∑
k=0
f [x0, x1, ..., xk]
k−1∏
j=0
(x− xj)
Obtenemos la fórmula de interpolación de Newton en diferencias divididas
La determinación de las diferencias divididas para puntos de datos tabulados se bosqueja en la tabla siguiente.
Se podrían encontrar una cuarta diferencia a partir de estos datos.
x f(x) 1ras dif. Div. 2das dif. Div. 3ras dif. Div. 4tas dif. Div.
f[x0, x1, x2, x3, x4] =
f[x1,x2,x3,x4]−f[x0,x1,x2,x3]
x4−x0
Ejemplo 4:
Dada la tabla del ejemplo 1 , encontrar los coe�cientes de la fórmula de las diferencias divididas progresivas
del polinomio interpolante de Newton y aproximar f(1.5) .
Solución:
1. Mediante lenguaje de programación (ForTran 77) se construye la tabla de diferencias divididas de orden
4 . Ver �gura 24.
Figure 24:
2. Como se puede observar los coe�cientes se encuentran a lo largo de la diagonal de la tabla, por lo tanto:
25
Calculo Numérico - UNJu - Facultad de Ingeniería �San Salvador de Jujuy
P4(x) = 0.7651977− 0.4837057(x− 1.0)− 0.1087337(x− 1.0)(x− 1.3) + 0.0660631(x− 1.0)(x− 1.3)(x−
1.6) + 0.0012088(x− 1.0)(x− 1.3)(x− 1.6)(x− 1.9)
3. Aproximación de f(1, 5):
P4(1, 5) = 0, 511787
4. Programa tablade diferencias divididas.
Program tablaDiferenciasDivididas
real F(0:10 ,0:10) , x(0:10)
print *
print *, "TABLA DE DIFERENCIAS DIVIDIDAS"
data NI/4/
data (x(I), I=0 ,4)/1.0 , 1.3, 1.6, 1.9, 2.2/
data (F(I,0),I=0 ,4)/0.7651977 , 0.620086 , 0.4554022 , 0.2818486
+ , 0.1103623/
do K=1, NI
J=NI - K
do I=0, J
F(I,K)=(F(I+1,K-1) -F(I,K-1) ) / (x(I+K) -x(I))
end do
end do
print *
print *, " I X F[x0] F[x0,x1] ..."
do I=0, NI
J=NI -I
write (*,440), I, x(I), (F(I,K),K=0,J)
end do
print *
440 format (1X,I2 ,8F12 .7)
stop
end
26
Calculo Numérico - UNJu - Facultad de Ingeniería �San Salvador de Jujuy
TECNOLOGÍAS LIBRES UTILIZADAS
Plataforma base: GNU/Linux � Debian �wheezy�
Procesador de texto: Lyx/Latex
Lenguaje de Programación: Fortran G77
Generador de Grá�cos: GNUPlot
Editor de Ecuaciones: Lyx/Latex
Editor de Diagramas: Dia
Planilla de Calculo: LibreO�ce Calc
BIBLIOGRAFÍA
Chapra Steven C., Canale Raymond P.; MÉTODOS NUMÉRICOS PARA INGENIEROS. Con apli-
caciones en computadoras personales, 1996, McGraw � Hill/Interamericana de México.
Burden Richard L., Faires J.Douglas; ANÁLISIS NUMÉRICO, 1996, Grupo Editorial Iberoamérica.
Sadosky Manuel, CÁLCULO NUMÉRICO Y GRÁFICO, 1971, Ediciones Librería del Colegio.
Nakamura Shoichiro, MÉTODOS NUMÉRICOS APLICADOS CON SOFTWARE, 1992, Pearson Edu-
cación
27