牛顿法的推导

最新推荐文章于 2023-12-04 19:00:33 发布

西蒙先生、

最新推荐文章于 2023-12-04 19:00:33 发布

阅读量2.2k

点赞数 2

分类专栏：机器学习文章标签：机器学习数学

本文链接：https://blog.csdn.net/qq_35661436/article/details/103798255

版权

机器学习专栏收录该内容

3 篇文章 0 订阅

订阅专栏

一、使用牛顿法可以用来求解的根。

具体原理如下：

首先初始化一个点为 (x_0,y_0) ，然后求解 f(x) 在 x_0 处的切线，设切线为 h_0(x) ， x_0 处的斜率为 f{}'(x) ，则切线为 $h_0(x) = f{}'(x)(x-x_0)+y_0$ 。

此切线与x轴交于 (x_1,0) ，此时在x轴上的点为 (x_1,y_1) ，然后在处求切线 $h_1(x) = f{}'(x)(x-x_1)+y_1$ ，交x轴为 (x_2,0) 。

不断重复此步骤最终得到 f(x)=0 的解。

因为我们想要求解 f(x)=0 的根。因此我们让 $h_(_i_)(x) = 0$ 。可以得到得到此公式 $x_{i+1} = x_{(i)} - \frac{f(x_{i})}{f{}'(x_{i})}$

由上式可以得到简化的迭代公式 $x_{i+1} := x_{(i)} - \frac{f(x_{i})}{f{}'(x_{i})}$ ，使用此公式不断的进行迭代。然后通过设置阈值得到数值解。

求此函数 $f(x) = x^{2}+2x+0.2$

具体代码如下：

#实现f = x**2+2*x - 0.2的根
import numpy as np

#定义fx
def f(x):
    f = x**2+2*x - 0.2
    return f
#定义导数
def df(x):
    df = 2*x+2
    return df
# f:函数
# df:导数
# x：zip变量
# err：误差，用来判断是否满足求解要求
def newtons(f, df, x, err):
    loop = 1
    x_i = float(x)
    e_tmp = err + 1
    while e_tmp > err:
        #超出循环次数
        if loop > 1000:
            print('cycles exceeded')
            break
        print('######loop'+str(loop))
        f_tmp = f(x_i)
        df_tmp = df(x_i)
        print('xi = ' + str(x_i) + ',f = ' + str(f_tmp) + ',df = ' + str(df_tmp))
        x_i = x_i - (f_tmp/df_tmp)
        e_tmp = abs(0 - f(x_i))
        print('err = ' + str(e_tmp))
        loop = loop + 1
    return x_i
x = newtons(f, df, 3, 0.01)
print(x)

二、牛顿法也被用于求函数的极值。

由于函数取极值的点处的导数值为零，故可用牛顿法求导函数的零点，其迭代式为

$x_{i+1}=x_{i}-{\frac {f^{\prime }(x_{i})}{f^{\prime \prime }(x_{i})}}.$

具体代码如下：

#实现f = x**2+2*x - 0.2的最小值
import numpy as np

#定义fx
def f(x):
    f = x**2+2*x - 0.2
    return f
#定义导数
def df(x):
    df = 2*x+2
    return df
#定义二阶导数
def ddf(x):
    df = 2
    return df
# f:函数
# df:导数
# x：zip变量
# err：误差，用来判断是否满足求解要求
def newtons(df, ddf, x, err):
    loop = 1
    x_i = float(x)
    e_tmp = err + 1
    while e_tmp > err:
        #超出循环次数
        if loop > 1000:
            print('cycles exceeded')
            break
        print('######loop'+str(loop))
        df_tmp = df(x_i)
        ddf_tmp = ddf(x_i)
        print('xi = ' + str(x_i) + ',df = ' + str(df_tmp) + ',ddf = ' + str(ddf_tmp))
        x_i = x_i - (df_tmp/ddf_tmp)
        e_tmp = abs(df(x_i) - df_tmp)
        print('err = ' + str(e_tmp))
        loop = loop + 1
    return x_i
x = newtons(df, ddf, 3, 0.01)
print(x)

多元函数：

import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
import numpy as np
 
def Fun(x,y):#原函数
    return x-y+2*x*x+2*x*y+y*y
 
def PxFun(x,y):#偏x导
    return 1+4*x+2*y
 
def PyFun(x,y):#偏y导
    return -1+2*x+2*y
#二次偏x导
def PPxFun(x,y):
    return 4
#二次偏y导
def PPyFun(x,y):
    return 2
 
#初始化
fig=plt.figure()#figure对象
ax=Axes3D(fig)#Axes3D对象
X,Y=np.mgrid[-2:2:40j,-2:2:40j]#取样并作满射联合
Z=Fun(X,Y)#取样点Z坐标打表
ax.plot_surface(X,Y,Z,rstride=1,cstride=1,cmap="rainbow")
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_zlabel('z')
 
#梯度下降
x=0
y=0#初始选取一个点
tag_x=[x]
tag_y=[y]
tag_z=[Fun(x,y)]#三个坐标分别打入表中，该表用于绘制点
new_x=x
new_y=y
Over=False
while Over==False:
    new_x-=PxFun(x,y)/PPxFun(x,y)
    new_y-=PyFun(x,y)/PPyFun(x,y)#分别作梯度下降
    if Fun(x,y)-Fun(new_x,new_y)<7e-9:#精度
        Over=True
    x=new_x
    y=new_y#更新旧点
    tag_x.append(x)
    tag_y.append(y)
    tag_z.append(Fun(x,y))#新点三个坐标打入表中
 
#绘制点/输出坐标
ax.plot(tag_x,tag_y,tag_z,'r.')
plt.title('(x,y)~('+str(x)+","+str(y)+')')
plt.show()