Wikilibros

Implementación de algoritmos/Matemáticas/Método de Newton

< Implementación de algoritmos‎ | Matemáticas

El algoritmo del Método de Newton consiste en una relación de recurrencia.

Sea f: [a, b] -> R función derivable definida en el intervalo real [a, b]. Empezamos con un valor inicial x0 y definimos para cada número natural n

Donde f ' denota la derivada de f.

Ejemplo editar

Consideremos el problema de encontrar un número positivo x tal que cos(x) = x3. Podríamos tratar de encontrar el cero de f(x) = cos(x) - x3.

Sabemos que f '(x) = -sin(x) - 3x2. Ya que cos(x) ≤ 1 para todo x y x3 > 1 para x>1, deducimos que nuestro cero está entre 0 y 1. Comenzaremos probando con el valor inicial x0 = 0,5

 

Los dígitos correctos están subrayados. En particular, x6 es correcto para el número de decimales pedidos. Podemos ver que el número de dígitos correctos después de la coma se incrementa desde 2 (para x3) a 5 y 10, ilustando la convergencia cuadrática.

En pseudocódigo, esto es: 

function newtonIterationFunction(x) {
    return  x - (cos(x) - x^3) / (-sin(x) - 3*x^2)     
}

var x := 0,5

for i from 0 to 99 {
    print "Iteraciones: " + i
    print "Valor aproximado: " + x
    xold := x
    x := newtonIterationFunction(x) 
    if x = xold {
        print "Solución encontrada!"
        break
    }
}

Código en Python editar

Programa escrito en Python correspondiente al ejemplo f(x) = cos(x) - x3. 

import math

def newtonIterationFunction(x): 
    return  x - (math.cos(x) - x**3) / (-math.sin(x) - 3*(x**2))     

 
x = 0.5
 
for i in range(100):
    print "Iteraciones: ",str(i),"Valor aproximado: ", str(x)
    xold = x
    x = newtonIterationFunction(x) 
    if x == xold:
        print "Solución encontrada!"
        break

Para ejecutar escribe en una terminal: 

$ python nombre-del-programa.py


Código en C editar

Programa escrito en C correspondiente al ejemplo f(x) = cos(x) - x3. 

#include <stdio.h>
#include <math.h>

double FuncionIteracion(double x);

int main() {
   int i;
   double x, xvieja;

   xvieja = 0.5; // valor inicial

   printf("Iteracion      valor\n");
   for(i=0; i<8; i++){
      x = FuncionIteracion(xvieja);
      xvieja = x;
      printf("%5.0d  %20.12f\n", i+1, x);
   }
   
   return 0;
}

double FuncionIteracion(double x)
{
   double valorFuncion, funcionOriginal, funcionDerivada;

   funcionOriginal = cos(x) - pow(x, 3);
   funcionDerivada = -sin(x) - 3*pow(x, 2);
   valorFuncion = x - funcionOriginal / funcionDerivada;

   return valorFuncion;
}

Para compilar en GNU/Linux con compilador de GNU, se escribe en una terminal: 

$ gcc programa.c -lm -o programa

donde la opción -lm le indica al compilador que utilice la biblioteca de matemáticas.

Código en Matlab editar

Programa escrito en Matlab para hallar las raíces usando el método de NEWTON-RAPHSON 

% Al escribir la función, usar x como variable.
x0=input('Ingrese el valor inicial: ');
tol=input('Ingrese el porcentaje de error: ');
f=input('Ingrese la función: ');
i=1;
fx(i)=x0;
syms x;
f1=subs(f,x,fx(i));
z=diff(f);
d=subs(z,x,fx(i));
ea(1)=100;
while abs(ea(i))>=tol;
    fx(i+1)=fx(i)-f1/d; f1=subs(f,x,fx(i+1)); d=subs(z,x,fx(i+1));
    ea(i+1)=abs((fx(i+1)-fx(i))/fx(i+1)*100);
    i=i+1;
end
fprintf('i fx(i) Error aprox (i) \n');
for j=1:i;
    fprintf('%2d \t %11.7f \t %7.3f \n',j-1,fx(j),ea(j));
end

El programa siguiente hace el cálculo para una superficie. 

syms x1
syms x2
syms x3
V=['sin(x1)+2^x2+log(x3)-7';'3*x1+2*x2-x3^3+1 ';'x1+x2+x3-5 '];
%se calcula el jacobiano:
DV(1,:)=[diff(V(1,:),x1), diff(V(1,:),x2),diff(V(1,:),x3)];
DV(2,:)=[diff(V(2,:),x1), diff(V(2,:),x2),diff(V(2,:),x3)];
DV(3,:)=[diff(V(3,:),x1), diff(V(3,:),x2),diff(V(3,:),x3)];
%se da el valor de partida:
x1=0;
x2=4;
x3=2;
x_1o=[x1;x2;x3];
%se calcula H en ese punto
Vo(1,:)=eval(V(1,:));
Vo(2,:)=eval(V(2,:));
Vo(3,:)=eval(V(3,:));
%Se calcula el Hessiano en ese punto
DV1=eval(DV);
%se calcula la Inversa del Hessiano en ese punto
DV_1=DV1^-1;
%se calcula el siguiente valor de iteración
x_1=[x1;x2;x3]-DV_1*Vo;
%cantidad de iteraciones maxima:
n=50; 
%se define a = n, si se cumple condicion de error antes, cambia
a=n; 

for i=1:n 
%error relativo entre aproximaciones sucecivas
    er=norm(x_1-x_1o)/norm(x_1);
    if er<.0001
        a=i;
       'break;' end
   
    x1=x_1(1);
    x2=x_1(2);
    x3=x_1(3);
    x_1o=[x1;x2;x3];
    Vo(1,:)=eval(V(1,:));
    Vo(2,:)=eval(V(2,:));
    Vo(3,:)=eval(V(3,:));
    DV1=eval(DV);
    DV_1=DV1^-1;
    x_1=[x1;x2;x3]-DV_1*Vo;  
    
end
    a
    x_1
Obtenido de «https://es.wikibooks.org/w/index.php?title=Implementación_de_algoritmos/Matemáticas/Método_de_Newton&oldid=369294»