In [35]:
def halleys(x0, f, J, H, abs_tol=1e-15, nl_max_its=50):
    x = x0
    fe = f(x)
    it = 0    
    print("Iteration: %d, x: %.1e, f: %.1e" % (it, x, fe))
    

    while abs(fe) > abs_tol and it < nl_max_its:
        Je = J(x)
        He = H(x)
        x = x - fe / Je * 1 / (1 - fe / Je * (He / (2 * Je)))
        fe = f(x)
        it += 1
        print("Iteration: %d, x: %.1e, f: %.1e" % (it, x, fe))

In [36]:
def newtons(x0, f, J, abs_tol=1e-15, nl_max_its=50, m=1):
    x = x0
    fe = f(x)
    it = 0    
    print("Iteration: %d, x: %.1e, f: %.1e" % (it, x, fe))
    

    while abs(fe) > abs_tol and it < nl_max_its:
        Je = J(x)
        x = x - m * fe / Je
        fe = f(x)
        it += 1
        print("Iteration: %d, x: %.1e, f: %.1e" % (it, x, fe))

In [25]:
def f(x):
    return x**2
def J(x):
    return 2*x
def H(x):
    return 2

# from math import *
# def f(x):
#     return sin(x)
# def J(x):
#     return cos(x)
# def H(x):
#     return -sin(x)

In [38]:
halleys(1, f, J, H)


Iteration: 0, x: 1.0e+00, f: 1.0e+00
Iteration: 1, x: 3.3e-01, f: 1.1e-01
Iteration: 2, x: 1.1e-01, f: 1.2e-02
Iteration: 3, x: 3.7e-02, f: 1.4e-03
Iteration: 4, x: 1.2e-02, f: 1.5e-04
Iteration: 5, x: 4.1e-03, f: 1.7e-05
Iteration: 6, x: 1.4e-03, f: 1.9e-06
Iteration: 7, x: 4.6e-04, f: 2.1e-07
Iteration: 8, x: 1.5e-04, f: 2.3e-08
Iteration: 9, x: 5.1e-05, f: 2.6e-09
Iteration: 10, x: 1.7e-05, f: 2.9e-10
Iteration: 11, x: 5.6e-06, f: 3.2e-11
Iteration: 12, x: 1.9e-06, f: 3.5e-12
Iteration: 13, x: 6.3e-07, f: 3.9e-13
Iteration: 14, x: 2.1e-07, f: 4.4e-14
Iteration: 15, x: 7.0e-08, f: 4.9e-15
Iteration: 16, x: 2.3e-08, f: 5.4e-16

In [39]:
newtons(1, f, J, m=2)


Iteration: 0, x: 1.0e+00, f: 1.0e+00
Iteration: 1, x: 0.0e+00, f: 0.0e+00

In [ ]: