Newton-Raphson法

Newton-Raphson法

用于求解方程的根的非常重要且经典的数值计算方法

1. 从例子体会牛顿法的精髓

以$f(x)=x^2-4$ 为例，函数曲线如下图

1）从$x_0=4$ 开始迭代计算，此时曲线上的坐标为 $(4,12)$。

在$(4,12)$点处，该函数的切线为 $f(x)=8x-20$

这时，切线$f(x)=8x-20$ 与 $x$轴的交点坐标为$(2.5,0)$。

2）令$x_1=2.5$，此时曲线上的坐标为 $(2.5,2.25)$。

在$(2.5,2.25)$点处，该函数的切线为 $f(x)=5x-10.25$

这时，切线$f(x)=5x-10.25$与 $x$轴的交点坐标为$(2.05,0)$。

3）$\cdots$

直到 $|x_n-x_{n-1}|<ϵ$，停止迭代（也可以是相对距离度量）。

2. 抽象牛顿法本质

该方法使用函数$f(x)$的泰勒级数的前2项近似$f(x)$求解$f(x)=0$的根。将$f(x)$函数在点$x_0$（初始点或已知点）的某邻域内展开成$n$ 阶泰勒公式，如下
$$f(x)=f(x_0)+f^{\prime }(x_0)(x-x_0)+\frac{f^{\prime \prime }(x_0)}{2!}(x-x_0{)}^2+\cdots +R_n(x)$$
其中 $R_n(x)$ 为$n$ 阶泰勒余项。

因此，对任意方程$f(x)=0$，牛顿法的迭代过程可以描述为
$$x_{n+1}=x_n-\frac{f(x_n)}{f^{\prime }(x_n)}$$

3. Python 代码实现

import matplotlib.pyplot as plt
import numpy as np


def function(x):
    """
    Function of the function.

        Parameters:
        -----------
            x : int | float
                var

        Returns:
        --------
            value : int | float
                calculated value of the function
    """
    value = x ** 2 - 4
    return value  # The main function


def derivative(x):
    """
        Derivative of the function.

            Parameters:
            -----------
                x : int | float 
                    var

            Returns:
            --------
                value : int | float
                    calculated value of the derivative of the function
        """
    return 2 * x  # The derivative of the main function

def newton(function, derivative, x0, tol, max_iter=100):
    """
    Calculate and plots the values througth newton-raphson method.

        Parameters:
        -----------
            function : function
                function
            derivative : function
                derivative
            x0: int | float
                x0
            tol: float
                tolerance
            max_iter: int
                number of the maximal iterations

        Returns:
        --------
            x1 : int | float
                calculated value of the function
    """
    x1 = 0

    if abs(x0-x1)<= tol and abs((x0-x1)/x0)<= tol:
        return x0

    print("k\t x0\t\t function(x0)")
    k = 1

    while k <= max_iter:
        x1 = x0 - (function(x0)/derivative(x0))
        print("x%d\t%e\t%e"%(k,x1,function(x1)))

        if abs(x0-x1)<= tol and abs((x0-x1)/x0)<= tol:
            plt.plot(x0, function(x0), '.')
            return x1

        x0 = x1
        k = k + 1
        plt.plot(x0, function(x0), '.')

        # Stops the method if it hits the number of maximum iterations
        if k > max_iter:
            print("ERROR: Exceeded max number of iterations")

    return x1  # Returns the value

sqrt = newton(function, derivative,9, 0.0000001)
print("The approximate value of x is: "+str(sqrt))

# Plotting configuration
u = np.arange(-10, 10, 0.0001) # Setting up values for x in the plot
w = u**2 - 4 # Define the main function again

plt.plot(u, w)
plt.axhline(y=0.0, color='black', linestyle='-')
plt.title('Newton-Raphson Graphics for' + ' y = x^2 - 2')
plt.xlabel('X')
plt.ylabel('Y')
plt.grid(True)
plt.legend(['Xn'], loc='upper left')
plt.show()

参考资料

[1] 机器学习 | Newton-Raphson法

[2] Newton-Raphson Method

[3] Toktom/Newton-Raphson-Method: Newton-Raphson Method Algorithm in Python