El Método de Newton-Raphson.
Si tenemos
una función f(x) continua y cerca de
una raíz p. Si la derivada f’(x)
existe, entonces puede utilizarse para desarrollar algoritmos que produzcan
sucesiones {pk} que converjan a p más rápidamente que los algoritmos de
bisección.
Consideremos
el caso de una función como la que se muestra en la figura.
Dos iteraciones del método de
Newton.
En estas
figuras se muestra dos iteraciones del método. En la figura de la izquierda
mostramos la línea recta que es tangente a la función f(x) (en negro), note que la línea recta cruza el eje x en x = 1. A la derecha tenemos la segunda
iteración tomando como valor inicial x=1.
Note como poco a poco se acerca a la solución.
La sucesión
que nos lleva a la solución esta dada por los puntos {p0, p1, p2, …, pk}. La
pendiente de la línea recta es
m = (0 – f(p0))/(p1
– p0)
Por otro lado sabemos, del cálculo diferencial,
que la pendiente de la línea tangente a una función es la primer derivada
valuada en ese punto. Así:
m = f’(p0)
Uniendo las ecuaciones tenemos
f’(p0) = (0
– f(p0))/(p1 – p0)
p1= p0 - f(p0)/ f’(p0)
De manera iterativa podemos hacer
p2= p1 - f(p1)/ f’(p1)
p3= p2 - f(p2)/ f’(p2)
pk+1= pk - f(pk)/ f’(pk)
Ejemplo.
Hacer un algoritmo iterativo que permita hacer el cálculo de la raíz cuadrada de A.
Para este caso nuestra función a
resolver es f(x) = x2-A.
La solución cuando f(x)=0 es x = A-0.5.
f(x) = x2-A
f’(x) =2 x
pk+1= pk - f(pk)/ f’(pk)
pk+1= pk - (pk-2-A))/(2 pk)
pk+1=( pk + A/ pk)/2
Los cálculos numéricos suponiendo A=5 son:
|
k |
pk |
|
0 |
2.0000 |
|
1 |
2.2500 |
|
2 |
2.2361 |
|
3 |
2.2361 |
|
4 |
2.2361 |
|
5 |
2.2361 |
|
6 |
2.2361 |
|
7 |
2.2361 |
|
8 |
2.2361 |
|
9 |
2.2361 |
|
10 |
2.2361 |
Ejemplo
Calcular los ceros de la función x-cos(x) utilizando el algoritmo de regula falsi en el intervalo [0,1].
|
k |
pk |
|
0 |
0.0000 |
|
1 |
1.0000 |
|
2 |
0.7504 |
|
3 |
0.7391 |
|
4 |
0.7391 |
|
5 |
0.7391 |
|
6 |
0.7391 |
|
7 |
0.7391 |
Método de Newton Raphson para sistemas no lineales
Consideremos el sistema
Utilizando la serie de Taylor podemos hacer una aproximación lineal
Si escribimos el sistema original como una función vectorial V=F(X), entonces la matriz jacobiana J(x,y) es el análogo bidimensional de la derivada. La aproximación lineal queda como:
donde
Entonces nuestra formulación bidimensional queda como:
y la actualización de la variable la hacemos:
Ejemplo.
Resolver el siguiente sistema de ecuaciones dado por
f1(x,y)
= x2 – 2x – y +0.5
f2(x,y)
= x2 + 4y2 – 4
El jacobiano es
Considerando como valores iniciales [0, 1] tenemos:
Primer iteración
cuya solución es Dx= 0.25 y Dy =0
Segunda iteración
cuya solución es Dx= -0.0274
y Dy =0.0061
Las demás iteraciones la resumimos es:
|
k |
x |
y |
|
0 |
0.0000 |
1.0000 |
|
1 |
-0.2500 |
1.0000 |
|
2 |
-0.2260 |
0.9939 |
|
3 |
-0.2222 |
0.9938 |
|
4 |
-0.2222 |
0.9938 |
|
5 |
-0.2222 |
0.9938 |
|
6 |
-0.2222 |
0.9938 |
|
7 |
-0.2222 |
0.9938 |
La implementación en Matlab es:
function Newton_Raphson(x0)
for it=1:10
dx = inv(J(x0))*f(x0);
x = x0 - dx;
x0 = x;
end;
function y = f(x)
y(1) = x(1)^2 - 2*x(1) - x(2) +0.5;
y(2) = x(1)^2 + 4*x(2)^2 -4;
y =
y';
end;
function y = J(x)
y(1,1) = 2*x(1) - 2;
y(1,2) = -1;
y(2,1) = 2*x(1);
y(2,2) = 8*x(2);
end;