Método de Newton para calcular raíces

José A. Alonso

04-11-2023

Los ceros de una función pueden calcularse mediante el método de Newton basándose en las siguientes propiedades:

Si \(b\) es una aproximación para el punto cero de \(f\), entonces \[ b-\frac{f(b)}{f'(b)} \] donde \(f'\) es la derivada de \(f\), es una mejor aproximación.
el límite de la sucesión \(x_n\) definida por \begin{align} x_0 &= 1 \newline x_{n+1} &= x_n-\frac{f(x_n)}{f'(x_n)} \end{align} es un cero de \(f\).

Definir la función

   puntoCero :: (Double -> Double) -> Double

tal que puntoCero f es un cero de la función f calculado usando la propiedad anterior. Por ejemplo,

   puntoCero cos  ==  1.5707963267949576

Soluciones

A continuación se muestran las soluciones en Haskell y las soluciones en Python.

Soluciones en Haskell

import Test.Hspec (Spec, hspec, it, shouldBe)

-- 1ª solución
-- ===========

puntoCero :: (Double -> Double) -> Double
puntoCero f = puntoCeroAux f 1
  where puntoCeroAux f' x | aceptable x = x
                          | otherwise   = puntoCeroAux f' (mejora x)
        aceptable b = abs (f b) < 0.00001
        mejora b    = b - f b / derivada f b

-- (derivada f x) es el valor de la derivada de la función f en el punto
-- x con aproximación 0.0001. Por ejemplo,
--    derivada sin pi == -0.9999999983354435
--    derivada cos pi == 4.999999969612645e-5
derivada :: (Double -> Double) -> Double -> Double
derivada f x = (f (x+a) - f x)/a
  where a = 0.0001

-- 2ª solución
-- ===========

puntoCero2 :: (Double -> Double) -> Double
puntoCero2 f = until aceptable mejora 1
  where aceptable b = abs (f b) < 0.00001
        mejora b    = b - f b / derivada f b

-- Verificación
-- ============

verifica :: IO ()
verifica = hspec spec

spec :: Spec
spec = do
  it "e1" $
    puntoCero cos `shouldBe` 1.5707963267949576
  it "e2" $
    puntoCero2 cos `shouldBe` 1.5707963267949576

-- La verificación es
--    λ> verifica
--
--    e1
--    e2
--
--    Finished in 0.0002 seconds
--    2 examples, 0 failures

Soluciones en Python

from math import cos, pi, sin
from typing import Callable

# 1ª solución
# ===========

# derivada(f, x) es el valor de la derivada de la función f en el punto
# x con aproximación 0.0001. Por ejemplo,
#    derivada(sin, pi) == -0.9999999983354435
#    derivada(cos, pi) == 4.999999969612645e-5
def derivada(f: Callable[[float], float], x: float) -> float:
    a = 0.0001
    return (f(x+a) - f(x)) / a

def puntoCero(f: Callable[[float], float]) -> float:
    def aceptable(b: float) -> bool:
        return abs(f(b)) < 0.00001
    def mejora(b: float) -> float:
        return b - f(b) / derivada(f, b)
    def aux(g: Callable[[float], float], x: float) -> float:
        if aceptable(x):
            return x
        return aux(g, mejora(x))
    return aux(f, 1)

# 2ª solución
# ===========

def puntoCero2(f: Callable[[float], float]) -> float:
    def aceptable(b: float) -> bool:
        return abs(f(b)) < 0.00001
    def mejora(b: float) -> float:
        return b - f(b) / derivada(f, b)
    y = 1.0
    while not aceptable(y):
        y = mejora(y)
    return y

# Verificación
# ============

def test_puntoCero () -> None:
    assert puntoCero(cos) == 1.5707963267949576
    assert puntoCero(cos) - pi/2 == 6.106226635438361e-14
    assert puntoCero(sin) == -5.8094940533562345e-15
    assert puntoCero2(cos) == 1.5707963267949576
    assert puntoCero2(cos) - pi/2 == 6.106226635438361e-14
    assert puntoCero2(sin) == -5.8094940533562345e-15
    print("Verificado")

# La comprobación es
#    >>> test_puntoCero()
#    Verificado