# Newton's method # Second order method

newton_method.py

import numpy as np
# ################################################################################
# Newton's method
# Second order method
def Newtons_method(max_iters, threshold, XY_init, func, grad_func, second_grad, learning_rate=0.05):
    X, Y = XY_init
    w = np.array([X, Y])
    w_history_N = X, Y
    f_history_N = func(X, Y, extra_param)
    delta_w = np.zeros(XY_init.shape)
    i = 0
    # start diff_f |f2 - f1| for stop criteria
    diff_f = 1.0e10
    eps_history_N = np.array([0.0])
    eps_history_xy_N = np.array([0.0 , 0.0])

    while i < max_iters and diff_f >= threshold:
        # print(grad_func(w[0], w[1], extra_param))
        # print(np.linalg.inv(second_grad(w[0], w[1], extra_param)))
        # Inverse Hessian matrix of second derivatives.
        inverseHessian = np.linalg.inv(second_grad(w[0], w[1], extra_param))
        # Hessian step calculation, it is not define descent direction.
        # The Hessian matrix characterizes the convexity or concavity of a
        # surface and can change the sign of the determinant at the inflection point.
        inverseHessian = np.abs(inverseHessian)
        delta_w = - learning_rate * np.dot(grad_func(w[0], w[1], extra_param), inverseHessian)

        w = w + delta_w
        # store the history of w and f
        w_history_N = np.vstack((w_history_N, w))
        f_history_N = np.vstack((f_history_N, func(w[0], w[1], extra_param)))

        i += 1
        diff_f = np.absolute(func(3.0, 2.0,extra_param) - f_history_N[-1])
        diff_xy = np.absolute(w_history_N[-1] - w_history_N[-2])

        eps_history_N = np.vstack((eps_history_N, diff_f))
        eps_history_xy_N= np.vstack((eps_history_xy_N, diff_xy))
        
    return w_history_N, f_history_N, eps_history_N, eps_history_xy_N
# ################################################################################

• set as public gist