22、DFP法

最新推荐文章于 2023-03-31 10:41:35 发布

chaoyuzhang

最新推荐文章于 2023-03-31 10:41:35 发布

阅读量404

点赞数

分类专栏：运筹优化、机器学习基础算法实现

本文链接：https://blog.csdn.net/chaoyuzhang/article/details/88747886

版权

运筹优化、机器学习基础算法实现专栏收录该内容

42 篇文章 9 订阅

订阅专栏

# -*- coding: utf-8 -*-
"""
Created on Fri Mar 22 17:07:02 2019

@author: zhangchaoyu
"""
import math
import copy

"""人为修改函数，梯度函数、步长迭代函数"""
"""
1、f(x,y) = 1.5x*x + 0.5y*y - xy - 2x
2、f(x,y) = 4(x-5)(x-5) + (y-6)(y-6)
"""
def gradient(X):
    x = X[0]
    y = X[1]
    #return [3*x-y-2, y-x]
    return [8*x-40, 2*y-12]

def step_best(X, P):
    px = P[0]
    py = P[1]
    x = X[0]
    y = X[1]
    #return -(px*(3*x-y-2)+py*(y-x))/(3*px*px-2*px*py+py*py)
    return -(px*(4*x-20)+py*(y-6))/(4*px*px+py*py)
def multiplication_M_and_V(H, G): 
    P = []
    for i in range(len(H)):
        P.append(-sum([H[i][j]*G[j] for j in range(len(H[i]))]))
    return P

def multiplication_M_and_M(A, B): 
    C = []
    for i in range(len(A)):
        temp = []
        for j in range(len(B[0])):
            temp.append(sum([A[i][k]*B[k][j] for k in range(len(A[i]))]))
        C.append(copy.deepcopy(temp))
    return C

def trans(A):
    return [[A[j][i] for j in range(len(A))] for i in range(len(A[0]))]

def H_iteration(H, dX, dG):
    H1 = copy.deepcopy(H)
    
    H2 = multiplication_M_and_M(trans([dX]), [dX])
    c2 = sum([dX[i]*dG[i] for i in range(len(dX))])
    H2 = [[H2[i][j]/c2 for j in range(len(H2[i]))] for i in range(len(H2))]
    
    c3 = multiplication_M_and_M(multiplication_M_and_M([dG], H), trans([dG]))[0][0]
    H3 = multiplication_M_and_M(multiplication_M_and_M(H, trans([dG])), trans(multiplication_M_and_M(H, trans([dG]))))
    H3 = [[H3[i][j]/c3 for j in range(len(H3[i]))] for i in range(len(H3))]
    
    H4 = [[H1[i][j]+H2[i][j]-H3[i][j] for j in range(len(H1[i]))] for i in range(len(H1))]
    
    return H4

def Newton(X0):    
    
    X = X0
    G = gradient(X)
    H = [[0 for j in range(len(X))] for i in range(len(X))]
    for i in range(len(H)):
        H[i][i] = 1
    P = multiplication_M_and_V(H, G)
    
    while math.sqrt(G[0]*G[0]+G[1]*G[1]) > 1e-5:
        
        step = step_best(X, P)
        dX = [step*P[i] for i in range(len(X))]
        
        X = [X[i]+dX[i] for i in range(len(X))]
        G1 = gradient(X)
        dG = [G1[i]-G[i] for i in range(len(G))]
        H = H_iteration(H, dX, dG)
        G = copy.deepcopy(G1)
        P = multiplication_M_and_V(H, G)
        print(X)
    
    return X