domingo, 11 de noviembre de 2012

Método de Newton-Raphson

Tal vez, dentro de los métodos para determinar el valor de las raíces de una función ya sea lineal o no lineal, el método de Newton es uno de los más utilizados. Si el valor inicial con el que se aproxima  la raíz es $x_i$, entonces se puede extender una tangente desde el punto $[x_i,f(x_i)]$ en donde $i=0,1,2,..,n$. De esta forma el punto de la tangente que corta al eje $x$ representa una aproximación mejorada de la raíz (ver la figura).

Figura: Ilustración del modo en que el método de Newton-Raphson se aproxima a la solución.


 El método de Newton-Raphson se puede deducir sobre la base de una interpretación geométrica (un método alterno basado en la serie de Taylor). Como podemos apreciar en la figura, la derivada de $x$ es equivalente a:
$$f'(x_i)=\frac{f(x_i)-0}{x_i-x_{i+1}}$$

Que se puede ordenar para obtener:

$$x_{i+1}=x_i-\frac{f(x_i)}{f'(x_i)}$$

Esta es la formula utilizada en el método de Newton-Raphson

Estimación del error

Para poder especificar un criterio de paro en el método de Newton-Raphson se suele utilizar la siguiente expresión para calcular el error en cada iteración:
$$error=\left | \frac{x_{i+1}-x_i}{x_{i+1}}\right | $$

Por otra parte si desarrollamos el modelo por medio de la serie de Taylor, esto nos proporciona una base teórica firme para poder estimar la velocidad de variación del error en cada iteración y así poder verificar la convergencia del método. Veamos:

La expansión de la serie de taylor se expresa de la siguiente manera:
$$f(x_{i+1})=f(x_i)+f'(x_i)(x_{i+1}-x_i)+\frac{f''(\xi)}{2!}(x_{i+1}-x_i)^2+... +\frac{f^n(\xi)}{n!}(x_{i+1}-x_i)^n\label{ref3}$$

En donde $\xi$ se encuentra en alguna parte en el intervalo entre $x_{i+1}$ y $x_i$. Si truncamos la serie de Taylor después del segundo termino derivado obtenemos la siguiente aproximación:

$$f(x_{i+1}) \cong f(x_i)+f'(x_i)(x_{i+1}-x_i)$$

En a intersección con el eje $x$ debemos tener $f(x_{i+1}))=0$, por lo tanto:

$$0 = f(x_i)+f'(x_i)(x_{i+1}-x_i)$$

Por lo que obtenemos la formula utilizada en el Método de Newton-Raphson
\begin{equation}
x_{i+1}=x_i-\frac{f(x_i)}{f'(x_i)} \label{ref1}
\end{equation}
Ahora si analizamos el error. Para esto hacemos $x_{i+1}=x_r$ en la ecuación \eqref{ref3}, en donde $x_r$ representa el valor real de la raiz y a la vez $f(x_r)=0$ :

$$0=f(x_i)+f'(x_i)(x_r-x_i)+\frac{f''(\xi)}{2!}(x_r-x_i)^2 \label{ref4}$$

y truncando  la ecuación \eqref{ref3} hasta el tercer termino:

$$f(x_{i+1}) \cong f(x_i)+f'(x_i)(x_{i+1}-x_i)+\frac{f''(\xi)}{2!}(x_{i+1}-x_i)^2 \label{ref2}$$

Luego restamos \eqref{ref4} de \eqref{ref2}:


$$0=f'(x_i)(x_r-x_{i+1})+\frac{f''(\xi)}{2!}(x_r-x_i)^2 \label{ref5}$$

Ahora sabiendo que el error entre la aproximación del método y el valor real es igual a la direfencia entre $x_{i+1}$ y $x_r$:

$$E_{t,i+1}=x_r-x_{i+1}$$
entonces:

$$0=f'(x_i)E_{t,i+1}+\frac{f''(\xi)}{2!}E_{t,i}^2 \label{ref6}$$

Si se supone que el método efectivamente converge a una solución, se podría decir entonces que tanto $x_i, \xi \longrightarrow x_r$, de esta manera la ecuación \eqref{ref6} quedaría así:

$$E_{t,i+1}=\frac{-f''(x_r)}{2f'(x_r)}E_{t.i}^2 \label{ref7}$$



Por lo que podemos ver en la ecuación \eqref{ref7} el error es proporcional al cuadrado del error anterior, este es uno de los motivos por los cuales en el caso de que la solución converja lo hace de una manera bastante rápida, sólo después de un par de iteraciones. A este comportamiento se le conoce como convergencia cuadrática.

A continuación les dejo un código que hice en MATLAB con la implementación del método.

Function Newton-Raphson

Function Newton-Raphson


%%%/////////////////////////////////////////////////
% /          Función Newton-Raphson            ///
%/        Desarrollada por Jorge De La Cruz     ///
%/                                             ///
%////////////////////////////////////////////////
% La función 'fun' debe ser introducida como un string
% x: representa el valor inicial necesario para iniciar
% las iteraciones.
% tol: La tolerancia deseada.
% imax: máximo número de iteraciones, evita que el
% programa se quede en un bucle infinito en caso de que
% el método no converja.
function [xr]=newton(fun,x,tol,imax)
err=100;
n=1;
while (tol<=err)&&(n<imax)
    xr = x-eval(fun)/(eval(diff(sym(fun))));
    err=abs((xr-x)/xr);
    error=err*100;
    x=xr;
    fprintf('Iteración: %d, Error: %f, Raiz: %f\n',n,error,xr)
    n=n+1;
end
x=-2:0.001:xr+1;
y=eval(fun);
plot(x,y)
title('METODO DE NEWTON-RAPHSON')
xlabel('$x$','interpreter','latex','fontsize',18)
ylabel('$y$','interpreter','latex','fontsize',18)
grid on
end
Y el código en Scilab sería:
0001  //Función a la cual se le desea determinar la raiz, fun.sci
0002  function y=fun(x)
0003      y=exp(-x)-x;
0004  endfunction

0001  //FUnción Newton-Raphson en otro archivo llamado newton.sci
0002  function xr=newton(d, tol, imax)
0003      err=100;
0004      n=1;
0005      while (tol<=err)&(imax>n)
0006          xr=d-(fun(d)/(derivative(fun,d)));
0007          n=n+1;
0008          err=abs((xr-d)/xr);
0009          e=err*100;
0010          d=xr;
0011          printf("iter: %i, error: %f, raiz: %f\n",n,e,xr)
0012      end
0013      x=xr-2:0.001:xr+2;
0014      y=fun(x);
0015      plot(x,y)
0016  endfunction

Licencia de Creative Commons
Newton-Raphson by Jorge De La Cruz is licensed under a Creative Commons Reconocimiento 3.0 Unported License.

11 comentarios:

harald dijo...

en que sectores del area de la ingenieria mecanica el metodo se hace presente? cual es su importancia para el area?

Unknown dijo...

Hola Herald, en muchos sectores. Es un método muy utilizado (y sus variantes) en la solución numérica de ecuaciones en el área de dinámica. En la solución de sistemas de ecuaciones diferenciales. En donde se buscan las raíces de un sistema de ecuaciones para saber los polos por ejemplo. Determinar puntos de bifurcación, etc etc

Anónimo dijo...

Buenas noches, Como sé por ejemplo que imax me está cumpliendo la funcion de parar las iteraciones cuando se presente algun bucle sin salida?
Gracias

Anónimo dijo...

no funciona en matlab :( me has decepcionado george.

Anónimo dijo...

a mi tambien me has decepcionado :(

Anónimo dijo...

Estamos todos tristes ya que su codigo no ha funcionado en matlab en nuestro grupo de trabajo. que dios no le bendiga hasta comprobar su error.

Anónimo dijo...

Dependíamos de su tutorial. Nos has traicionado.....

Anónimo dijo...

calmensen, jorge enrike es dios y nadie le dice qe su codig o es divino y majestuoso. su codigo de matlab s el pan q komemos tdos los dias. qien se mta con jorge, se mete con migo. pinches webones.

Anónimo dijo...

my nuded teets in my bio
http://www.elrincondelcura.es/

Anónimo dijo...

clique aquí para el código correcto: http://www.pornogay.com

Jorge Enrique dijo...

Hola Chicos!
Disculpen por no haberles respondido antes, lo que pasa es que ya no le estoy dando mantenimiento a este blog. Con referencia al problema que están teniendo, es porque en las versiones más actuales de Matlab la nomenclatura para declarar las variables simbólicas a cambiado. He corregido el código para que se pueda ejecutar en Matlab R2017b.

mi identidad la pueden verificar en: https://keybase.io/jdelacruz26
y el nuevo sitio en donde he empezado a alojar entradas de blogs sobre distintos temas es en: https://jdelacruz26.github.io/blog/

tengo mucho material que añadir pero poco tiempo, si desean algo en especial y si les puedo ayudar, con mucho gusto lo haré.

P.D: el nuevo código lo pueden encontrar aquí: https://jdelacruz26.github.io/code/2013/11/01/first-blog.html

Saludos!