El método de Newton e implementación en Python

El método de Newton es un algoritmo numérico mediante el cual se puede buscar las raíces de las funciones que emplea la derivada para iterar hacia la solución. Obteniendo generalmente el resultado en menos pasos que los métodos basados en la bisección o la secante. La implementación del método de Newton requiere conocer la derivada de la función, por lo que esta ha de ser derivable. Por otro lado, a diferencia de los métodos de la bisección y la secante, no requiere conocer un intervalo donde la función cambia de signo.

Formulación del método de Newton

Supongamos que se tiene una función f(x) de la que se conoce su derivada f'(x) en el entorno de un punto x_0. Bajo estos supuestos se puede aproximar el valor de la función en torno a este punto con la tangente. Lo que lleva a una ecuación lineal en x_0 tal como la siguiente

y = f(x_0) + f'(x0) (x - x_0)

Ahora, al igualar esta ecuación a cero, se puede obtener una primera aproximación de la posición de la raíz de la función despejando el valor dex. Esto es, partiendo de

0 = f(x_0) + f'(x0) (x_1 - x_0)

Se puede obtener

En Analytics Lane

La caverna del consumo, o cómo Greenspan miraba calzoncillos para hacer política monetaria – El bestiario de los indicadores económicos absurdos (parte 1)

x_1 = x_0 - \frac{f(x_0)}{f'(x_0)}

Al ser la primera aproximación de la raíz se ha denotado este valor con x_1. Si esta solución se usa en un proceso iterativo entonces se puede obtener una serie que se aproxima a la raíz de la ecuación

x_{n+1} = x_n - \frac{f(x_n)}{f'(x_n)}

Implementación del método de Newton

Los pasos para la implementación la búsqueda de la raíz de una función mediante el método de Newton se puede resumir en los siguientes pasos:

1. Evaluar la función f(x) en el punto x_n, si el valor absoluto es menor que valor \epsilon, este valor se puede considerar como una solución aproximada.
2. Evaluar la derivada de la función f'(x) en el punto x_n. Cuando el valor es igual a cero, nos encontramos en un punto de inflexión, no se puede continuar ya que la tangente es horizontal. Por lo que nos encontraremos en un punto donde el algoritmo no puede converger.
3. Obtener el valor x_{n+1} mediante la expresión x_{n+1} = x_n - \frac{f(x_n)}{f'(x_n)}.
4. Volver al punto 1 hasta que el algoritmo converja o se alcance un límite de iteraciones fijado.

Implementación en Python

En base a los pasos de la sección anterior se puede implementar este algoritmo en Python usando un código como el siguiente.

def newton(fun, der, x_n, epsilon=1e-6, steps=50):   
    for n in range(steps + 1):
        # Evaluación de la función para ver si el resultado es válido
        f_x = fun(x_n)
        if abs(f_x) < epsilon:
            return x_n
        
        # Evaluación de la derivada
        d_f = der(x_n)
        if d_f == 0:
            print('Error la derivada es cero')
            return None
        
        # Estimación del siguiente punto
        x_n = x_n - f_x / d_f
    
    print('Se ha alcanzado el límite de iteraciones')
    return None

Siendo esta una función que tiene cinco parámetros de entrada de los que solamente tres son obligatorios. En concreto necesita la función a evaluar fun, su derivada der y un punto de inicio x_n. Además, cuenta con dos parámetros opciones mediante los que se pueden modificar los criterios de parada del algoritmo. epsilon permite fijar el error máximo de la raíz y steps el número de iteraciones máximas.

Validación de la implementación

Para validar la implementación se puede usar la función basada en un polinomio con la que se evaluó el método de la secante. Pudiéndose comprobar que obtiene más fácilmente el resultado correcto.

fun = lambda x: x**2 + 2*x - 8
der = lambda x: 2*x + 2

newton(fun, der, 10, steps=10) # 2.0000001002152237

El punto de inicio es un valor importante, en funciones que tiene más de una raíz generalmente se va a localizar la más cercana al punto de inicio. Por ejemplo, una función basada en el seno como \sin\left(\frac{\pi}{2}\right) que tiene las raíces en 2 n \pi, siendo n un número entero.

import math

fun = lambda x: math.sin(x / 2)
der = lambda x: math.cos(x / 2) / 2

newton(fun, der, 2, steps=10) # 5.847132402824612e-13
newton(fun, der, 5, steps=10) # 6.28318530721986

En donde se puede apreciar que, en el primero de los casos, cuando el punto de inicio es 2, se obtiene como resultado un valor próximo a 0. Por otro lado, al cambiar el punto de inicio a 5 el resultado que se obtiene es 2 \pi .

Conclusiones

La implementación del método de Newton en Python es algo relativamente sencilla si se conoce la derivada de la función a evaluar. Ofreciendo grandes ventajas frente a otros métodos. Por ejemplo, respecto a los basados en la bisección o la secante no requiere conocer previamente un intervalo en el que se encuentra una raíz.

Imagen de 849356 en Pixabay