10.4: Método de Newton

En esta sección necesitaremos lo siguiente:

using Pkg; Pkg.add("Plots");

using Plots

Resolver una ecuación es una parte muy importante de las matemáticas y otros campos científicos. Sin embargo, encontrar soluciones generales a las ecuaciones no es muy fácil. Considera una ecuación cúbica como

\[15x^3 −143x^2 +226x+280=0 \label{cubic1} \]

En el espíritu de la fórmula cuadrática, hay una fórmula cúbica. Gran parte de la página de wikipedia pasa tiempo resolviendo el cúbico con todas las posibilidades. En definitiva, no es muy fácil. En lugar de usar tal fórmula, un enfoque más robusto es resolverlo numéricamente.

Esta ecuación cúbica en el lado izquierdo de\ ref {cubic1} en realidad factores, pero encontrar esos factores es bastante difícil de hacer en general. Una gráfica de la función se da con

f(x) = 15x^3-143x^2 +226x+280
plot([f,x->0],-2.5,10,label=["f(x)" "x-axis"])

Los tres puntos de intersección entre la curva roja (\(x\)eje -eje) y la línea azul (la función,\( y = f (x)\)) son las tres raíces.

El método de Newton comienza con una “conjetura” en la raíz y luego la refina. Deja\(x_0\) ser la conjetura, entonces

Luego nos fijamos en usar el método de Newton para resolver esto. Lo siguiente devuelve otro valor para la raíz.

\[x_1 = x_0 - \frac{f(x_0)}{f'(x_0)} \label{newton} \]

Cuando el método de Newton funciona (y a menudo lo hace bien), el valor de\(x_1\) está más cerca de la raíz verdadera que\(x_0\). Tenemos la función\(f(x)\) definida en julia arriba, y podemos definir la derivada con

df(x)=45x^2-286x+226

Waiting for kernel...

\[ df(x) = 45x^2-286x+226 \label{diff:cubic1} \]

Entonces podemos hacer un paso del método de Newton con

x0=0
x1 = x0 - f(x0)/df(x0)

Luego podemos hacer otro paso del método de Newton usando x1 como valor de entrada y

x2 = x1 - f(x1)/df(x1)

y señalar que parece que los puntos se están acercando a la raíz (ver la trama anterior) que está entre 0 y -1. Algunas iteraciones más de esto deberían obtener un valor muy preciso de la raíz.

Ejercicio

Realiza dos pasos más del método de Newton. Eso es encontrar x3 y x4.

# insert your code here.

Una función del método de Newton

Debería parecer natural tratar de escribir una función para realizar el método de Newton durante algún número de pasos. Caminaremos por el proceso de pensamiento en el desarrollo de la función. Comenzamos con la plantilla para esta función. Las cosas que se necesitan para encontrar una raíz usando el método de Newton es 1) la función 2) la derivada y 3) la suposición inicial. Por lo tanto comenzaremos con

function newton(f::Function,df::Function,x0::Number) 

end

Aquí hay consideraciones al construir la función:

• Necesitarás hacer los dos pasos anteriores muchas veces por lo que necesitarás un bucle. Como no sabes cuántas veces necesitarás correr por el bucle usa un bucle while y tu condición debe ser que los dos pasos x0 y x1 estén separados entre sí.

• Comprobando si dos números de punto flotante son iguales generalmente no es una buena idea, porque tienen que ser iguales a todos los bits, así que en su lugar ejecutaremos el bucle while mientras que la diferencia en los números x0 y x1 son mayores que algunos por defecto (como 10−6).

Aquí hay más un marco de la función:

function newton(f::Function, df::Function, x0::Number)
  x1 = x0 - f(x0)/df(x0)
  while abs(x1-x0)>1e-6 # while the two numbers are larger than 10^(-6)
    x0 = x1
    x1 = x0 - f(x0)/df(x0) 
  end
  x1
end

Waiting for kernel...

Usando esto ahora podemos llamar a esta función como

newton(f,df,0)

Solo para simplificar esto, definiremos dx como -f (x0) /df (x0), que es la distancia entre valores x sucesivos, así que usaremos esto para determinar cuándo detener el bucle:

function newton(f::Function, df::Function, x0::Real) 
  local dx=f(x0)/df(x0)
  while abs(dx)>1e-6
    x0 = x0-dx
    dx = f(x0)/df(x0)
  end
  x0
end

Y solo para asegurar que la función siga devolviendo el mismo resultado,

newton(f,df,0)

 

Uso de la diferenciación automática

Si ha utilizado Sistemas de Álgebra Computacional como Maple o Mathematica, sabe que las computadoras tienen la capacidad de diferenciarse. Julia no tiene esta capacidad incorporada, aunque tiene alguna capacidad de hacer algunas de las características de estos programas. Existe un sistema llamado diferenciación automática que calculará la derivada exacta a una función en un punto dado. Es decir, si tienes una función f (x) y un número a, te dará f′ (a). El paquete se llama ForwardDiff y es posible que deba agregarlo y luego

using Pkg; Pkg.add("ForwardDiff");
using ForwardDiff

Por ejemplo, si define:

g(x) =x^2

g (generic function with 1 method)

ForwardDiff.derivative(g,3), ForwardDiff.derivative(g,-1) 

Método de Newton con diferenciación automática

Usando el paquete ForwardDiff, podemos simplificar el método de Newton y no requerir que necesitemos encontrar el derivado al llamarlo. La siguiente función

function newton(f::Function, x0::Real)
  local dx = f(x0)/ForwardDiff.derivative(f,x0) 
  while abs(dx)>1e-6
    x0 = x0 - dx
    dx = f(x0)/ForwardDiff.derivative(f,x0)
  end
  x0
end

Waiting for kernel...