优化方法之最速下降法(Python实现)

最速下降法是优化方法中很重要的一类基础搜素算法,其具体算法如下:
在这里插入图片描述
在接下来的算法实现中,我们使用二阶微分的方法求解最优步长(当然也可以用一维搜索算法进行实现),直接取函数的一阶偏导作为方向向量.
具体实现代码如下:

from sympy import*
import math
import matplotlib.pyplot as plt
import numpy as np

def ObjFunction(x1,x2):
    value =2*x1*x1+x2*x2
    return value
def Square(x1,x2):
    return math.sqrt(math.pow(x1,2)+math.pow(x2,2))
x1,x2=sy.symbols('x1,x2')
Obj=2*x1*x1+x2*x2
z=diff(Obj,x1,x1)
print(z)
result=z.subs({x1:1,x2:3})#替换数据
print(result)

#指定初始值
Epsilon=0.1
x=np.array([1,1])
df1=diff(Obj,x1)
df2=diff(Obj,x2)
Jacobian=np.array([df1,df2])
Hessian=([[0,0],[0,0]])
for i in range (0,2,1): #求解Hessian矩阵
    for j in range(0,2,1):
        if i==0:
            Hessian[i][j]=diff(Jacobian[i],x1)
        else:
            Hessian[i][j]=diff(Jacobian[i],x2)
Direction=np.array([-Jacobian[0].subs({x1:x[0],x2:x[1]}),-Jacobian[1].subs({x1:x[0],x2:x[1]})]) #求解方向向量的初始向量
print(Direction)
#计算Lambda的初始值
Temp_Jacobian=np.array([0,0])
for i in range(0,2,1):
    Temp_Jacobian[i]=Jacobian[i].subs({x1:x[0],x2:x[1]})
Temp_Hessian=np.array(([0,0],[0,0]))
for i in range(0,2,1): #求解Hessian矩阵的实例化
    for j in range(0,2,1):
        Temp_Hessian[i][j]=Hessian[i][j].subs({x1:x[0],x2:x[1]})
Temp=np.array([0,0])
Temp[0]=Temp_Jacobian[0]*Hessian[0][0]+Temp_Jacobian[1]*Hessian[1][0]
Temp[1]=Temp_Jacobian[1]*Hessian[0][1]+Temp_Jacobian[1]*Hessian[1][1]
Lambda=np.dot(Temp_Jacobian.T,Temp_Jacobian)/np.dot(Temp,Temp_Jacobian) #求解Lambda示例的初始化

#使用最速下降法求解目标函数最优值
result=[] #存储计算的函数值
Minimum=[0,0] #求解最优值
while Square(Direction[0],Direction[1])<Epsilon:
    #更新x的值
    result.append(ObjFunction(x[0],x[1]))
    Minimum=x[0]
    Minimum=x[1]
    print(Minimum)
    x[0]=x[0]+Direction[0]*Lambda
    x[1]=x[1]+Direction[1]*Lambda
    
    #更新Direction的值
    Direction=np.array([-Jacobian[0].subs({x1:x[0],x2:x[1]}),-Jacobian[1].subs({x1:x[0],x2:x[1]})]) 
    
    #更新Lambda
    for i in range(0,2,1):
        Temp_Jacobian[i]=Jacobian[i].subs({x1:x[0],x2:x[1]})
        Temp_Hessian=np.array(([0,0],[0,0]))
    for i in range(0,2,1): #求解Hessian矩阵的实例化
        for j in range(0,2,1):
            Temp_Hessian[i][j]=Hessian[i][j].subs({x1:x[0],x2:x[1]})
    Temp=np.array([0,0])
    Temp[0]=Temp_Jacobian[0]*Hessian[0][0]+Temp_Jacobian[1]*Hessian[1][0]
    Temp[1]=Temp_Jacobian[1]*Hessian[0][1]+Temp_Jacobian[1]*Hessian[1][1]
    temp1=(Temp_Jacobian[0]*Temp_Jacobian[0]+Temp_Jacobian[1]*Temp_Jacobian[1])
    temp2=(Temp[0]*Temp_Jacobian[0]+Temp[1]*Temp_Jacobian[1])
    Lambda=temp1/temp2
print(Lambda)
print(ObjFunction(x[0],x[1]))

    


在这里插入图片描述
运行结果如上所示。

评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值